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, PetscInt 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 = PetscInfo(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=PetscCookieRegister("TAO DA Application",&DAAPP_COOKIE);CHKERRQ(info);
56: }
57: info=PetscObjectGetComm((PetscObject)daapplication,&comm); CHKERRQ(info);
58: info = PetscHeaderCreate(daapp,_p_DA_APPLICATION,PetscInt,DAAPP_COOKIE,-1,"DA Application",comm,DAAppDestroy,0); CHKERRQ(info);
59:
60: daapp->nda=0;
61: daapp->ndamax=PETSCDAAPPMAXGRIDS;
62: daapp->currentlevel=0;
63: daapp->IsComplementarity=PETSC_FALSE;
64: daapp->kspflag=SAME_NONZERO_PATTERN;
65: daapp->computedafunction=0;
66: daapp->computedagradient=0;
67: daapp->computedafunctiongradient=0;
68: daapp->computedahessian=0;
69: daapp->usrdafctx=0; daapp->usrdafgctx=0; daapp->usrdagctx=0; daapp->usrdahctx=0;
70: daapp->computedabounds=0;
71: daapp->bounddactx=0;
72: daapp->nbeforemonitors=0;
73: daapp->naftermonitors=0;
74: PetscInt ii;
75: info=TaoAppAddObject(daapplication,PPetscDAAppXXX,(void*)daapp,&ii); CHKERRQ(info);
76: info=TaoAppSetDestroyRoutine(daapplication,(int(*)(void*))(TaoAppDestroyDAApp), (void*)daapplication); CHKERRQ(info);
77: info=TaoAppSetOptionsRoutine(daapplication,DAAppSetOptions); CHKERRQ(info);
79: return(0);
80: }
84: /*@
85: TaoAppSetDAApp - Extend the functionality of a TaoApplication by
86: adding support for applications defined
87: a PETSc Distributed Array (DA) object. The routines used to evaluate
88: the objective function, constraint, and derivative information are
89: designed specifically for these problems.
91: Input Parameters:
92: + daappliation - an existing application object
93: . tda - an array of Distributed Array objects
94: - nnda - the number of Distibuted Array objects
96: .seealso TaoDAAppSolve(), TaoApplicationCreate();
98: Level: beginner
100: .keywords: Application, DA
101: @*/
102: int TaoAppSetDAApp(TAO_APPLICATION daapplication, DA* tda, PetscInt nnda){
103: int info;
104: PetscInt i;
105: DA_APPLICATION daapp;
106: PetscScalar zero=0.0;
108: if (nnda > PETSCDAAPPMAXGRIDS){
109: SETERRQ1(1,"Number of grids cannot exceed %d.\n",PETSCDAAPPMAXGRIDS); }
110: info = DAAppExtend(daapplication); CHKERRQ(info);
111: info = TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
112: info = PetscInfo(daapplication,"TAO Extending TAO_APPLICATION for DA_APPLICATION object.\n"); CHKERRQ(info);
113: info = DAAppSetMatType(daapplication,(const MatType)MATAIJ); CHKERRQ(info);
114: daapp->nda=nnda;
115: daapp->currentlevel=0;
116: for (i=0;i<nnda; i++){
118: info = PetscObjectReference((PetscObject)(tda[i])); CHKERRQ(info);
119: daapp->grid[i].da=tda[i];
120: info = DACreateGlobalVector(daapp->grid[i].da,&daapp->grid[i].X); CHKERRQ(info);
121: info = VecDuplicate(daapp->grid[i].X,&daapp->grid[i].R); CHKERRQ(info);
122: info = VecDuplicate(daapp->grid[i].X,&daapp->grid[i].RHS); CHKERRQ(info);
123: daapp->grid[i].XL=0;
124: daapp->grid[i].XU=0;
125: daapp->grid[i].H=0;
126: daapp->grid[i].Interpolate=0;
127: daapp->grid[i].CScale=0;
128: daapp->grid[i].mgrid=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 = PetscInfo(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: PetscInt i,nnda=daapp->nda;
176: int info;
181: if (--((PetscObject)daapp)->refct > 0) return(0);
183: info = PetscInfo(daapp,"Destroying work vectors for DA_APPLICATION object.\n"); CHKERRQ(info);
184: for (i=0;i<nnda; i++){
185: if (daapp->grid[i].coloring){
186: info = ISColoringDestroy(daapp->grid[i].coloring); CHKERRQ(info); daapp->grid[i].coloring=0;
187: }
188: if (daapp->grid[i].H){
189: info = MatDestroy(daapp->grid[i].H); CHKERRQ(info);daapp->grid[i].H=0;
190: }
191: if (daapp->grid[i].XL && daapp->grid[i].XU){
192: info = VecDestroy(daapp->grid[i].XL); CHKERRQ(info);daapp->grid[i].XL=0;
193: info = VecDestroy(daapp->grid[i].XU); CHKERRQ(info);daapp->grid[i].XU=0;
194: }
195: if (i>0){
196: info = MatDestroy(daapp->grid[i].Interpolate); CHKERRQ(info);daapp->grid[i].Interpolate=0;
197: info = VecDestroy(daapp->grid[i].CScale); CHKERRQ(info);daapp->grid[i].CScale=0;
198: }
199: info = VecDestroy(daapp->grid[i].RHS); CHKERRQ(info);daapp->grid[i].RHS=0;
200: info = VecDestroy(daapp->grid[i].R); CHKERRQ(info);daapp->grid[i].R=0;
201: info = VecDestroy(daapp->grid[i].X); CHKERRQ(info);daapp->grid[i].X=0;
202: info = DADestroy(daapp->grid[i].da); CHKERRQ(info);daapp->grid[i].da=0;
203: }
204: daapp->nda=0;
205: daapp->currentlevel=0;
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, PetscInt 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, PetscInt *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, PetscInt *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;
297: PetscInt clevel;
298: Mat H=*M;
299: DA_APPLICATION daapp=(DA_APPLICATION)ctx;
304: PetscStackPush("TAO User DA Hessian");
305: clevel=daapp->currentlevel;
306:
307: // daapp->FDGrad=daapp->grid[i].RHS; /* Needed for finite differences gradient vector */
308: info = PetscInfo1(daapp,"TAO hessian evaluation at DA_APPLICATION object, level %d.\n",clevel); CHKERRQ(info);
309: info = (*daapp->computedahessian)(daapplication,daapp->grid[clevel].da,X,H,daapp->usrdahctx);CHKERRQ(info);
310: *flag=daapp->kspflag;
312: PetscStackPop;
314: return(0);
315: }
320: /*@
321: DAAppSetHessianMat - Creates the matrices used to store the Hessian at finer grid levels
323: Collective on TAO_APPLICATION
325: Input Parameters:
326: . daapplication - the DA Application object
328: Level: advanced
331: .keywords: DA, hessian
333: .seealso: DAAppSetHessianRoutine(), DAAppSetMatType();
335: @*/
336: int DAAppSetHessianMat(TAO_APPLICATION daapplication){
337: int info;
338: PetscInt i,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_NEW_NONZERO_LOCATION_ERR, PETSC_TRUE); CHKERRQ(info);
348: info = MatSetOption(daapp->grid[i].H,MAT_STRUCTURALLY_SYMMETRIC, PETSC_TRUE); CHKERRQ(info);
349: info = MatSetOption(daapp->grid[i].H,MAT_SYMMETRIC, PETSC_TRUE); CHKERRQ(info);
350: info = PetscInfo1(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_GLOBAL,daapp->HessianMatrixType,&(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 = PetscInfo1(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 = PetscInfo(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 = PetscInfo1(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 = PetscInfo(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 = PetscInfo1(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 = PetscInfo(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: PetscInt level=daapp->currentlevel;
617: int info;
621: info = PetscInfo1(daapp->grid[level].da,"Compute variable bounds in level %d of DA_APPLICATION object.\n",level); CHKERRQ(info);
622: info = PetscObjectReference((PetscObject)XL);
623: info = PetscObjectReference((PetscObject)XU);
624: daapp->grid[level].XL=XL;
625: daapp->grid[level].XU=XU;
626: if (daapp->computedabounds){
627: info = (*daapp->computedabounds)(daapplication,daapp->grid[level].da,XL,XU,daapp->bounddactx); CHKERRQ(info);
628: }
629: return(0);
630: }
635: /*@C
636: DAAppSetVariableBoundsRoutine - Set a routine that will evaluate the bounds
637: on the variables.
639: Collective on TAO_APPLICATION
641: Input Parameters:
642: + daapplication - the DA Application object
643: . bound - the function pointer for the bound evaluation routine
644: - ctx - the monitor context
646: Calling sequence of bound:
647: $ bound(TAO_APPLICATION daapplication, DA da,Vec xl, Vec xu, void *ctx);
649: + daapplication - the TAO_APPLICATION context
650: . da - the Distributed Array
651: . xl - vector of upper bounds
652: . xu - vector of lower bounds
653: - ctx - user-defined monitor context set from DAAppSetVariableBoundsRoutine()
655: Level: beginner
657: .keywords: DA, bounds
659: .seealso: TaoDAAppSolve();
661: @*/
662: int DAAppSetVariableBoundsRoutine(TAO_APPLICATION daapplication, int (*bound)(TAO_APPLICATION,DA,Vec,Vec, void*),void *ctx){
663: int info;
664: DA_APPLICATION daapp;
666: info = TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
667: info=TaoAppSetVariableBoundsRoutine(daapplication,DAApplicationEvaluateVariableBounds,(void*)daapp);CHKERRQ(info);
668: daapp->computedabounds=bound;
669: daapp->bounddactx=ctx;
670: info = PetscInfo(daapp,"Set variable bounds in DA_APPLICATION object.\n"); CHKERRQ(info);
671: return(0);
672: }
678: /*@C
679: DAAppSetConstraintRoutine - Set a routine that will evaluate the constraint
680: function on the given DA at the given point.
682: Collective on TAO_APPLICATION
684: Input Parameters:
685: + daapplication - the TAO Application object
686: . grad - the function pointer for the gradient evaluation routine
687: - ctx - the gradient context
689: Calling sequence of grad:
690: $ f(TAO_APPLICATION daapplication,DA da, Vec x,Vec r,void *ctx);
692: + daapplication - the DA_APPLICATION context
693: . da - the Distributed Array
694: . x - input vector
695: . r - constraint vector
696: - ctx - user defined gradient context set from DAAppSetGradientRoutine()
698: Level: intermediate
700: Options Database Key:
701: . -tao_view - view the gradient after each evaluation using PETSC_VIEWER_STDOUT_SELF
703: .keywords: DA, gradient
705: .seealso: DAAppSetJacobianRoutine();
706: @*/
707: int DAAppSetConstraintRoutine(TAO_APPLICATION daapplication, int (*f)(TAO_APPLICATION,DA,Vec,Vec, void*),void *ctx){
708: int info;
709: DA_APPLICATION daapp;
711: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
712: if (f){
713: info=TaoAppSetConstraintRoutine(daapplication, DAAppPetscGradient, (void*)daapp);CHKERRQ(info);
714: } else {
715: info=TaoAppSetConstraintRoutine(daapplication,0,0);CHKERRQ(info);
716: }
717: daapp->IsComplementarity=PETSC_TRUE;
718: info=DAAppSetGradientRoutine(daapplication,f,ctx); CHKERRQ(info);
719: daapp->computedagradient=f;
720: daapp->usrdagctx=ctx;
721: daapp->usrdafgctx=ctx;
723: info = PetscInfo(daapp,"Set constraint function pointer for DA_APPLICATION object.\n"); CHKERRQ(info);
724: return(0);
725: }
730: /*@C
731: DAAppSetJacobianRoutine - Set a routine that will evaluate the Jacobian
732: on the given DA at the given point.
734: Collective on TAO_APPLICATION
736: Input Parameters:
737: + daapplication - the DA Application object
738: . jac - the function pointer for the Jacobian evaluation routine
739: - ctx - the jacobian context
741: Calling sequence of hess:
742: $ jac(TAO_APPLICATION daapplication, DA da, Vec x,Mat J,void *ctx);
744: + daapplication - the TAO_APPLICATION context
745: . da - the Distributed Array
746: . x - input vector
747: . J - Jacobian matrix
748: - ctx - user-defined hessian context set from DAAppSetHessianRoutine()
750: Level: intermediate
752: Options Database Key:
753: . -tao_view_jacobian - view the jacobian after each evaluation using PETSC_VIEWER_STDOUT_WORLD
755: .keywords: DA, hessian
757: .seealso: DAAppSetConstraintRoutine();
759: @*/
760: int DAAppSetJacobianRoutine(TAO_APPLICATION daapplication, int (*jac)(TAO_APPLICATION,DA,Vec,Mat,void*),void *ctx){
761: int info;
762: PetscInt i,nnda;
763: DA_APPLICATION daapp;
765: info = TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
766: if (!jac){ SETERRQ(1,"No DA Hessian routine provided\n");}
767: daapp->computedahessian=jac;
768: daapp->usrdahctx=ctx;
769: nnda=daapp->nda;;
770: for (i=0;i<nnda; i++){
771: if (daapp->grid[i].H==0){
772: info = DAGetMatrix(daapp->grid[i].da,daapp->HessianMatrixType,&daapp->grid[i].H); CHKERRQ(info);
773: info = MatSetOption(daapp->grid[i].H,MAT_NEW_NONZERO_LOCATION_ERR,PETSC_TRUE); CHKERRQ(info);
774: info = MatSetOption(daapp->grid[i].H,MAT_STRUCTURALLY_SYMMETRIC,PETSC_TRUE); CHKERRQ(info);
775: info = PetscInfo1(daapp->grid[i].da,"TAO create Jacobian for DA_APPLICATION object, level %d.\n",i); CHKERRQ(info);
776: }
777: }
778:
779: info = TaoAppSetJacobianRoutine(daapplication,DAAppPetscHessian,(void*)daapp);CHKERRQ(info);
780: info = TaoAppSetJacobianMat(daapplication,daapp->grid[0].H,daapp->grid[0].H);CHKERRQ(info);
781:
782: return(0);
783: }
787: /*@C
788: DAAppSetBeforeMonitor - Set a routine that will be called before the
789: optimization on each grid.
791: Collective on TAO_APPLICATION
793: Input Parameters:
794: + daapplication - the DA Application object
795: . beforemonitor - a monitor routine called before solving the application on each DA
796: - ctx - the monitor context
798: Calling sequence of monitor:
799: $ beforemonitor(TAO_APPLICATION daapplication, DA da,int level, void *ctx);
801: + daapplication - this TAO_APPLICATION context
802: . da - the Distributed Array
803: . level - the grid level that will be solved next (level 0 is coarsest)
804: - ctx - user-defined function context set from DAAppSetBeforeMonitor()
806: Note:
807: These monitors are different from the monitors that can be called after
808: each iteration of the optimization algorithm.
810: Note:
811: The beforemonitor and aftermonitor are good for setting up and destroying the application data.
813: Level: intermediate
815: .keywords: DA, monitor
817: .seealso: DAAppSetAfterMonitor(), TaoSetMonitor(), TaoAppSetMonitor();
819: @*/
820: int DAAppSetBeforeMonitor(TAO_APPLICATION daapplication, int (*beforemonitor)(TAO_APPLICATION,DA,PetscInt, void*), void *ctx){
821: PetscInt n;
822: int info;
823: DA_APPLICATION daapp;
825: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
826: if (beforemonitor){
827: n=daapp->nbeforemonitors;
828: if (n>=MAX_DAAP_MONITORS){
829: SETERRQ(1,"TAO ERROR: Too many beforemonitors set in DAAPP");
830: }
831: daapp->beforemonitor[n]=beforemonitor;
832: daapp->beforemonitorctx[n]=ctx;
833: daapp->nbeforemonitors++;
834: info = PetscInfo(daapplication,"TAO: Set user beforemonitor for DA_APPLICATION object.\n"); CHKERRQ(info);
835: }
836: return(0);
837: }
841: /*@C
842: DAAppSetAfterMonitor - Set a routine that will be called after the
843: optimization on each grid.
845: Collective on TAO_APPLICATION
847: Input Parameters:
848: + daapplication - the DA Application object
849: . aftermonitor - a monitor routine called after solving the application on each DA
850: - ctx - the monitor context
852: Calling sequence of monitor:
853: $ aftermonitor(TAO_APPLICATION daapplication, DA da, int level, void *ctx);
855: + daapplication - this TAO_APPLICATION context
856: . da - the Distributed Array
857: . level - the grid level that will be solved next (level 0 is coarsest)
858: - ctx - user-defined function context set from DAAppSetAfterMonitor()
860: Note:
861: These monitors are different from the monitors that can be called after
862: each iteration of the optimization algorithm.
864: Note:
865: The beforemonitor and aftermonitor are good for setting up and destroying the application data.
867: Level: intermediate
869: .keywords: DA, monitor
871: .seealso: DAAppSetBeforeMonitor(), TaoSetMonitor(), TaoAppSetMonitor();
873: @*/
874: int DAAppSetAfterMonitor(TAO_APPLICATION daapplication, int (*aftermonitor)(TAO_APPLICATION,DA,PetscInt, void*), void *ctx){
875: int info;
876: PetscInt n;
877: DA_APPLICATION daapp;
879: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
880: if (aftermonitor){
881: n=daapp->naftermonitors;
882: if (n>=MAX_DAAP_MONITORS){
883: SETERRQ(1,"TAO ERROR: Too many aftermonitors set in DAAPP");
884: }
885: daapp->aftermonitor[n]=aftermonitor;
886: daapp->aftermonitorctx[n]=ctx;
887: daapp->naftermonitors++;
888: info = PetscInfo(daapp,"TAO: Set user aftermonitor for DA_APPLICATION object.\n"); CHKERRQ(info);
889: }
890: return(0);
891: }
897: /*@
898: DAAppGetSolution - Sets a pointer to an existing vector that contains
899: the solution.
901: Collective on TAO_APPLICATION
903: Input Parameters:
904: + daapplication - the DA Application object
905: - level - the DA on which the current solution is desired (Level 0 is coarsest)
907: Output Parameters:
908: . X - the current solution
910: Level: beginner
912: .seealso: DAAppSetBeforeMonitor(), DAAppSetAfterMonitor();
914: .keywords: Solution, DA
915: @*/
916: int DAAppGetSolution(TAO_APPLICATION daapplication, PetscInt level, Vec *X){
917: int info;
918: DA_APPLICATION daapp;
920: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
921: if (level<0 || level>=daapp->nda){
922: SETERRQ(1,"No Solution available. This grid does not exist");
923: } else if (X && daapp && daapp->grid && level<daapp->nda){
924: *X=daapp->grid[level].X;
925: }
926: return(0);
927: }
932: /*@C
933: DAAppGetHessianMat - Sets a pointer to an existing matrix that contains
934: the Hessian at the specified level.
936: Collective on TAO_APPLICATION
938: Input Parameters:
939: + daapplication - the DA Application object
940: - level - the DA on which the current Hessian is desired (Level 0 is coarsest)
942: Output Parameters:
943: . H - the current Hessian
945: Level: intermediate
947: .seealso: DAAppSetHessianRoutine;
949: .keywords: Hessian, DA
950: @*/
951: int DAAppGetHessianMat(TAO_APPLICATION daapplication, PetscInt level, Mat *H){
952: int info;
953: DA_APPLICATION daapp;
955: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
956: if (level<0 || level>=daapp->nda){
957: SETERRQ(1,"No Hessian available. This grid does not exist");
958: } else if (H && daapp && daapp->grid && level<daapp->nda){
959: *H=daapp->grid[level].H;
960: }
961: return(0);
962: }
966: /*@
967: DAAppGetInterpolationMatrix - Sets pointers to existing vectors
968: containg upper and lower bounds on the variables.
970: Collective on TAO_APPLICATION
972: Input Parameters:
973: + daapplication - the DA Application object
974: - level - the DA on which the current solution is desired (Level 0 is coarsest)
976: Output Parameters:
977: + Interpolate - the interpolating matrix
978: - CScale - the scaling matrix for restriction Matrix
980: Level: intermediate
982: .seealso: DAAppSetBeforeMonitor(), DAAppSetAfterMonitor();
984: .keywords: Interpolate, DA
985: @*/
986: int DAAppGetInterpolationMatrix(TAO_APPLICATION daapplication, PetscInt level, Mat *Interpolate, Vec *CScale){
987: int info;
988: DA_APPLICATION daapp;
990: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
991: if ( !daapp || level<1 || level>=daapp->nda){
992: SETERRQ1(1,"No Interpolator available for level %d. This grid does not exist",level);
993: }
994: if (Interpolate){
995: *Interpolate=daapp->grid[level].Interpolate;
996: }
997: if (CScale){
998: *CScale=daapp->grid[level].CScale;
999: }
1000: return(0);
1001: }
1006: /*@
1007: DAAppGetVariableBounds - Sets pointers to existing vectors
1008: containg upper and lower bounds on the variables.
1010: Collective on TAO_APPLICATION
1012: Input Parameters:
1013: + daapplication - the DA Application object
1014: - level - the DA on which the current solution is desired (Level 0 is coarsest)
1016: Output Parameters:
1017: + XL - the current solution
1018: - XU - the current solution
1020: Level: intermediate
1022: .seealso: DAAppSetBeforeMonitor(), DAAppSetAfterMonitor();
1024: .keywords: Solution, DA
1025: @*/
1026: int DAAppGetVariableBounds(TAO_APPLICATION daapplication, PetscInt level, Vec *XL, Vec *XU){
1027: int info;
1028: DA_APPLICATION daapp;
1030: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
1031: if ( !daapp || level<0 || level>=daapp->nda){
1032: SETERRQ(1,"No Solution available. This grid does not exist.");
1033: }
1034: if (XL){
1035: *XL=daapp->grid[level].XL;
1036: }
1037: if (XU){
1038: *XU=daapp->grid[level].XU;
1039: }
1040: return(0);
1041: }
1047: /*@
1048: DAAppSetInitialSolution - Sets the initial solution for
1049: the grid application.
1051: Collective on TAO_APPLICATION
1053: Input Parameters:
1054: + daapplication - the DA Application object
1055: - X0 - the initial solution vector
1057: Level: beginner
1059: .seealso: DAAppGetSolution();
1061: .keywords: Initial Solution, DA
1062: @*/
1063: int DAAppSetInitialSolution(TAO_APPLICATION daapplication, Vec X0){
1064: int info;
1065: DA_APPLICATION daapp;
1068: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
1069: info=VecCopy(X0,daapp->grid[0].X); CHKERRQ(info);
1070: return(0);
1071: }
1077: /*
1078: DAAppSetMatStructure - Set the MatStructure flag to be used in KSPSetOperators:
1080: Collective on TAO_APPLICATION
1081:
1082: Input Parameters:
1083: + daapplication - the DA Application object
1084: - flag - indicates the KSP flag
1086: Level: intermediate
1088: Note:
1089: The default value is SAME_NONZERO_PATTERN
1091: .seealso: KSPSetOperators(), TaoAppGetKSP();
1093: .keywords: Linear Solver, Multigrid, DA, KSP
1095: @*/
1096: int DAAppSetMatStructure(TAO_APPLICATION daapplication, MatStructure pflag){
1097: int info;
1098: DA_APPLICATION daapp;
1101: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
1102: daapp->kspflag=pflag;
1103: info = PetscInfo(daapp,"Set MatStructure flag in TAO solver.\n"); CHKERRQ(info);
1105: return(0);
1106: }
1110: /*@
1111: DAAppSetMatType - set the matrix type to be used for the Hessian matrix
1113: Collective on TAO_APPLICATION
1114:
1115: Input Parameters:
1116: + daapplication - the DA Application object
1117: - mattype - the matrix type
1119: Level: intermediate
1121: Options Database Key:
1122: . -tao_da_mattype - select matrix type
1124: Note:
1125: The default value is MATMPIAIJ.
1127: Note:
1128: Call before DAAppSetHessianRoutine()
1130: .seealso: DAAppSetMatStructure(), DAAppSetHessianRoutine(), DAGetMatrix();
1132: .keywords: DA, Hessian
1134: @*/
1135: int DAAppSetMatType(TAO_APPLICATION daapplication, const MatType mattype){
1136: int info;
1137: DA_APPLICATION daapp;
1140: info=TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
1141: info=PetscStrcpy(daapp->HessianMatrixType,mattype); CHKERRQ(info);
1143: return(0);
1144: }
1149: /*@
1150: DAAppSetOptions - Sets various parameters to be used in this application
1151: and the TAO solver.
1153: Collective on TAO_APPLICATION
1155: Input Parameters:
1156: . daapplication - the DA Application object
1158: Level: intermediate
1160: Note: This routine will be called by the TaoAppSetFromOptions().
1162: .keywords: DA, options
1164: .seealso: TaoDAAppSolve();
1166: @*/
1167: int DAAppSetOptions(TAO_APPLICATION daapplication){
1168: int info;
1169: // char *prefix;
1170: PetscTruth flg1=PETSC_FALSE,flg2=PETSC_FALSE;
1171: DA_APPLICATION daapp;
1172: // MPI_Comm comm;
1176: info = PetscInfo(daapplication,"TaoDAAppSetOptions(): Reading command line for options\n"); CHKERRQ(info);
1177: info = TaoAppDAApp(daapplication,&daapp); CHKERRQ(info);
1178: /*
1179: info = PetscObjectGetOptionsPrefix((PetscObject)daapplication,&prefix); CHKERRQ(info);
1180: info = PetscObjectGetComm((PetscObject)daapplication,&comm); CHKERRQ(info);
1181: info = PetscOptionsBegin(comm,prefix,"TAO PETSC DA APPLICATIONS","solver");CHKERRQ(info);
1182: */
1183: info = DAAppSetMultiGridOptions(daapplication);CHKERRQ(info);
1186: flg1=PETSC_FALSE,flg2=PETSC_FALSE;
1187: info = PetscOptionsTruth("-tao_da_samepreconditioner","Use same preconditioner in KSP","DAAppSetMatStructure",
1188: PETSC_FALSE,&flg2,&flg1);CHKERRQ(info);
1189: if (flg1 && flg2) {
1190: info = DAAppSetMatStructure(daapplication, SAME_PRECONDITIONER);CHKERRQ(info);
1191: } else if (flg1){
1192: info = DAAppSetMatStructure(daapplication, SAME_NONZERO_PATTERN);CHKERRQ(info);
1193: }
1195: info = PetscOptionsString("-tao_da_mattype","the hessian should have matrix type","DAAppSetMatType",
1196: daapp->HessianMatrixType,daapp->HessianMatrixType,18,&flg1);CHKERRQ(info);
1198: /*
1199: info = PetscOptionsEnd();CHKERRQ(info);
1200: */
1201: return(0);
1202: }