Actual source code: daapp.c
1: #include "taodaapplication.h" /*I "taodaapplication.h" I*/
2: #include taoapp.h
3: #include "petsc.h"
4: #include "daapp_impl.h"
6: int DAAPP_COOKIE=0;
7: static char PPetscDAAppXXX[] = "PPetscDAApp";
9: int TaoAppDAApp(TAO_APPLICATION, DA_APPLICATION *);
10: int DAAppDestroy(DA_APPLICATION);
11: int DAAppExtend(TAO_APPLICATION);
15: /* @C
16: DAApplicationCreate - Creates a TaoApplication defined
17: a PETSc Distributed Array (DA) object. The routines used to evaluate
18: the objective function, constraint, and derivative information are
19: designed specifically for these problems.
21: Input Parameters:
22: + comm - an MPI communicator
23: . tda - an array of Distributed Array objects
24: - nnda - the number of Distibuted Array objects
26: Output Parameters:
27: . newdaapplication - the TaoApplication structure
29: .seealso TaoDAAppSolve(), TaoApplicationCreate(), TaoAppDestroy();
31: Level: beginner
33: .keywords: Application, DA
34: @ */
35: int DAApplicationCreate(MPI_Comm comm, DA* tda, int nnda, TAO_APPLICATION* newdaapplication){
36: int info;
37: TAO_APPLICATION ttapp;
39: info = TaoApplicationCreate(comm,&ttapp); CHKERRQ(info);
40: info = TaoAppSetDAApp(ttapp,tda,nnda); CHKERRQ(info);
41: info = PetscLogInfo((tda[0],"Creating a DA_APPLICATION object.\n")); CHKERRQ(info);
42: *newdaapplication=ttapp;
43: return(0);
44: }
48: int DAAppExtend(TAO_APPLICATION daapplication){
49: int info;
50: MPI_Comm comm;
51: DA_APPLICATION daapp;
54: if (DAAPP_COOKIE==0){
55: info=PetscLogClassRegister(&DAAPP_COOKIE,"TAO DA Application");CHKERRQ(info);
56: }
57: info=PetscObjectGetComm((PetscObject)daapplication,&comm); CHKERRQ(info);
58: info = PetscHeaderCreate(daapp,_p_DA_APPLICATION,int,DAAPP_COOKIE,-1,"DA Application",comm,DAAppDestroy,0); CHKERRQ(info);
59: info = PetscLogObjectCreate(daapp); CHKERRQ(info);
60: info = PetscLogObjectMemory(daapp, sizeof(struct _p_DA_APPLICATION)); CHKERRQ(info);
61:
62: daapp->nda=0;
63: daapp->ndamax=PETSCDAAPPMAXGRIDS;
64: daapp->currentlevel=0;
65: daapp->IsComplementarity=PETSC_FALSE;
66: daapp->kspflag=SAME_NONZERO_PATTERN;
67: daapp->computedafunction=0;
68: daapp->computedagradient=0;
69: daapp->computedafunctiongradient=0;
70: daapp->computedahessian=0;
71: daapp->usrdafctx=0; daapp->usrdafgctx=0; daapp->usrdagctx=0; daapp->usrdahctx=0;
72: daapp->computedabounds=0;
73: daapp->bounddactx=0;
74: daapp->nbeforemonitors=0;
75: daapp->naftermonitors=0;
76: int ii;
77: info=TaoAppAddObject(daapplication,PPetscDAAppXXX,(void*)daapp,&ii); CHKERRQ(info);
78: info=TaoAppSetDestroyRoutine(daapplication,(int(*)(void*))(TaoAppDestroyDAApp), (void*)daapplication); CHKERRQ(info);
79: info=TaoAppSetOptionsRoutine(daapplication,DAAppSetOptions); CHKERRQ(info);
81: return(0);
82: }
86: /*@
87: TaoAppSetDAApp - Extend the functionality of a TaoApplication by
88: adding support for applications defined
89: a PETSc Distributed Array (DA) object. The routines used to evaluate
90: the objective function, constraint, and derivative information are
91: designed specifically for these problems.
93: Input Parameters:
94: + daappliation - an existing application object
95: . tda - an array of Distributed Array objects
96: - nnda - the number of Distibuted Array objects
98: .seealso TaoDAAppSolve(), TaoApplicationCreate();
100: Level: beginner
102: .keywords: Application, DA
103: @*/
104: int TaoAppSetDAApp(TAO_APPLICATION daapplication, DA* tda, int nnda){
105: int i,info;
106: DA_APPLICATION daapp;
107: PetscScalar zero=0.0;
109: if (nnda > PETSCDAAPPMAXGRIDS){
110: SETERRQ1(1,"Number of grids cannot exceed %d.\n",PETSCDAAPPMAXGRIDS); }
111: info = DAAppExtend(daapplication); CHKERRQ(info);
112: info = TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
113: info = PetscLogInfo((daapplication,"TAO Extending TAO_APPLICATION for DA_APPLICATION object.\n")); CHKERRQ(info);
114: info = DAAppSetMatType(daapplication,MATAIJ); CHKERRQ(info);
115: daapp->nda=nnda;
116: daapp->currentlevel=0;
117: for (i=0;i<nnda; i++){
119: info = PetscObjectReference((PetscObject)(tda[i])); CHKERRQ(info);
120: daapp->grid[i].da=tda[i];
121: info = DACreateGlobalVector(daapp->grid[i].da,&daapp->grid[i].X); CHKERRQ(info);
122: info = VecDuplicate(daapp->grid[i].X,&daapp->grid[i].R); CHKERRQ(info);
123: info = VecDuplicate(daapp->grid[i].X,&daapp->grid[i].RHS); CHKERRQ(info);
124: daapp->grid[i].XL=0;
125: daapp->grid[i].XU=0;
126: daapp->grid[i].H=0;
127: daapp->grid[i].Interpolate=0;
128: daapp->grid[i].CScale=0;
129: if (i>0){
130: info = DAGetInterpolation(daapp->grid[i-1].da,daapp->grid[i].da,&daapp->grid[i].Interpolate,&daapp->grid[i].CScale); CHKERRQ(info);
131: }
132: }
133: info = VecSet(daapp->grid[0].X, zero); CHKERRQ(info);
134: info = TaoAppSetInitialSolutionVec(daapplication,daapp->grid[0].X); CHKERRQ(info);
135: info = PetscLogInfo((daapp,"Create work vectors for DA_APPLICATION object.\n")); CHKERRQ(info);
137: return(0);
138: }
142: int TaoAppDAApp(TAO_APPLICATION daapplication, DA_APPLICATION *daapp){
143: int info;
144: DA_APPLICATION daapp2;
147: info=TaoAppQueryForObject(daapplication,PPetscDAAppXXX,(void**)&daapp2); CHKERRQ(info);
148: if (daapp2){
150: } else {
151: SETERRQ(1,"TAO ERROR: Must call TaoAppSetDAApp() first");
152: }
153: *daapp=daapp2;
154: return(0);
155: }
159: int TaoAppDestroyDAApp(TAO_APPLICATION daapplication){
160: int info;
161: DA_APPLICATION daapp;
163: info=TaoAppQueryForObject(daapplication,PPetscDAAppXXX,(void**)&daapp); CHKERRQ(info);
164: if (daapp){
165: info=DAAppDestroy(daapp); CHKERRQ(info);
166: // info=TaoAppQueryRemoveObject(daapplication,PPetscDAAppXXX); CHKERRQ(info);
167: }
169: return(0);
170: }
174: int DAAppDestroy(DA_APPLICATION daapp){
175: int i,info,nnda=daapp->nda;
180: if (--daapp->refct > 0) return(0);
182: info = PetscLogInfo((daapp,"Destroying work vectors for DA_APPLICATION object.\n")); CHKERRQ(info);
183: for (i=0;i<nnda; i++){
184: if (daapp->grid[i].coloring){
185: info = ISColoringDestroy(daapp->grid[i].coloring); CHKERRQ(info) daapp->grid[i].coloring=0;
186: }
187: if (daapp->grid[i].H){
188: info = MatDestroy(daapp->grid[i].H); CHKERRQ(info);daapp->grid[i].H=0;
189: }
190: if (daapp->grid[i].XL && daapp->grid[i].XU){
191: info = VecDestroy(daapp->grid[i].XL); CHKERRQ(info);daapp->grid[i].XL=0;
192: info = VecDestroy(daapp->grid[i].XU); CHKERRQ(info);daapp->grid[i].XU=0;
193: }
194: if (i>0){
195: info = MatDestroy(daapp->grid[i].Interpolate); CHKERRQ(info);daapp->grid[i].Interpolate=0;
196: info = VecDestroy(daapp->grid[i].CScale); CHKERRQ(info);daapp->grid[i].CScale=0;
197: }
198: info = VecDestroy(daapp->grid[i].RHS); CHKERRQ(info);daapp->grid[i].RHS=0;
199: info = VecDestroy(daapp->grid[i].R); CHKERRQ(info);daapp->grid[i].R=0;
200: info = VecDestroy(daapp->grid[i].X); CHKERRQ(info);daapp->grid[i].X=0;
201: info = DADestroy(daapp->grid[i].da); CHKERRQ(info);daapp->grid[i].da=0;
202: }
203: daapp->nda=0;
204: daapp->currentlevel=0;
205: PetscLogObjectDestroy(daapp);
206: PetscHeaderDestroy(daapp);
207: return(0);
208: }
212: /*@
213: DAAppGetDA - Get the DA on a the specified level.
214:
216: Input Parameters:
217: + daapplication - the TAO DAApplication structure
218: - n - the level
220: Output Parameter:
221: . da - address of the pointer to the DA.
224: .seealso TaoAppSetDAApp();
226: Level: intermediate
228: .keywords: Application, DA
229: @*/
230: int DAAppGetDA(TAO_APPLICATION daapplication, int n, DA *da){
231: int info;
232: DA_APPLICATION daapp;
234: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
235: *da=daapp->grid[n].da;
236: return(0);
237: }
241: /*@
242: DAAppGetNumberOfDAGrids - Get the number of grids specified on the application
243:
245: Input Parameters:
246: . daapplication - the TAO DAApplication structure
248: Output Parameter:
249: . n - number of DA grids specified
251: .seealso TaoAppSetDAApp();
253: Level: intermediate
255: .keywords: Application, DA
256: @*/
257: int DAAppGetNumberOfDAGrids(TAO_APPLICATION daapplication, int *n){
258: int info;
259: DA_APPLICATION daapp;
261: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
262: *n=daapp->nda;
263: return(0);
264: }
268: /*@
269: DAAppGetCurrentLevel - Get the number of the current DA grid
271: Input Parameters:
272: . daapplication - the TAO DAApplication structure
274: Output Parameter:
275: . n - number of DA grids specified
277: .seealso TaoAppSetDAApp();
279: Level: intermediate
281: .keywords: Application, DA
282: @*/
283: int DAAppGetCurrentLevel(TAO_APPLICATION daapplication, int *n){
284: int info;
285: DA_APPLICATION daapp;
287: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
288: *n=daapp->currentlevel;
289: return(0);
290: }
295: static int DAAppPetscHessian(TAO_APPLICATION daapplication , Vec X , Mat *M, Mat *MPre, MatStructure *flag, void*ctx){
296: int info, clevel;
297: Mat H=*M;
298: DA_APPLICATION daapp=(DA_APPLICATION)ctx;
303: PetscStackPush("TAO User DA Hessian");
304: clevel=daapp->currentlevel;
305:
306: // daapp->FDGrad=daapp->grid[i].RHS; /* Needed for finite differences gradient vector */
307: info = PetscLogInfo((daapp,"TAO hessian evaluation at DA_APPLICATION object, level %d.\n",clevel)); CHKERRQ(info);
308: info = (*daapp->computedahessian)(daapplication,daapp->grid[clevel].da,X,H,daapp->usrdahctx);CHKERRQ(info);
309: *flag=daapp->kspflag;
311: PetscStackPop;
313: return(0);
314: }
319: /*@
320: DAAppSetHessianMat - Sets the matrices used to store the Hessian
322: Collective on TAO_APPLICATION
324: Input Parameters:
325: . daapplication - the DA Application object
327: Level: intermediate
329: Note: This routine should be used when finite differences are used to compute
330: the Hessian.
332: .keywords: DA, hessian
334: .seealso: DAAppSetHessianRoutine(), DAAppSetMatType();
336: @*/
337: int DAAppSetHessianMat(TAO_APPLICATION daapplication){
338: int i,info,nnda;
339: DA_APPLICATION daapp;
342: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
343: nnda=daapp->nda;
344: for (i=0;i<nnda; i++){
345: if (daapp->grid[i].H==0){
346: info = DAGetMatrix(daapp->grid[i].da,daapp->HessianMatrixType,&daapp->grid[i].H); CHKERRQ(info);
347: info = MatSetOption(daapp->grid[i].H,MAT_NO_NEW_NONZERO_LOCATIONS); CHKERRQ(info);
348: info = MatSetOption(daapp->grid[i].H,MAT_STRUCTURALLY_SYMMETRIC); CHKERRQ(info);
349: info = MatSetOption(daapp->grid[i].H,MAT_SYMMETRIC); CHKERRQ(info);
350: info = PetscLogInfo((daapp->grid[i].da,"TAO create Hessian for DA_APPLICATION object, level %d.\n",i)); CHKERRQ(info);
351: info = DAGetColoring(daapp->grid[i].da,IS_COLORING_LOCAL,&(daapp->grid[i].coloring));CHKERRQ(info);
352: }
353: }
354: info = TaoAppSetHessianMat(daapplication,daapp->grid[0].H,daapp->grid[0].H);CHKERRQ(info);
355:
356: return(0);
357: }
360: /*@C
361: DAAppSetHessianRoutine - Set a routine that will evaluate the hessian
362: on the given DA at the given point.
364: Collective on TAO_APPLICATION
366: Input Parameters:
367: + daapplication - the DA Application object
368: . hess - the function pointer for the hessian evaluation routine
369: - ctx - the hessian context
371: Calling sequence of hess:
372: $ hess(TAO_APPLICATION daapplication,DA da, Vec x,Mat H,void *ctx);
374: + daapplication - the TAO_APPLICATION context
375: . da - the Distributed Array
376: . x - input vector
377: . H - hessian matrix
378: - ctx - user-defined hessian context set from DAAppSetHessianRoutine()
380: Level: beginner
382: Options Database Key:
383: . -tao_view_hessian - view the hessian after each evaluation using PETSC_VIEWER_STDOUT_WORLD
385: .keywords: DA, hessian
387: .seealso: DAAppSetObjectiveAndGradientRoutine();
388: @*/
389: int DAAppSetHessianRoutine(TAO_APPLICATION daapplication, int (*hess)(TAO_APPLICATION,DA,Vec,Mat,void*),void *ctx){
390: int info;
391: DA_APPLICATION daapp;
394: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
395: daapp->computedahessian=hess;
396: daapp->usrdahctx=ctx;
397: if (hess){
398: info = TaoAppSetHessianRoutine(daapplication,DAAppPetscHessian,(void*)daapp);CHKERRQ(info);
399: } else {
400: info = TaoAppSetHessianRoutine(daapplication,0,0);CHKERRQ(info);
401: }
402: info = DAAppSetHessianMat(daapplication);CHKERRQ(info);
403: return(0);
404: }
408: static int DAAppPetscFunctionGradient(TAO_APPLICATION daapplication , Vec X ,double *ff, Vec G, void*ctx){
409: int info;
410: DA_APPLICATION daapp=(DA_APPLICATION)ctx;
415: PetscStackPush("TAO User DA Objective and gradient function");
416: info = (*daapp->computedafunctiongradient)(daapplication,daapp->grid[daapp->currentlevel].da,X,ff,G,daapp->usrdafgctx);
417: CHKERRQ(info);
418: info = PetscLogInfo((daapplication,"TAO function evaluation at DA_APPLICATION object, level %d.\n",daapp->currentlevel)); CHKERRQ(info);
419: PetscStackPop;
421: return(0);
422: }
427: /*@C
428: DAAppSetObjectiveAndGradientRoutine - Set a routine that will evaluate the objective and
429: gradient functions on the given DA at the given point.
431: Collective on TAO_APPLICATION
433: Input Parameters:
434: + daapplication - the DA Application object
435: . grad - the function pointer for the gradient evaluation routine
436: - ctx - the function-gradient context
438: Calling sequence of funcgrad:
439: $ funcgrad(TAO_APPLICATION daapplication,DA da, Vec x,double *f,Vec g,void *ctx);
441: + daapplication - the TAO_APPLICATION daapplication context
442: . da - the Distributed Array
443: . x - input vector
444: . f - objective value
445: . g - gradient vector
446: - ctx - user defined function-gradient context set from DAAppSetObjectiveAndGradientRoutine()
448: Fortran Note:
449: If your Fortran compiler does not recognize symbols over 31 characters in length, then
450: use the identical routine with the shortened name DAAppSetObjectiveAndGradientRou()
452: Level: beginner
454: Options Database Key:
455: . -tao_view_gradient - view the gradient after each evaluation using PETSC_VIEWER_STDOUT_SELF
457: .keywords: DA, Gradient, Objective Function
459: .seealso: DAAppSetObjectiveRoutine(), DAAppSetGradientRoutine();
461: @*/
462: int DAAppSetObjectiveAndGradientRoutine(TAO_APPLICATION daapplication, int (*funcgrad)(TAO_APPLICATION,DA,Vec,double*,Vec, void*),void *ctx){
463: int info;
464: DA_APPLICATION daapp;
466: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
467: if (funcgrad){
468: info=TaoAppSetObjectiveAndGradientRoutine(daapplication,
469: DAAppPetscFunctionGradient,(void*)daapp);CHKERRQ(info);
470: } else {
471: info=TaoAppSetObjectiveAndGradientRoutine(daapplication,0,0);CHKERRQ(info);
472: }
473: info=TaoAppSetInitialSolutionVec(daapplication,daapp->grid[0].X);CHKERRQ(info);
474: daapp->computedafunctiongradient=funcgrad;
475: daapp->usrdafgctx=ctx;
476: info = PetscLogInfo((daapp,"Set objective function pointer for DA_APPLICATION object.\n")); CHKERRQ(info);
477: return(0);
478: }
482: static int DAAppPetscGradient(TAO_APPLICATION daapplication , Vec X ,Vec G, void*ctx){
483: int info;
484: DA_APPLICATION daapp=(DA_APPLICATION)ctx;
489: PetscStackPush("TAO User DA Gradient Evaluation");
490: info = (*daapp->computedagradient)(daapplication,daapp->grid[daapp->currentlevel].da,X,G,daapp->usrdafgctx);
491: CHKERRQ(info);
492: info = PetscLogInfo((daapplication,"TAO gradient evaluation at DA_APPLICATION object, level %d.\n",daapp->currentlevel)); CHKERRQ(info);
493: PetscStackPop;
495: return(0);
496: }
501: /*@C
502: DAAppSetGradientRoutine - Set a routine that will evaluate the gradient
503: function on the given DA at the given point.
505: Collective on TAO_APPLICATION
507: Input Parameters:
508: + daapplication - the DA Application object
509: . grad - the function pointer for the gradient evaluation routine
510: - ctx - the gradient context
512: Calling sequence of grad:
513: $ grad(TAO_APPLICATION daapplication,DA da, Vec x,Vec g,void *ctx);
515: + daapplication - the TAO_APPLICATION context
516: . da - the Distributed Array
517: . x - input vector
518: . g - gradient vector
519: - ctx - user defined gradient context set from DAAppSetGradientRoutine()
521: Level: intermediate
523: Options Database Key:
524: . -tao_view_gradient - view the gradient after each evaluation using PETSC_VIEWER_STDOUT_SELF
526: .keywords: DA, gradient
528: .seealso: DAAppSetObjectiveRoutine(), DAAppSetObjectiveAndGradientRoutine();
530: @*/
531: int DAAppSetGradientRoutine(TAO_APPLICATION daapplication, int (*grad)(TAO_APPLICATION,DA,Vec,Vec, void*),void *ctx){
532: int info;
533: DA_APPLICATION daapp;
535: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
536: if (grad){
537: info=TaoAppSetGradientRoutine(daapplication, DAAppPetscGradient, (void*)daapp);CHKERRQ(info);
538: } else {
539: info=TaoAppSetGradientRoutine(daapplication,0,0);CHKERRQ(info);
540: }
541: daapp->computedagradient=grad;
542: daapp->usrdagctx=ctx;
543: info = PetscLogInfo((daapp,"Set gradient function pointer for DA_APPLICATION object.\n")); CHKERRQ(info);
544: return(0);
545: }
549: static int DAAppPetscObjectiveFunction(TAO_APPLICATION daapplication , Vec X ,double* f, void*ctx){
550: int info;
551: DA_APPLICATION daapp=(DA_APPLICATION)ctx;
556: PetscStackPush("TAO User DA Objective function");
557: info = (*daapp->computedafunction)(daapplication,daapp->grid[daapp->currentlevel].da,X,f,daapp->usrdafgctx);
558: CHKERRQ(info);
559: info = PetscLogInfo((daapp,"TAO gradient evaluation at DA_APPLICATION object, level %d.\n",daapp->currentlevel)); CHKERRQ(info);
560: PetscStackPop;
562: return(0);
563: }
568: /*@C
569: DAAppSetObjectiveRoutine - Set a routine that will evaluate the objective
570: function on the given DA at the given point.
572: Collective on TAO_APPLICATION
574: Input Parameters:
575: + daapplication - the DA Application object
576: . func - the function pointer for the objecive evaluation routine
577: - ctx - the monitor context
579: Calling sequence of func:
580: $ func(TAO_APPLICATION daapplication,DA da,Vec x,double *f,void *ctx);
582: + daapplication - the TAO_APPLICATION context
583: . da - the Distributed Array
584: . x - input vector
585: . f - application sets equal to the function value
586: - ctx - user-defined function context set from DAAppSetObjectiveRoutine()
588: Level: beginner
590: .keywords: DA, objective
592: .seealso: DAAppSetObjectiveAndGradientRoutine();
594: @*/
595: int DAAppSetObjectiveRoutine(TAO_APPLICATION daapplication, int (*func)(TAO_APPLICATION,DA,Vec,double*, void*),void *ctx){
596: int info;
597: DA_APPLICATION daapp;
599: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
600: if (func){
601: info=TaoAppSetObjectiveRoutine(daapplication, DAAppPetscObjectiveFunction, (void*)daapp);CHKERRQ(info);
602: } else {
603: info=TaoAppSetObjectiveRoutine(daapplication,0,0);CHKERRQ(info);
604: }
605: info=TaoAppSetInitialSolutionVec(daapplication, daapp->grid[0].X);CHKERRQ(info);
606: daapp->computedafunction=func;
607: daapp->usrdafctx=ctx;
608: info = PetscLogInfo((daapp,"Set objective function pointer for DA_APPLICATION object.\n")); CHKERRQ(info);
609: return(0);
610: }
614: int DAApplicationEvaluateVariableBounds(TAO_APPLICATION daapplication, Vec XL, Vec XU, void*ctx){
615: DA_APPLICATION daapp=(DA_APPLICATION)ctx;
616: int level=daapp->currentlevel, info;
620: info = PetscLogInfo((daapp->grid[level].da,"Compute variable bounds in level %d of DA_APPLICATION object.\n",level)); CHKERRQ(info);
621: info = PetscObjectReference((PetscObject)XL);
622: info = PetscObjectReference((PetscObject)XU);
623: daapp->grid[level].XL=XL;
624: daapp->grid[level].XU=XU;
625: if (daapp->computedabounds){
626: info = (*daapp->computedabounds)(daapplication,daapp->grid[level].da,XL,XU,daapp->bounddactx); CHKERRQ(info);
627: }
628: return(0);
629: }
634: /*@C
635: DAAppSetVariableBoundsRoutine - Set a routine that will evaluate the bounds
636: on the variables.
638: Collective on TAO_APPLICATION
640: Input Parameters:
641: + daapplication - the DA Application object
642: . bound - the function pointer for the bound evaluation routine
643: - ctx - the monitor context
645: Calling sequence of bound:
646: $ bound(TAO_APPLICATION daapplication, DA da,Vec xl, Vec xu, void *ctx);
648: + daapplication - the TAO_APPLICATION context
649: . da - the Distributed Array
650: . xl - vector of upper bounds
651: . xu - vector of lower bounds
652: - ctx - user-defined monitor context set from DAAppSetVariableBoundsRoutine()
654: Level: beginner
656: .keywords: DA, bounds
658: .seealso: TaoDAAppSolve();
660: @*/
661: int DAAppSetVariableBoundsRoutine(TAO_APPLICATION daapplication, int (*bound)(TAO_APPLICATION,DA,Vec,Vec, void*),void *ctx){
662: int info;
663: DA_APPLICATION daapp;
665: info = TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
666: info=TaoAppSetVariableBoundsRoutine(daapplication,DAApplicationEvaluateVariableBounds,(void*)daapp);CHKERRQ(info);
667: daapp->computedabounds=bound;
668: daapp->bounddactx=ctx;
669: info = PetscLogInfo((daapp,"Set variable bounds in DA_APPLICATION object.\n")); CHKERRQ(info);
670: return(0);
671: }
677: /*@C
678: DAAppSetConstraintRoutine - Set a routine that will evaluate the constraint
679: function on the given DA at the given point.
681: Collective on TAO_APPLICATION
683: Input Parameters:
684: + daapplication - the TAO Application object
685: . grad - the function pointer for the gradient evaluation routine
686: - ctx - the gradient context
688: Calling sequence of grad:
689: $ f(TAO_APPLICATION daapplication,DA da, Vec x,Vec r,void *ctx);
691: + daapplication - the DA_APPLICATION context
692: . da - the Distributed Array
693: . x - input vector
694: . r - constraint vector
695: - ctx - user defined gradient context set from DAAppSetGradientRoutine()
697: Level: intermediate
699: Options Database Key:
700: . -tao_view - view the gradient after each evaluation using PETSC_VIEWER_STDOUT_SELF
702: .keywords: DA, gradient
704: .seealso: DAAppSetJacobianRoutine();
705: @*/
706: int DAAppSetConstraintRoutine(TAO_APPLICATION daapplication, int (*f)(TAO_APPLICATION,DA,Vec,Vec, void*),void *ctx){
707: int info;
708: DA_APPLICATION daapp;
710: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
711: if (f){
712: info=TaoAppSetConstraintRoutine(daapplication, DAAppPetscGradient, (void*)daapp);CHKERRQ(info);
713: } else {
714: info=TaoAppSetConstraintRoutine(daapplication,0,0);CHKERRQ(info);
715: }
716: daapp->IsComplementarity=PETSC_TRUE;
717: info=DAAppSetGradientRoutine(daapplication,f,ctx); CHKERRQ(info);
718: daapp->computedagradient=f;
719: daapp->usrdagctx=ctx;
720: daapp->usrdafgctx=ctx;
722: info = PetscLogInfo((daapp,"Set constraint function pointer for DA_APPLICATION object.\n")); CHKERRQ(info);
723: return(0);
724: }
729: /*@C
730: DAAppSetJacobianRoutine - Set a routine that will evaluate the Jacobian
731: on the given DA at the given point.
733: Collective on TAO_APPLICATION
735: Input Parameters:
736: + daapplication - the DA Application object
737: . jac - the function pointer for the Jacobian evaluation routine
738: - ctx - the jacobian context
740: Calling sequence of hess:
741: $ jac(TAO_APPLICATION daapplication, DA da, Vec x,Mat J,void *ctx);
743: + daapplication - the TAO_APPLICATION context
744: . da - the Distributed Array
745: . x - input vector
746: . J - Jacobian matrix
747: - ctx - user-defined hessian context set from DAAppSetHessianRoutine()
749: Level: intermediate
751: Options Database Key:
752: . -tao_view_jacobian - view the jacobian after each evaluation using PETSC_VIEWER_STDOUT_WORLD
754: .keywords: DA, hessian
756: .seealso: DAAppSetConstraintRoutine();
758: @*/
759: int DAAppSetJacobianRoutine(TAO_APPLICATION daapplication, int (*jac)(TAO_APPLICATION,DA,Vec,Mat,void*),void *ctx){
760: int i,info,nnda;
761: DA_APPLICATION daapp;
763: info = TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
764: if (!jac){ SETERRQ(1,"No DA Hessian routine provided\n");}
765: daapp->computedahessian=jac;
766: daapp->usrdahctx=ctx;
767: nnda=daapp->nda;;
768: for (i=0;i<nnda; i++){
769: if (daapp->grid[i].H==0){
770: info = DAGetMatrix(daapp->grid[i].da,daapp->HessianMatrixType,&daapp->grid[i].H); CHKERRQ(info);
771: info = MatSetOption(daapp->grid[i].H,MAT_NO_NEW_NONZERO_LOCATIONS); CHKERRQ(info);
772: info = MatSetOption(daapp->grid[i].H,MAT_STRUCTURALLY_SYMMETRIC); CHKERRQ(info);
773: info = PetscLogInfo((daapp->grid[i].da,"TAO create Jacobian for DA_APPLICATION object, level %d.\n",i)); CHKERRQ(info);
774: }
775: }
776:
777: info = TaoAppSetJacobianRoutine(daapplication,DAAppPetscHessian,(void*)daapp);CHKERRQ(info);
778: info = TaoAppSetJacobianMat(daapplication,daapp->grid[0].H,daapp->grid[0].H);CHKERRQ(info);
779:
780: return(0);
781: }
785: /*@C
786: DAAppSetBeforeMonitor - Set a routine that will be called before the
787: optimization on each grid.
789: Collective on TAO_APPLICATION
791: Input Parameters:
792: + daapplication - the DA Application object
793: . beforemonitor - a monitor routine called before solving the application on each DA
794: - ctx - the monitor context
796: Calling sequence of monitor:
797: $ beforemonitor(TAO_APPLICATION daapplication, DA da,int level, void *ctx);
799: + daapplication - this TAO_APPLICATION context
800: . da - the Distributed Array
801: . level - the grid level that will be solved next (level 0 is coarsest)
802: - ctx - user-defined function context set from DAAppSetBeforeMonitor()
804: Note:
805: These monitors are different from the monitors that can be called after
806: each iteration of the optimization algorithm.
808: Note:
809: The beforemonitor and aftermonitor are good for setting up and destroying the application data.
811: Level: intermediate
813: .keywords: DA, monitor
815: .seealso: DAAppSetAfterMonitor(), TaoSetMonitor(), TaoAppSetMonitor();
817: @*/
818: int DAAppSetBeforeMonitor(TAO_APPLICATION daapplication, int (*beforemonitor)(TAO_APPLICATION,DA,int, void*), void *ctx){
819: int n,info;
820: DA_APPLICATION daapp;
822: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
823: if (beforemonitor){
824: n=daapp->nbeforemonitors;
825: if (n>=MAX_DAAP_MONITORS){
826: SETERRQ(1,"TAO ERROR: Too many beforemonitors set in DAAPP");
827: }
828: daapp->beforemonitor[n]=beforemonitor;
829: daapp->beforemonitorctx[n]=ctx;
830: daapp->nbeforemonitors++;
831: info = PetscLogInfo((daapplication,"TAO: Set user beforemonitor for DA_APPLICATION object.\n")); CHKERRQ(info);
832: }
833: return(0);
834: }
838: /*@C
839: DAAppSetAfterMonitor - Set a routine that will be called after the
840: optimization on each grid.
842: Collective on TAO_APPLICATION
844: Input Parameters:
845: + daapplication - the DA Application object
846: . aftermonitor - a monitor routine called after solving the application on each DA
847: - ctx - the monitor context
849: Calling sequence of monitor:
850: $ aftermonitor(TAO_APPLICATION daapplication, DA da, int level, void *ctx);
852: + daapplication - this TAO_APPLICATION context
853: . da - the Distributed Array
854: . level - the grid level that will be solved next (level 0 is coarsest)
855: - ctx - user-defined function context set from DAAppSetAfterMonitor()
857: Note:
858: These monitors are different from the monitors that can be called after
859: each iteration of the optimization algorithm.
861: Note:
862: The beforemonitor and aftermonitor are good for setting up and destroying the application data.
864: Level: intermediate
866: .keywords: DA, monitor
868: .seealso: DAAppSetBeforeMonitor(), TaoSetMonitor(), TaoAppSetMonitor();
870: @*/
871: int DAAppSetAfterMonitor(TAO_APPLICATION daapplication, int (*aftermonitor)(TAO_APPLICATION,DA,int, void*), void *ctx){
872: int info,n;
873: DA_APPLICATION daapp;
875: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
876: if (aftermonitor){
877: n=daapp->naftermonitors;
878: if (n>=MAX_DAAP_MONITORS){
879: SETERRQ(1,"TAO ERROR: Too many aftermonitors set in DAAPP");
880: }
881: daapp->aftermonitor[n]=aftermonitor;
882: daapp->aftermonitorctx[n]=ctx;
883: daapp->naftermonitors++;
884: info = PetscLogInfo((daapp,"TAO: Set user aftermonitor for DA_APPLICATION object.\n")); CHKERRQ(info);
885: }
886: return(0);
887: }
893: /*@
894: DAAppGetSolution - Sets a pointer to an existing vector that contains
895: the solution.
897: Collective on TAO_APPLICATION
899: Input Parameters:
900: + daapplication - the DA Application object
901: - level - the DA on which the current solution is desired (Level 0 is coarsest)
903: Output Parameters:
904: . X - the current solution
906: Level: beginner
908: .seealso: DAAppSetBeforeMonitor(), DAAppSetAfterMonitor();
910: .keywords: Solution, DA
911: @*/
912: int DAAppGetSolution(TAO_APPLICATION daapplication, int level, Vec *X){
913: int info;
914: DA_APPLICATION daapp;
916: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
917: if (level<0 || level>=daapp->nda){
918: SETERRQ(1,"No Solution available. This grid does not exist");
919: } else if (X && daapp && daapp->grid && level<daapp->nda){
920: *X=daapp->grid[level].X;
921: }
922: return(0);
923: }
928: /*@C
929: DAAppGetHessianMat - Sets a pointer to an existing matrix that contains
930: the Hessian at the specified level.
932: Collective on TAO_APPLICATION
934: Input Parameters:
935: + daapplication - the DA Application object
936: - level - the DA on which the current Hessian is desired (Level 0 is coarsest)
938: Output Parameters:
939: . H - the current Hessian
941: Level: intermediate
943: .seealso: DAAppSetHessianRoutine;
945: .keywords: Hessian, DA
946: @*/
947: int DAAppGetHessianMat(TAO_APPLICATION daapplication, int level, Mat *H){
948: int info;
949: DA_APPLICATION daapp;
951: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
952: if (level<0 || level>=daapp->nda){
953: SETERRQ(1,"No Hessian available. This grid does not exist");
954: } else if (H && daapp && daapp->grid && level<daapp->nda){
955: *H=daapp->grid[level].H;
956: }
957: return(0);
958: }
962: /*@
963: DAAppGetInterpolationMatrix - Sets pointers to existing vectors
964: containg upper and lower bounds on the variables.
966: Collective on TAO_APPLICATION
968: Input Parameters:
969: + daapplication - the DA Application object
970: - level - the DA on which the current solution is desired (Level 0 is coarsest)
972: Output Parameters:
973: + Interpolate - the interpolating matrix
974: - CScale - the scaling matrix for restriction Matrix
976: Level: intermediate
978: .seealso: DAAppSetBeforeMonitor(), DAAppSetAfterMonitor();
980: .keywords: Interpolate, DA
981: @*/
982: int DAAppGetInterpolationMatrix(TAO_APPLICATION daapplication, int level, Mat *Interpolate, Vec *CScale){
983: int info;
984: DA_APPLICATION daapp;
986: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
987: if ( !daapp || level<1 || level>=daapp->nda){
988: SETERRQ1(1,"No Interpolator available for level %d. This grid does not exist",level);
989: }
990: if (Interpolate){
991: *Interpolate=daapp->grid[level].Interpolate;
992: }
993: if (CScale){
994: *CScale=daapp->grid[level].CScale;
995: }
996: return(0);
997: }
1002: /*@
1003: DAAppGetVariableBounds - Sets pointers to existing vectors
1004: containg upper and lower bounds on the variables.
1006: Collective on TAO_APPLICATION
1008: Input Parameters:
1009: + daapplication - the DA Application object
1010: - level - the DA on which the current solution is desired (Level 0 is coarsest)
1012: Output Parameters:
1013: + XL - the current solution
1014: - XU - the current solution
1016: Level: intermediate
1018: .seealso: DAAppSetBeforeMonitor(), DAAppSetAfterMonitor();
1020: .keywords: Solution, DA
1021: @*/
1022: int DAAppGetVariableBounds(TAO_APPLICATION daapplication, int level, Vec *XL, Vec *XU){
1023: int info;
1024: DA_APPLICATION daapp;
1026: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
1027: if ( !daapp || level<0 || level>=daapp->nda){
1028: SETERRQ(1,"No Solution available. This grid does not exist.");
1029: }
1030: if (XL){
1031: *XL=daapp->grid[level].XL;
1032: }
1033: if (XU){
1034: *XU=daapp->grid[level].XU;
1035: }
1036: return(0);
1037: }
1043: /*@
1044: DAAppSetInitialSolution - Sets the initial solution for
1045: the grid application.
1047: Collective on TAO_APPLICATION
1049: Input Parameters:
1050: + daapplication - the DA Application object
1051: - X0 - the initial solution vector
1053: Level: beginner
1055: .seealso: DAAppGetSolution();
1057: .keywords: Initial Solution, DA
1058: @*/
1059: int DAAppSetInitialSolution(TAO_APPLICATION daapplication, Vec X0){
1060: int info;
1061: DA_APPLICATION daapp;
1064: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
1065: info=VecCopy(X0,daapp->grid[0].X); CHKERRQ(info);
1066: return(0);
1067: }
1073: /*
1074: DAAppSetMatStructure - Set the MatStructure flag to be used in KSPSetOperators:
1076: Collective on TAO_APPLICATION
1077:
1078: Input Parameters:
1079: + daapplication - the DA Application object
1080: - flag - indicates the KSP flag
1082: Level: intermediate
1084: Note:
1085: The default value is SAME_NONZERO_PATTERN
1087: .seealso: KSPSetOperators(), TaoGetKSP();
1089: .keywords: Linear Solver, Multigrid, DA, KSP
1091: @*/
1092: int DAAppSetMatStructure(TAO_APPLICATION daapplication, MatStructure pflag){
1093: int info;
1094: DA_APPLICATION daapp;
1096: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
1097: daapp->kspflag=pflag;
1098: info = PetscLogInfo((daapp,"Set MatStructure flag in TAO solver.\n")); CHKERRQ(info);
1099: return(0);
1100: }
1104: /*@
1105: DAAppSetMatType - set the matrix type to be used for the Hessian matrix
1107: Collective on TAO_APPLICATION
1108:
1109: Input Parameters:
1110: + daapplication - the DA Application object
1111: - mattype - the matrix type
1113: Level: intermediate
1115: Options Database Key:
1116: . -tao_da_mattype - select matrix type
1118: Note:
1119: The default value is MATMPIAIJ.
1121: Note:
1122: Call before DAAppSetHessianRoutine()
1124: .seealso: DAAppSetMatStructure(), DAAppSetHessianRoutine(), DAGetMatrix();
1126: .keywords: DA, Hessian
1128: @*/
1129: int DAAppSetMatType(TAO_APPLICATION daapplication, MatType mattype){
1130: int info;
1131: DA_APPLICATION daapp;
1133: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
1134: info=PetscStrcpy(daapp->HessianMatrixType,mattype); CHKERRQ(info);
1135: return(0);
1136: }
1141: /*@
1142: DAAppSetOptions - Sets various parameters to be used in this application
1143: and the TAO solver.
1145: Collective on TAO_APPLICATION
1147: Input Parameters:
1148: . daapplication - the DA Application object
1150: Level: intermediate
1152: Note: This routine will be called by the TaoAppSetFromOptions().
1154: .keywords: DA, options
1156: .seealso: TaoDAAppSolve();
1158: @*/
1159: int DAAppSetOptions(TAO_APPLICATION daapplication){
1160: int info;
1161: // char *prefix;
1162: PetscTruth flg1=PETSC_FALSE,flg2=PETSC_FALSE;
1163: DA_APPLICATION daapp;
1164: // MPI_Comm comm;
1168: info = PetscLogInfo((daapplication,"TaoDAAppSetOptions(): Reading command line for options\n")); CHKERRQ(info);
1169: info = TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
1170: /*
1171: info = PetscObjectGetOptionsPrefix((PetscObject)daapplication,&prefix); CHKERRQ(info);
1172: info = PetscObjectGetComm((PetscObject)daapplication,&comm); CHKERRQ(info);
1173: info = PetscOptionsBegin(comm,prefix,"TAO PETSC DA APPLICATIONS","solver");CHKERRQ(info);
1174: */
1175: info = DAAppSetMultiGridOptions(daapplication);CHKERRQ(info);
1177: flg1=PETSC_FALSE,flg2=PETSC_FALSE;
1178: info = PetscOptionsTruth("-tao_da_samepreconditioner","Use same preconditioner in KSP","DAAppSetMatStructure",
1179: PETSC_FALSE,&flg2,&flg1);CHKERRQ(info);
1180: if (flg1 && flg2) {
1181: info = DAAppSetMatStructure(daapplication, SAME_PRECONDITIONER);CHKERRQ(info);
1182: } else if (flg1){
1183: info = DAAppSetMatStructure(daapplication, SAME_NONZERO_PATTERN);CHKERRQ(info);
1184: }
1186: info = PetscOptionsString("-tao_da_mattype","the hessian should have matrix type","DAAppSetMatType",
1187: daapp->HessianMatrixType,daapp->HessianMatrixType,18,&flg1);CHKERRQ(info);
1189: /*
1190: info = PetscOptionsEnd();CHKERRQ(info);
1191: */
1192: return(0);
1193: }