Actual source code: petsckernel.c

  1: #include "tao_general.h"  /*I "tao_general.h"  I*/
  2: #ifdef TAO_USE_PETSC

  4: #include "src/tao_impl.h"

  6: TaoTruth TaoBeganPetsc = TAO_FALSE;


 10: #include "petscsys.h"

 12: static PetscFList TaoList = 0;

 16: /*
 17:    TaoPrintVersion - Prints TAO version info.

 19:    Collective on MPI_Comm
 20: */
 21: int TaoPrintVersion(MPI_Comm comm)
 22: {
 23:   int info=0;
 24: 

 27:   info = PetscHelpPrintf(comm,"--------------------------------------------\
 28: ------------------------------\n"); CHKERRQ(info);
 29:   info = PetscHelpPrintf(comm,"\t   %s\n",TAO_VERSION_NUMBER); CHKERRQ(info);
 30:   info = PetscHelpPrintf(comm,"%s",TAO_AUTHOR_INFO); CHKERRQ(info);
 31:   info = PetscHelpPrintf(comm,"See docs/manualpages/index.html for help. \n"); CHKERRQ(info);
 32: #if !defined(PARCH_win32)
 33:   info = PetscHelpPrintf(comm,"TAO libraries linked from %s\n",TAO_LIB_DIR); CHKERRQ(info);
 34: #endif
 35:   info = PetscHelpPrintf(comm,"--------------------------------------------\
 36: ------------------------------\n"); CHKERRQ(info);

 38:   PetscFunctionReturn(info);
 39: }

 43: /*
 44:    TaoPrintHelpIntro - Prints introductory TAO help info.

 46:    Collective on MPI_Comm
 47: */
 48: int TaoPrintHelpIntro(MPI_Comm comm)
 49: {
 50:   int info=0;
 51: 

 54:   info = PetscHelpPrintf(comm,"--------------------------------------------\
 55: ------------------------------\n"); CHKERRQ(info);
 56:   info = PetscHelpPrintf(comm,"TAO help information includes that for the PETSc libraries, which provide\n"); CHKERRQ(info);
 57:   info = PetscHelpPrintf(comm,"low-level system infrastructure and linear algebra tools.\n"); CHKERRQ(info);
 58:   info = PetscHelpPrintf(comm,"--------------------------------------------\
 59: ------------------------------\n"); CHKERRQ(info);

 61:   PetscFunctionReturn(info);
 62: }

 66: int TaoPrintStatement(TAO_SOLVER tao, const char *statement){
 69:   PetscPrintf(tao->comm,statement);
 70:   return(0);
 71: }

 75: int TaoPrintInt(TAO_SOLVER tao, const char *statement, int n){
 78:   PetscPrintf(tao->comm,statement,n);
 79:   return(0);
 80: }

 84: int TaoPrintDouble(TAO_SOLVER tao, const char *statement,double dd){
 87:   PetscPrintf(tao->comm,statement,dd);
 88:   return(0);
 89: }

 93: int TaoPrintString(TAO_SOLVER tao, const char *statement,const char *str){
 96:   PetscPrintf(tao->comm,statement,str);
 97:   return(0);
 98: }

102: int TaoOptionsHead(const char *heading){
103:   int info;
105:   info = PetscOptionsHead(heading);CHKERRQ(info);
106:   return(0);
107: }

111: int TaoOptionsTail(){
112:   int info;
114:   info = PetscOptionsTail();CHKERRQ(info);
115:   return(0);
116: }

120: int TaoOptionInt(const char *opt,const char *text,const char *man,int defaultv,int *value,TaoTruth *set){
121:   int info;
122:   PetscTruth flg=PETSC_FALSE;
124:   info = PetscOptionsInt(opt,text,man,defaultv,value,&flg);CHKERRQ(info);
125:   if (set){
126:     if (flg==PETSC_TRUE) *set=TAO_TRUE; else *set=TAO_FALSE;
127:   }
128:   return(0);
129: }

133: int TaoOptionDouble(const char *opt,const char *text,const char *man,double defaultv,double *value,TaoTruth *set){
134:   int info;
135:   PetscReal pdv=(PetscReal)defaultv, pv;
136:   PetscTruth flg=PETSC_FALSE;
138:   info = PetscOptionsReal(opt,text,man,pdv,&pv,&flg);CHKERRQ(info);
139:   if (set){
140:     if (flg==PETSC_TRUE) *set=TAO_TRUE; else *set=TAO_FALSE;
141:   }
142:   if (value&&flg==PETSC_TRUE){
143:     *value=(double)pv;
144:   }

146:   return(0);
147: }

151: int TaoOptionString(const char *opt,const char *text,const char *man,const char* defaultv,char *value, int len, TaoTruth *set){
152:   int info;
153:   PetscTruth flg=PETSC_FALSE;
155:   info = PetscOptionsString(opt,text,man,defaultv,value,len,&flg);CHKERRQ(info);
156:   if (set){
157:     if (flg==PETSC_TRUE) *set=TAO_TRUE; else *set=TAO_FALSE;
158:   }
159:   return(0);
160: }

164: int TaoOptionName(const char *opt,const char *text,const char *man,TaoTruth *set){
165:   int info;
166:   PetscTruth flg=PETSC_FALSE;
168:   info = PetscOptionsName(opt,text,man,&flg);CHKERRQ(info);
169:   if (set){
170:     if (flg==PETSC_TRUE) *set=TAO_TRUE; else *set=TAO_FALSE;
171:   }
172:   return(0);
173: }

177: int TaoMethodsList(const char *opt,const char *ltext,const char *man,const char *defaultv,char *value,int len,TaoTruth *set){
178:   int info;
179:   PetscTruth flg=PETSC_FALSE;
181:   info = PetscOptionsList(opt,ltext,man,TaoList,defaultv,value,len,&flg); CHKERRQ(info);
182:   if (set){
183:     if (flg==PETSC_TRUE) *set=TAO_TRUE; else *set=TAO_FALSE;
184:   }
185:   return(0);
186: }


191: int TaoFindSolver(TAO_SOLVER tao, TaoMethod type,  int (**r)(TAO_SOLVER) ){
192:   int info;
194:   info = PetscFListFind(tao->comm,TaoList,type,(void (**)(void))r);CHKERRQ(info);
195:   if (*r){
196:     info = PetscObjectChangeTypeName((PetscObject)tao,type);CHKERRQ(info);
197:   }
198:   return(0);
199: }

203: /*@C
204:    TaoRegisterDestroy - Frees the list of minimization solvers that were
205:    registered by TaoRegisterDynamic().

207:    Not Collective

209:    Level: advanced

211: .keywords: TAO_SOLVER, register, destroy

213: .seealso: TaoRegisterAll()
214: @*/
215: int TaoRegisterDestroy(){
216:   int info;
218:   if (TaoList) {
219:     info=PetscFListDestroy(&TaoList);CHKERRQ(info);
220:     TaoList=0;
221:   }
222:   if (TaoBeganPetsc) {
223:     info = PetscFinalize();CHKERRQ(info);
224:   }

226:   return(0);
227: }

231: /* --------------------------------------------------------------------- */
232: /*MC
233:    TaoRegisterDynamic - Adds a method to the TAO_SOLVER package for unconstrained minimization.

235:    Synopsis:
236:    TaoRegisterDynamic(char *name_solver,char *path,char *name_Create,int (*routine_Create)(TAO_SOLVER))

238:    Not collective

240:    Input Parameters:
241: +  name_solver - name of a new user-defined solver
242: .  path - path (either absolute or relative) the library containing this solver
243: .  name_Create - name of routine to Create method context
244: -  routine_Create - routine to Create method context

246:    Notes:
247:    TaoRegisterDynamic() may be called multiple times to add several user-defined solvers.

249:    If dynamic libraries are used, then the fourth input argument (routine_Create)
250:    is ignored.

252:    Environmental variables such as ${TAO_DIR}, ${PETSC_ARCH}, ${PETSC_DIR}, 
253:    and others of the form ${any_environmental_variable} occuring in pathname will be 
254:    replaced with the appropriate values used when PETSc and TAO were compiled.

256:    Sample usage:
257: .vb
258:    TaoRegisterDynamic("my_solver","/home/username/my_lib/lib/libg_c++/solaris/mylib.a",
259:                 "MySolverCreate",MySolverCreate);
260: .ve

262:    Then, your solver can be chosen with the procedural interface via
263: $     TaoSetMethod(solver,"my_solver")
264:    or at runtime via the option
265: $     -tao_method my_solver

267:    Level: advanced

269: .keywords: TAO_SOLVER, register

271: .seealso: TaoRegisterAll(), TaoRegisterDestroy()
272: M*/
276: int TaoRegisterDynamic(const char *sname,const char *path,const char *name,int (*function)(TAO_SOLVER))
277: {
278:   char fullname[256];
279:   int  info;

282:   info = PetscFListConcat(path,name,fullname); CHKERRQ(info);
283: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
284:   info = PetscFListAdd(&TaoList,sname,fullname, 0);CHKERRQ(info);
285: #else
286:   info = PetscFListAdd(&TaoList,sname,fullname, (void (*)())function);CHKERRQ(info);
287: #endif
288:   return(0);
289: }

294: /*@C
295:    TaoCompareMethod - Determines whether the TAO_SOLVER structure is of a
296:    specified type.

298:    Not Collective

300:    Input Parameter:
301: .  tao - the TAO_SOLVER solver context
302: .  type - a TAO_SOLVER solver method

304:    Output Parameter:
305: .  same - TAO_TRUE if 'tao' is of method 'type'

307:    Level: developer

309: .keywords: method
310: @*/
311: int TaoCompareMethod(TAO_SOLVER tao, TaoMethod type,TaoTruth *issame){
312:   int info;
313:   PetscTruth flag;
314:   TaoMethod currenttype;
317:   info = TaoGetMethod(tao,&currenttype);CHKERRQ(info);
318:   info = PetscStrcmp(type,currenttype,&flag);CHKERRQ(info);
319:   if (issame){
320:     if (flag==PETSC_TRUE) *issame=TAO_TRUE; else *issame=TAO_FALSE;
321:   }
322:   return(0);
323: }

327: int TaoStrcmp(const char *str1,const char *str2,TaoTruth *flag){
328:   int info;
329:   PetscTruth flg;
331:   info = PetscStrcmp(str1,str2,&flg);CHKERRQ(info);
332:   if (flg==PETSC_TRUE) *flag=TAO_TRUE; else *flag=TAO_FALSE;
333:   return(0);
334: }

338: int TaoStrcpy(char* str1,const char*str2){
339:   int info;
341:   info = PetscStrcpy(str1,str2);CHKERRQ(info);
342:   return(0);
343: }

345: #
348: static int TaoPublish_Petsc(PetscObject obj)
349: {
350: #if defined(PETSC_HAVE_AMS)
351:   TAO_SOLVER       v = (TAO_SOLVER) obj;
352:   int          info;
353: #endif

355:   TaoFunctionBegin;

357: #if defined(PETSC_HAVE_AMS)
358:   /* if it is already published then return */
359:   if (v->amem >=0 ) TaoFunctionReturn(0);

361:   info = PetscObjectPublishBaseBegin(obj);CHKERRQ(info);
362:   info = AMS_Memory_add_field((AMS_Memory)v->amem,"Iteration",&v->iter,1,AMS_INT,AMS_READ,AMS_COMMON,AMS_REDUCT_UNDEF);CHKERRQ(info);
363:   info = AMS_Memory_add_field((AMS_Memory)v->amem,"Residual",&v->fc,1,AMS_DOUBLE,AMS_READ,AMS_COMMON,AMS_REDUCT_UNDEF);CHKERRQ(info);
364:   info = AMS_Memory_add_field((AMS_Memory)v->amem,"Residual",&v->norm,1,AMS_DOUBLE,AMS_READ,AMS_COMMON,AMS_REDUCT_UNDEF);CHKERRQ(info);
365:   info = AMS_Memory_add_field((AMS_Memory)v->amem,"Iteration",(int*)&v->reason,1,AMS_INT,AMS_READ,AMS_COMMON,AMS_REDUCT_UNDEF);CHKERRQ(info);
366:   info = PetscObjectPublishBaseEnd(obj);CHKERRQ(info);
367: #endif
368:   TaoFunctionReturn(0);
369: }


374: int TaoObjectCreate( TAO_SOLVER *newsolver, MPI_Comm comm){
375:   TAO_SOLVER solver;
376:   int info;

379:   info = PetscHeaderCreate(solver,_p_TAO_SOLVER,int,TAO_COOKIE,-1,"TAO Solver",comm,TaoDestroy,0); CHKERRQ(info);
380:   info = PetscLogObjectCreate(solver); CHKERRQ(info);
381:   info = PetscLogObjectMemory(solver, sizeof(struct _p_TAO_SOLVER)); CHKERRQ(info);
382:   solver->bops->publish      = TaoPublish_Petsc;
383:   info=PetscPublishAll(solver);CHKERRQ(info);
384:   *newsolver = solver;
385:   return(0);
386: }

390: int TaoObjectDestroy( TAO_SOLVER solver){
391:   int info;
393:   /* if memory was published with AMS then destroy it */

395:   info = PetscObjectDepublish(solver);CHKERRQ(info);
396:   PetscLogObjectDestroy(solver);
397:   PetscHeaderDestroy(solver);

399:   return(0);
400: }
401: static int one = 1;
402: static char *executable = (char *)"tao";
403: static char **executablePtr = &executable;

407: int TaoLogClassRegister(int*cookie,const char*objectname){
408:   int info;
409:   int argc;
410:   char **args;

413:   info = TaoGetArgs(&argc,&args); CHKERRQ(info);

415: #if !defined(PARCH_t3d)
416:   info = PetscSetHelpVersionFunctions(TaoPrintHelpIntro,TaoPrintVersion); CHKERRQ(info);
417: #endif

419:   if (!PetscInitializeCalled) {
420:     if (argc&&args){
421:       info = PetscInitialize(&argc,&args,0,0); CHKERRQ(info);
422:     } else {
423:       info = PetscInitialize(&one,&executablePtr,0,0); CHKERRQ(info);
424:     }
425:     TaoBeganPetsc = TAO_TRUE;
426:   }
427:   info=PetscLogClassRegister(cookie,objectname);CHKERRQ(info);
428:   return(0);
429: }

431: #endif

433: /* The PetscObject is reduced and microkernal capabilities are absent
434: #define PetscObjectComposeFunctionDynamic(a,b,c,d) 0
435: #define PetscObjectQueryFunction(a,b,c) 0
436:   PetscLogInfo((void *,const char*,...){ return 0);}

438: */