Actual source code: taoapp_util.c
1: #include "taoapp_petsc.h" /*I "tao.h" I*/
2: #include tao.h
4: static int Tao_Solve=0;
9: /*@
10: TaoSolveApplication - Find a solution to the application using the set TAO solver.
12: Collective on TAO_APPLICATION
14: Input Parameters:
15: + taoapp - the TAO_APPLICATION context
16: - tao - the TAO_SOLVER context
18: Level: beginner
20: .keywords: solve
22: .seealso: TaoAppGetSolutionVec()
24: @*/
25: int TaoSolveApplication(TAO_APPLICATION taoapp, TAO_SOLVER tao){
26: int info;
29: info = TaoSetupApplicationSolver(taoapp,tao);CHKERRQ(info);
30: info = PetscLogEventBegin(Tao_Solve,tao,0,0,0);CHKERRQ(info);
31: info = TaoSolve(tao); CHKERRQ(info);
32: info = PetscLogEventEnd(Tao_Solve,tao,0,0,0);CHKERRQ(info);
33: // info = TaoGetSolutionStatus(tao,&taoapp->iter,&taoapp->fval,&taoapp->residual,0,0,&taoapp->reason);
34: CHKERRQ(info);
35: return(0);
36: }
40: /*@
41: TaoSetupApplicationSolver - This routine creates the vectors,
42: matrices, linear solvers, and other data structures used in
43: the during the optimization process. The application provides
44: the solver with an objective function, constraints, derivative
45: information, and application data structures. These structures
46: include a vector of variables, and Hessian matrix.
48: Collective on TAO_SOLVER
50: Input Parameters:
51: + myapp - user application context
52: - tao - the TAO_SOLVER solver context
54: Level: intermediate
56: Note:
57: This routine should be called before TaoGetKSP(), but after
58: TaoAppSetInitialSolutionVec() and after TaoAppSetHessianMat() (when Newton solvers are used).
60: Note:
61: This method is called during TaoSetOptions() and TaoSolveApplication()
62:
63: .keywords: application, context
65: @*/
66: int TaoSetupApplicationSolver(TAO_APPLICATION myapp, TAO_SOLVER tao ){
67: int info;
68: TaoPetscApplication* taopetscapp;
71: if (Tao_Solve==0){
72: info = PetscLogEventRegister(&Tao_Solve,"TaoSolve",TAO_APP_COOKIE); CHKERRQ(info);
73: }
74: info = TaoAppGetTaoPetscApp(myapp,&taopetscapp);
75: info = TaoSetApplication(tao,taopetscapp);CHKERRQ(info);
76: taopetscapp->tao=tao;
77: return(0);
78: }
82: /*@
83: TaoSetOptions - Sets various TAO parameters from user options
85: Collective on TAO_APPLICATION
87: Input Parameters:
88: + taoapp - the TAO Application (optional)
89: - tao - the TAO optimization solver (optional)
90: Level: beginner
92: Note:
93: This routine should be called after TaoSetupApplicationSolver()
95: Note:
96: This routine must be called if there are multiple processors involved and
97: the MPI Communicator is different than MPI_COMM_WORLD.
99: .keywords: options
101: .seealso: TaoSolveApplication()
103: @*/
104: int TaoSetOptions(TAO_APPLICATION taoapp, TAO_SOLVER tao){
105: int info;
106: const char *prefix=0;
107: PetscTruth flg;
108: KSP ksp;
109: MPI_Comm comm=MPI_COMM_WORLD;
113: if (tao){
115: info = PetscObjectGetOptionsPrefix((PetscObject)tao,&prefix); CHKERRQ(info);
116: info = PetscObjectGetComm((PetscObject)tao,&comm);CHKERRQ(info);
117: info = PetscOptionsBegin(comm,prefix,"TAO PETSC APPLICATIONS ","solver");CHKERRQ(info);
119: info = TaoSetFromOptions(tao); CHKERRQ(info);
121: flg=PETSC_FALSE;
122: info = PetscOptionsName("-tao_xmonitor","Use graphics convergence","TaoPetscXMonitor",&flg);CHKERRQ(info);
123: if (flg){
124: info = TaoSetPetscXMonitor(tao); CHKERRQ(info);
125: }
127: info = PetscOptionsEnd();CHKERRQ(info);
128: }
130: if (taoapp){
131: info = TaoAppSetFromOptions(taoapp); CHKERRQ(info);
132: }
134: if (tao && taoapp){
135: info = TaoSetupApplicationSolver(taoapp,tao);CHKERRQ(info);
136: info = TaoGetKSP(tao,&ksp);CHKERRQ(info);
137: if (ksp){
138: info = KSPSetFromOptions(ksp);CHKERRQ(info);
139: }
140: info = PetscOptionsName("-tao_lmvmh","User supplies approximate hessian for LMVM solvers","TaoLMVMSetH0",&flg);
141: if (flg){
142: info=TaoBLMVMSetH0(tao,TAO_TRUE);CHKERRQ(info);
143: info=TaoLMVMSetH0(tao,TAO_TRUE);CHKERRQ(info);
144: }
145: }
146:
147: return(0);
148: }
149:
152: static char TaoPetscAppXXX[] = "TaoPetscApp";
156: static int TaoAppDestroyPetscAppXXX(void*ctx){
157: TaoPetscApplication *mctx=(TaoPetscApplication*)ctx;
159: if (mctx){
160: delete mctx;
161: }
162: return(0);
163: }
167: int TaoAppGetTaoPetscApp(TAO_APPLICATION taoapp,TaoPetscApplication**tpapp){
168: int ii,info;
169: MPI_Comm comm;
170: TaoPetscApplication*ttpapp;
174: info=TaoAppQueryForObject(taoapp,TaoPetscAppXXX,(void**)&ttpapp); CHKERRQ(info);
175: if (ttpapp==0){
176: info=PetscObjectGetComm((PetscObject)taoapp,&comm); CHKERRQ(info);
177: ttpapp=new TaoPetscApplication(comm);
178: info=TaoAppAddObject(taoapp,TaoPetscAppXXX,(void*)ttpapp,&ii); CHKERRQ(info);
179: info=TaoAppSetDestroyRoutine(taoapp,TaoAppDestroyPetscAppXXX, (void*)ttpapp); CHKERRQ(info);
180: }
181: ttpapp->papp=taoapp;
182: *tpapp=ttpapp;
183: return(0);
184: }
188: /*@C
189: TaoGetKSP - Gets the linear solver used by the optimization solver
191: Input Parameters:
192: . tao - the TAO solver
194: Output Parameters:
195: . ksp - the KSP linear solver used in the optimization solver
197: Level: intermediate
199: .keywords: Application
200: @*/
201: int TaoGetKSP(TAO_SOLVER tao, KSP *ksp){
202: int info;
203: KSP newksp;
204: TaoLinearSolver *tsolver;
205:
208: if (ksp){
209: *ksp=0;
210: info = TaoGetLinearSolver(tao,&tsolver); CHKERRQ(info);
211: if (tsolver){
212: newksp=(KSP)(tsolver->LinearSolverObject);
213: *ksp=newksp;
214: }
215: }
216: return(0);
217: }
222: int TaoAppGetTaoSolver(TAO_APPLICATION taoapp, TAO_SOLVER *tao){
223: int info;
224: TaoPetscApplication* taopetscapp;
228: info = TaoAppGetTaoPetscApp(taoapp,&taopetscapp); CHKERRQ(info);
229: if (tao){ *tao=taopetscapp->tao; }
230: return(0);
231: }
235: int TaoAppGetKSP(TAO_APPLICATION taoapp, KSP *ksp){
236: int info;
237: TAO_SOLVER tao;
240: info=TaoAppGetTaoSolver(taoapp,&tao); CHKERRQ(info);
241: if (tao){
242: info=TaoGetKSP(tao,ksp); CHKERRQ(info);
243: } else if (ksp){
244: *ksp=0;
245: }
246: return(0);
247: }
251: /*@C
252: TaoGetVariableBoundVecs - Get the vectors with the
253: lower and upper bounds in current solver.
255: Input Parameters:
256: . tao - the TAO solver
258: Output Parameter:
259: + XL - the lower bounds
260: - XU - the upper bounds
262: Level: intermediate
264: Note:
265: These vectors should not be destroyed.
267: .keywords: Application, bounds
269: .seealso: TaoAppGetGradientVec(), TaoAppGetSolutionVec(), TaoAppSetVariableBoundsRoutine()
270: @*/
271: int TaoGetVariableBoundVecs(TAO_SOLVER tao, Vec *XL, Vec *XU){
272: int info;
273: TaoVec *XXLL,*XXUU;
276: info = TaoGetVariableBounds(tao,&XXLL,&XXUU); CHKERRQ(info);
277: if (XL){
278: *XL=0;
279: if (XXLL){
280: *XL=(Vec)(XXLL->VecObject);
281: }
282: }
283: if (XU){
284: *XU=0;
285: if (XXUU){
286: *XU=(Vec)(XXUU->VecObject);
287: }
288: }
289: return(0);
290: }
294: /*@C
295: TaoAppGetGradientVec - Get the vector with the
296: gradient in the current application.
298: Input Parameters:
299: . tao - the solver
301: Output Parameter:
302: . G - the gradient vector
304: Level: intermediate
306: Note:
307: This vector should not be destroyed.
309: .keywords: Application, gradient
311: .seealso: TaoAppGetSolutionVec()
312: @*/
313: int TaoAppGetGradientVec(TAO_SOLVER tao, Vec *G){
314: int info;
315: TaoVec* GG;
318: info = TaoGetGradient(tao,&GG); CHKERRQ(info);
319: if (G&&GG) { *G=(Vec)GG->VecObject;}
320: return(0);
321: }
326: /*@
327: TaoCopyDualsOfVariableBounds - Copy the current dual variables
328: corresponding the lower and upper bounds on the variables.
329:
330: Input Parameters:
331: . tao - the solver
333: Output Parameter:
334: + DXL - the lower bounds
335: - DXU - the upper bounds
337: Level: intermediate
339: Note:
340: Existing vectors of commensurate distribution to the
341: variable bounds should be passed into this routine.
343: Note:
344: These variables may have to be computed. It may not be efficient
345: to call this routine in a Monitor.
347: Note:
348: These variables can be interpreted as the sensitivity of
349: the objective value with respect to the bounds on the variables.
351: .keywords: Application, bounds, duals
353: .seealso: TaoAppGetGradientVec(), TaoAppGetSolutionVec(), TaoAppSetVariableBoundsRoutine()
354: @*/
355: int TaoCopyDualsOfVariableBounds(TAO_SOLVER tao, Vec DXL, Vec DXU){
356: int info;
357: TaoVecPetsc *ddxl,*ddxu;
360: info = TaoWrapPetscVec(DXL,&ddxl); CHKERRQ(info);
361: info = TaoWrapPetscVec(DXU,&ddxu); CHKERRQ(info);
362: info = TaoGetDualVariables(tao, ddxl, ddxu); CHKERRQ(info);
363: info = TaoVecDestroy(ddxl); CHKERRQ(info);
364: info = TaoVecDestroy(ddxu); CHKERRQ(info);
365: return(0);
366: }
372: /*@C
373: TaoSetInequality Constraints - Set inequality constraints for OOQP
375: Collective on TAO_SOLVER
377: Input Parameters:
378: + taoapp - the TAO_APPLICATION context
379: . ll - vector to store lower bounds
380: . uu - vector to store upper bounds
381: - AA - matrix to store linear inequality constraints
383: Level: intermediate
385: .keywords: TAO_SOLVER, inequalities
387: @*/
388: int TaoSetInequalityConstraints(TAO_APPLICATION taoapp, Vec ll, Mat A, Vec uu){
390: int info;
391: TaoPetscApplication* taopetscapp;
393: info = TaoAppGetTaoPetscApp(taoapp,&taopetscapp);
394: info = TaoVecDestroy(taopetscapp->taofvll); CHKERRQ(info); taopetscapp->taofvll=0;
395: info = TaoWrapPetscVec(ll,&taopetscapp->taofvll); CHKERRQ(info);
396: info = TaoVecDestroy(taopetscapp->taofvuu); CHKERRQ(info); taopetscapp->taofvuu=0;
397: info = TaoWrapPetscVec(uu,&taopetscapp->taofvuu); CHKERRQ(info);
398: info = TaoMatDestroy(taopetscapp->taoAA); CHKERRQ(info); taopetscapp->taoAA=0;
399: info = TaoWrapPetscMat(A,&taopetscapp->taoAA); CHKERRQ(info);
400: return(0);
401: }