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: /*@
 67:   TaoPrintStatement - prints a character string to stdout.

 69:   Not Collective
 70:   
 71:   Input Parameters:
 72: + tao - the TAO_SOLVER solver context
 73: - statement - the string to print

 75:   Level: beginner
 76: @*/
 77: int TaoPrintStatement(TAO_SOLVER tao, const char *statement){
 80:   PetscPrintf(tao->comm,statement);
 81:   return(0);
 82: }

 86: int TaoPrintInt(TAO_SOLVER tao, const char *statement, int n){
 89:   PetscPrintf(tao->comm,statement,n);
 90:   return(0);
 91: }

 95: int TaoPrintDouble(TAO_SOLVER tao, const char *statement,double dd){
 98:   PetscPrintf(tao->comm,statement,dd);
 99:   return(0);
100: }

104: int TaoPrintString(TAO_SOLVER tao, const char *statement,const char *str){
107:   PetscPrintf(tao->comm,statement,str);
108:   return(0);
109: }

113: int TaoOptionsHead(const char *heading){
114:   int info;
116:   info = PetscOptionsHead(heading);CHKERRQ(info);
117:   return(0);
118: }

122: int TaoOptionsTail(){
123:   int info;
125:   info = PetscOptionsTail();CHKERRQ(info);
126:   return(0);
127: }

131: int TaoOptionInt(const char *opt,const char *text,const char *man,int defaultv,int *value,TaoTruth *set){
132:   int info;
133:   PetscTruth flg=PETSC_FALSE;
135:   info = PetscOptionsInt(opt,text,man,defaultv,value,&flg);CHKERRQ(info);
136:   if (set){
137:     if (flg==PETSC_TRUE) *set=TAO_TRUE; else *set=TAO_FALSE;
138:   }
139:   return(0);
140: }

144: int TaoOptionDouble(const char *opt,const char *text,const char *man,double defaultv,double *value,TaoTruth *set){
145:   int info;
146:   PetscReal pdv=(PetscReal)defaultv, pv;
147:   PetscTruth flg=PETSC_FALSE;
149:   info = PetscOptionsReal(opt,text,man,pdv,&pv,&flg);CHKERRQ(info);
150:   if (set){
151:     if (flg==PETSC_TRUE) *set=TAO_TRUE; else *set=TAO_FALSE;
152:   }
153:   if (value&&flg==PETSC_TRUE){
154:     *value=(double)pv;
155:   }

157:   return(0);
158: }

162: int TaoOptionString(const char *opt,const char *text,const char *man,const char* defaultv,char *value, int len, TaoTruth *set){
163:   int info;
164:   PetscTruth flg=PETSC_FALSE;
166:   info = PetscOptionsString(opt,text,man,defaultv,value,len,&flg);CHKERRQ(info);
167:   if (set){
168:     if (flg==PETSC_TRUE) *set=TAO_TRUE; else *set=TAO_FALSE;
169:   }
170:   return(0);
171: }

175: int TaoOptionName(const char *opt,const char *text,const char *man,TaoTruth *set){
176:   int info;
177:   PetscTruth flg=PETSC_FALSE;
179:   info = PetscOptionsName(opt,text,man,&flg);CHKERRQ(info);
180:   if (set){
181:     if (flg==PETSC_TRUE) *set=TAO_TRUE; else *set=TAO_FALSE;
182:   }
183:   return(0);
184: }

188: int TaoOptionList(const char *opt, const char *ltext, const char *man, 
189:                   const char **list, int nlist, const char *defaultv, 
190:                   int *value, TaoTruth *set)
191: {
192:   int info;
193:   PetscTruth flg=PETSC_FALSE;

196:   info = PetscOptionsEList(opt, ltext, man, list, nlist, defaultv, value, &flg); CHKERRQ(info);
197:   if (set) {
198:     if (PETSC_TRUE == flg) {
199:       *set=TAO_TRUE; 
200:     }
201:     else {
202:       *set=TAO_FALSE;
203:     }
204:   }
205:   return(0);
206: }

210: int TaoMethodsList(const char *opt,const char *ltext,const char *man,const char *defaultv,char *value,int len,TaoTruth *set){
211:   int info;
212:   PetscTruth flg=PETSC_FALSE;
214:   info = PetscOptionsList(opt,ltext,man,TaoList,defaultv,value,len,&flg); CHKERRQ(info);
215:   if (set){
216:     if (flg==PETSC_TRUE) *set=TAO_TRUE; else *set=TAO_FALSE;
217:   }
218:   return(0);
219: }


224: int TaoFindSolver(TAO_SOLVER tao, TaoMethod type,  int (**r)(TAO_SOLVER) ){
225:   int info;
227:   info = PetscFListFind(TaoList,tao->comm,type,(void (**)(void))r);CHKERRQ(info);
228:   if (*r){
229:     info = PetscObjectChangeTypeName((PetscObject)tao,type);CHKERRQ(info);
230:   }
231:   return(0);
232: }

236: /*@C
237:    TaoRegisterDestroy - Frees the list of minimization solvers that were
238:    registered by TaoRegisterDynamic().

240:    Not Collective

242:    Level: advanced

244: .keywords: TAO_SOLVER, register, destroy

246: .seealso: TaoRegisterAll()
247: @*/
248: int TaoRegisterDestroy(){
249:   int info;
251:   if (TaoList) {
252:     info=PetscFListDestroy(&TaoList);CHKERRQ(info);
253:     TaoList=0;
254:   }
255:   if (TaoBeganPetsc) {
256:     info = PetscFinalize();CHKERRQ(info);
257:   }

259:   return(0);
260: }

264: /* --------------------------------------------------------------------- */
265: /*MC
266:    TaoRegisterDynamic - Adds a method to the TAO_SOLVER package for unconstrained minimization.

268:    Synopsis:
269:    TaoRegisterDynamic(char *name_solver,char *path,char *name_Create,int (*routine_Create)(TAO_SOLVER))

271:    Not collective

273:    Input Parameters:
274: +  name_solver - name of a new user-defined solver
275: .  path - path (either absolute or relative) the library containing this solver
276: .  name_Create - name of routine to Create method context
277: -  routine_Create - routine to Create method context

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

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

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

289:    Sample usage:
290: .vb
291:    TaoRegisterDynamic("my_solver","/home/username/my_lib/lib/libg_c++/solaris/mylib.a",
292:                 "MySolverCreate",MySolverCreate);
293: .ve

295:    Then, your solver can be chosen with the procedural interface via
296: $     TaoSetMethod(solver,"my_solver")
297:    or at runtime via the option
298: $     -tao_method my_solver

300:    Level: advanced

302: .keywords: TAO_SOLVER, register

304: .seealso: TaoRegisterAll(), TaoRegisterDestroy()
305: M*/
309: int TaoRegisterDynamic(const char *sname,const char *path,const char *name,int (*function)(TAO_SOLVER))
310: {
311:   char fullname[256];
312:   int  info;

315:   info = PetscFListConcat(path,name,fullname); CHKERRQ(info);
316:   //#if defined(PETSC_USE_DYNAMIC_LIBRARIES)
317:   //info = PetscFListAddDynamic(&TaoList,sname,fullname, 0);CHKERRQ(info);
318:   //#else
319:   info = PetscFListAddDynamic(&TaoList,sname,fullname, (void (*)())function);CHKERRQ(info);
320:   //#endif
321:   return(0);
322: }

327: /*@C
328:    TaoCompareMethod - Determines whether the TAO_SOLVER structure is of a
329:    specified type.

331:    Not Collective

333:    Input Parameter:
334: .  tao - the TAO_SOLVER solver context
335: .  type - a TAO_SOLVER solver method

337:    Output Parameter:
338: .  same - TAO_TRUE if 'tao' is of method 'type'

340:    Level: developer

342: .keywords: method
343: @*/
344: int TaoCompareMethod(TAO_SOLVER tao, TaoMethod type,TaoTruth *issame){
345:   int info;
346:   PetscTruth flag;
347:   TaoMethod currenttype;
350:   info = TaoGetMethod(tao,&currenttype);CHKERRQ(info);
351:   info = PetscStrcmp(type,currenttype,&flag);CHKERRQ(info);
352:   if (issame){
353:     if (flag==PETSC_TRUE) *issame=TAO_TRUE; else *issame=TAO_FALSE;
354:   }
355:   return(0);
356: }

360: int TaoStrcmp(const char *str1,const char *str2,TaoTruth *flag){
361:   int info;
362:   PetscTruth flg;
364:   info = PetscStrcmp(str1,str2,&flg);CHKERRQ(info);
365:   if (flg==PETSC_TRUE) *flag=TAO_TRUE; else *flag=TAO_FALSE;
366:   return(0);
367: }

371: int TaoStrcpy(char* str1,const char*str2){
372:   int info;
374:   info = PetscStrcpy(str1,str2);CHKERRQ(info);
375:   return(0);
376: }


381: static int TaoPublish_Petsc(PetscObject obj)
382: {
383: #if defined(PETSC_HAVE_AMS)
384:   TAO_SOLVER       v = (TAO_SOLVER) obj;
385:   int          info;
386: #endif

388:   TaoFunctionBegin;

390: #if defined(PETSC_HAVE_AMS)
391:   /* if it is already published then return */
392:   if (v->amem >=0 ) TaoFunctionReturn(0);

394:   info = PetscObjectPublishBaseBegin(obj);CHKERRQ(info);
395:   info = AMS_Memory_add_field((AMS_Memory)v->amem,"Iteration",&v->iter,1,AMS_INT,AMS_READ,AMS_COMMON,AMS_REDUCT_UNDEF);CHKERRQ(info);
396:   info = AMS_Memory_add_field((AMS_Memory)v->amem,"Residual",&v->fc,1,AMS_DOUBLE,AMS_READ,AMS_COMMON,AMS_REDUCT_UNDEF);CHKERRQ(info);
397:   info = AMS_Memory_add_field((AMS_Memory)v->amem,"Residual",&v->norm,1,AMS_DOUBLE,AMS_READ,AMS_COMMON,AMS_REDUCT_UNDEF);CHKERRQ(info);
398:   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);
399:   info = PetscObjectPublishBaseEnd(obj);CHKERRQ(info);
400: #endif
401:   TaoFunctionReturn(0);
402: }


407: int TaoObjectCreate( TAO_SOLVER *newsolver, MPI_Comm comm){
408:   TAO_SOLVER solver;
409:   int info;

412:   info = PetscHeaderCreate(solver,_p_TAO_SOLVER,int,TAO_COOKIE,-1,"TAO Solver",comm,TaoDestroy,0); CHKERRQ(info);
413:   info = PetscLogObjectCreate(solver); CHKERRQ(info);
414:   info = PetscLogObjectMemory(solver, sizeof(struct _p_TAO_SOLVER)); CHKERRQ(info);
415:   solver->bops->publish      = TaoPublish_Petsc;
416:   info=PetscPublishAll(solver);CHKERRQ(info);
417:   *newsolver = solver;
418:   return(0);
419: }

423: int TaoObjectDestroy( TAO_SOLVER solver){
424:   int info;
426:   /* if memory was published with AMS then destroy it */

428:   info = PetscObjectDepublish(solver);CHKERRQ(info);
429:   PetscLogObjectDestroy(solver);
430:   PetscHeaderDestroy(solver); 

432:   return(0);
433: }
434: static int one = 1;
435: static char *executable = (char *)"tao";
436: static char **executablePtr = &executable;

440: int TaoLogClassRegister(int*cookie,const char*objectname){
441:   int info;
442:   int argc;
443:   char **args;

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

448: #if !defined(PARCH_t3d)
449:   info = PetscSetHelpVersionFunctions(TaoPrintHelpIntro,TaoPrintVersion); CHKERRQ(info);
450: #endif

452:   if (!PetscInitializeCalled) {
453:     if (argc&&args){
454:       info = PetscInitialize(&argc,&args,0,0); CHKERRQ(info);
455:     } else {
456:       info = PetscInitialize(&one,&executablePtr,0,0); CHKERRQ(info);
457:     }
458:     TaoBeganPetsc = TAO_TRUE;
459:   }
460:   info=PetscLogClassRegister(cookie,objectname);CHKERRQ(info);
461:   return(0);
462: }

464: #endif

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

471: */