/*File: lennard-jones_ga.c Purpose: to test Lennard-Jones problem using Tao/GA vectors using multiple processors. Date created: August 16, 2002 */ /************************************************************** Author: Limin Zhang, Ph.D. Mathematics Department Columbia Basin College Pasco, WA 99301 Limin.Zhang@cbc2.org Mentor: Jarek Naplocha, Ph.D. Environmental Molecular Science Laboratory Pacific Northwest National Laboratory Richland, WA 99352 Date: 4/22/2002 Purpose: to design and implement an application for between TAO and global arrays. **************************************************************/ /* Include MPI header files */ #include <mpi.h> /* Program usage: mpirun -machinefile m -np 1 lennard-jones_ga [-help] [all TAO options] */ /* Program usage: mpirun -machinefile m -np 4 lennard-jones_ga [-help] [all TAO options] */ /* Include tao header files so we can use TAO */ #include "tao.h" #include "taoapp_ga.h" /* Include ga header files so that we can use GA */ #include "ga.h" #include "taovec_ga.h" #include "macdecls.h" /*some standard header files */ #include <stdio.h> #include <math.h> #define NDIM 2 #define NATOMS 16 static char help[] = "This example demonstrates use of the TAO package to \n\ solve an unconstrained minimization problem on multiple processors. We \n\ minimize the extended Lennard-Jones function."; /*T Concepts: TAO - Solving an unconstrained minimization problem Routines: TaoInitialize(); TaoFinalize(); TaoSetFromOptions(); Routines: TaoGAApplicationCreate(); TaoGASetInitialSolutionVec(); Routines: TaoCreate(); TaoGASetObjectiveAndGradientRoutine(); Routines: TaoSetGAHessian(); TaoSetGAInitialVector(); Routines: TaoSolveGAApplication(); TaoDestroy(); TaoGAAppDestroy(); Routines: TaoGetTerminationReason(); TaoView(); Processors: 1 T*/ /* User-defined application context - contains data needed by the application-provided call-back routines that evaluate the function, gradient, and hessian. */ typedef struct { int n; /* dimension n = NDIM*NATOMS */ int ndim; int natoms; int memHandle; } AppCtx; /* -------------- User-defined routines ---------- */ int FormFunctionGradient (TAO_GA_APPLICATION gaapp, GAVec ga_X, double *f, GAVec ga_G, void *ptr); int InitializeVariables(GAVec ga_X, AppCtx *user); int main (int argc, char **argv) { double startTime; int info; /* used to check for functions returning nonzeros */ GAVec ga_x; /* solution vector */ TAO_SOLVER tao; /* TAO_SOLVER solver context */ TAO_GA_APPLICATION taoapp; /* TAO application context */ TaoTerminateReason reason; AppCtx user; /* user-defined application context */ /*initialize GA and MPI */ int heap = 4000, stack = 4000; MPI_Init (&argc, &argv); /* initialize MPI */ GA_Initialize (); /* initialize GA */ if (!MA_init(MT_F_DBL, stack, heap)) GA_Error("MA_init failed", stack+heap); /* Initialize TAO */ TaoInitialize (&argc, &argv, (char *) 0, help); startTime = MPI_Wtime(); /* Initialize problem parameters */ user.natoms = NATOMS; user.ndim = NDIM; user.n = user.natoms*user.ndim; /* Create working space */ if (MA_push_stack(C_DBL, 2*user.n, "Vector buffers", &user.memHandle) == MA_FALSE) ga_error("MAIN::ma_alloc_get failed",2*user.n); /* Allocate Global Array vector for the solution */ int dims[2]; dims[0] = user.n; ga_x = NGA_Create (C_DBL, 1, dims, "GA_X", NULL); if (!ga_x) GA_Error ("lennard-jones.main::NGA_Create ga_x", ga_x); /* The TAO code begins here */ /* Create TAO solver with desired solution method */ info = TaoCreate (MPI_COMM_WORLD, "tao_cg_fr", &tao); CHKERRQ(info); info = TaoGAApplicationCreate (MPI_COMM_WORLD, &taoapp); CHKERRQ(info); /* Set initial vector */ info = InitializeVariables(ga_x, &user); CHKERRQ(info); info = TaoGAAppSetInitialSolutionVec(taoapp, ga_x); CHKERRQ(info); /* Set routines for function, gradient */ info = TaoGAAppSetObjectiveAndGradientRoutine (taoapp, FormFunctionGradient, (void *) &user); CHKERRQ(info); /* Check for TAO command line options */ info = TaoSetFromOptions (tao); CHKERRQ(info); /* SOLVE THE APPLICATION */ info = TaoSolveGAApplication (taoapp, tao); CHKERRQ(info); /* To View TAO solver information use */ info = TaoView(tao); CHKERRQ(info); /* Get termination information */ info = TaoGetTerminationReason (tao, &reason); CHKERRQ(info); if (reason <= 0) printf("Try a different TAO method, adjust some parameters, or check the function evaluation routines\n"); printf("TIME TAKEN = %lf\n", MPI_Wtime()-startTime); /*output the solutions */ printf ("The solution is :\n"); GA_Print (ga_x); /* Free TAO data structures */ info = TaoDestroy (tao); CHKERRQ(info); info = TaoGAAppDestroy (taoapp); CHKERRQ(info); /* Free GA data structures */ GA_Destroy (ga_x); if (!MA_pop_stack(user.memHandle)) ga_error("Main::MA_pop_stack failed",0); /* Finalize TAO, GA, and MPI */ TaoFinalize (); GA_Terminate (); MPI_Finalize (); return 0; } /* -------------------------------------------------------------------- */ /* FormFunctionGradient - Evaluates the function, f(X), and gradient, G(X). Input Parameters: . tao - the TAO_SOLVER context . ga_X - input vector . ptr - optional user-defined context, as set by TaoSetFunctionGradient() Output Parameters: . ga_G - vector containing the newly evaluated gradient . f - function value Note: Some optimization methods ask for the function and the gradient evaluation at the same time. Evaluating both at once may be more efficient that evaluating each separately. */ int FormFunctionGradient (TAO_GA_APPLICATION gaapp, GAVec ga_X, double *f, GAVec ga_G, void *ptr) { int lo, hi; //the global coordinates AppCtx *user = (AppCtx *) ptr; int i,j; double *g, *x; double xx,yy,zz,temp,rij; MA_get_pointer(user->memHandle, &x); g = x + user->ndim*user->natoms; lo=0; hi=user->n-1; /* range of array indices */ NGA_Get(ga_X, &lo, &hi, x, &hi); *f = 0; for (i=0; i < user->n; i++) g[i] = 0.0; if (user->ndim == 2) { for (j=1; j < user->natoms; j++) { for (i=0; i<j; i++) { xx = x[2*j] - x[2*i]; yy = x[2*j+1] - x[2*i+1]; rij = xx*xx + yy*yy; temp = 1.0/rij/rij/rij; *f += temp*(temp-2.0); temp *= 12.0*(temp-1.0)/rij; g[2*j] -= xx*temp; g[2*j+1] -= yy*temp; g[2*i] += xx*temp; g[2*i+1] += yy*temp; } } } else if (user->ndim == 3) { for (j=1; j < user->natoms; j++) { for (i=0; i < j; i++) { xx = x[3*j] - x[3*i]; yy = x[3*j+1] - x[3*i+1]; zz = x[3*j+2] - x[3*i+2]; rij = xx*xx + yy*yy + zz*zz; temp = 1.0/rij/rij/rij; *f += temp*(temp-2.0); temp *= 12.0*(temp-1.0)/rij; g[3*j] -= xx*temp; g[3*j+1] -= yy*temp; g[3*j+2] -= zz*temp; g[3*i] += xx*temp; g[3*i+1] += yy*temp; g[3*i+2] += zz*temp; } } } NGA_Put(ga_G, &lo, &hi, g, &hi); return 0; } InitializeVariables">int InitializeVariables(GAVec ga_X, AppCtx *user) { double *x; double xx, yy, zz; int isqrtn, icrtn, left, i, j, k, il, jl, ctr; int lo, hi; lo = 0; hi = user->n - 1; /* range of array indices */ MA_get_pointer(user->memHandle, &x); NGA_Get (ga_X, &lo, &hi, x, &hi); if (user->ndim == 2) { isqrtn = (int) sqrt( (double) user->natoms); left = user->natoms - isqrtn * isqrtn; xx = 0.0; yy = 0.0; for (j=0; j<=isqrtn + left/isqrtn; j++) { for (i=0; i < TaoMin(isqrtn, user->natoms - j*isqrtn); i++) { ctr = j*isqrtn + i; x[2*ctr] = xx; x[2*ctr+1] = yy; xx += 1.0; } yy += 1.0; xx = 0.0; } } else if (user->ndim == 3) { icrtn = (int) pow((user->natoms + 0.5),1.0/3.0); left = user->natoms - icrtn * icrtn * icrtn; xx = yy = zz = 0.0; for (k=0; k <= icrtn + left; k++) { jl = TaoMin(icrtn, (user->natoms - k*icrtn*icrtn)/icrtn+1); for (j=0; j<jl; j++) { il = TaoMin(icrtn, user->natoms - k*icrtn*icrtn - j*icrtn); for (i=0; i<il; i++) { ctr = k*icrtn*icrtn + j*icrtn + i; x[3*ctr] = xx; x[3*ctr+1] = yy; x[3*ctr+2] = zz; xx += 1.0; } yy += 1.0; xx = 0.0; } zz += 1.0; yy = 0.0; } } NGA_Put(ga_X, &lo, &hi, x, &hi); return 0; }