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;

  7: int TaoAppGetTaoPetscApp(TAO_APPLICATION taoapp,TaoPetscApplication**tpapp);

 11: /*@
 12:   TaoSolveApplication - Find a solution to the application using the set TAO solver.

 14:    Collective on TAO_APPLICATION

 16:    Input Parameters:
 17: +  taoapp - the TAO_APPLICATION context
 18: -  tao - the TAO_SOLVER context

 20:    Level: beginner

 22: .keywords: solve

 24: .seealso: TaoAppGetSolutionVec()

 26: @*/
 27: int TaoSolveApplication(TAO_APPLICATION taoapp, TAO_SOLVER tao){
 28:   int info;

 31:   info = TaoSetupApplicationSolver(taoapp,tao);CHKERRQ(info);
 32:   info = PetscLogEventBegin(Tao_Solve,tao,0,0,0);CHKERRQ(info);
 33:   info = TaoSolve(tao); CHKERRQ(info);
 34:   info = PetscLogEventEnd(Tao_Solve,tao,0,0,0);CHKERRQ(info);
 35:   //  info = TaoGetSolutionStatus(tao,&taoapp->iter,&taoapp->fval,&taoapp->residual,0,0,&taoapp->reason); 
 36:   CHKERRQ(info);
 37:   return(0);
 38: }

 42: /*@
 43:    TaoSetupApplicationSolver - This routine creates the vectors,
 44:    matrices, linear solvers, and other data structures used in
 45:    the during the optimization process.  The application provides
 46:    the solver with an objective function, constraints, derivative 
 47:    information, and application data structures.  These structures
 48:    include a vector of variables, and Hessian matrix.

 50:    Collective on TAO_SOLVER

 52:    Input Parameters:
 53: +  myapp - user application context
 54: -  tao - the TAO_SOLVER solver context

 56:    Level: intermediate

 58:    Note: 
 59:    This routine should be called before TaoGetKSP(), but after 
 60:    TaoAppSetInitialSolutionVec() and after TaoAppSetHessianMat() (when Newton solvers are used). 

 62:    Note: 
 63:    This method is called during  TaoSetOptions() and TaoSolveApplication()
 64:    
 65: .keywords: application, context

 67: @*/
 68: int TaoSetupApplicationSolver(TAO_APPLICATION myapp, TAO_SOLVER tao ){
 69:   int info;
 70:   TaoPetscApplication* taopetscapp;
 73:   if (Tao_Solve==0){
 74:     info = PetscLogEventRegister(&Tao_Solve,"TaoSolve",TAO_APP_COOKIE); CHKERRQ(info);
 75:   }
 76:   info = TaoAppGetTaoPetscApp(myapp,&taopetscapp);
 77:   info = TaoSetApplication(tao,taopetscapp);CHKERRQ(info);
 78:   taopetscapp->tao=tao;
 79:   return(0);
 80: }

 84: /*@
 85:   TaoSetOptions - Sets various TAO parameters from user options

 87:    Collective on TAO_APPLICATION

 89:    Input Parameters:
 90: +  taoapp - the TAO Application (optional)
 91: -  tao - the TAO optimization solver (optional)
 92:    Level: beginner

 94:    Note: 
 95:    This routine should be called after TaoSetupApplicationSolver()

 97:    Note: 
 98:    This routine must be called if there are multiple processors involved and 
 99:    the MPI Communicator is different than MPI_COMM_WORLD.

101: .keywords:  options

103: .seealso: TaoSolveApplication()

105: @*/
106: int TaoSetOptions(TAO_APPLICATION taoapp, TAO_SOLVER tao){
107:   int info;
108:   const char *prefix=0;
109:   PetscTruth flg;
110:   MPI_Comm comm=MPI_COMM_WORLD;


114:   if (tao){
116:     info = PetscObjectGetOptionsPrefix((PetscObject)tao,&prefix); CHKERRQ(info);
117:     info = PetscObjectGetComm((PetscObject)tao,&comm);CHKERRQ(info);
118:     info = PetscOptionsBegin(comm,prefix,"TAO PETSC APPLICATIONS ","solver");CHKERRQ(info);

120:     info = TaoSetFromOptions(tao); CHKERRQ(info);

122:     flg=PETSC_FALSE;
123:     info = PetscOptionsName("-tao_xmonitor","Use graphics convergence","TaoPetscXMonitor",&flg);CHKERRQ(info);
124:     if (flg){
125:       info = TaoSetPetscXMonitor(tao); CHKERRQ(info);
126:     }

128:     info = PetscOptionsEnd();CHKERRQ(info);
129:   }

131:   if (taoapp){
132:     info = TaoAppSetFromOptions(taoapp); CHKERRQ(info);
133:   }

135:   if (tao && taoapp){
136:     info = TaoSetupApplicationSolver(taoapp,tao);CHKERRQ(info);
137:     info = PetscOptionsName("-tao_lmvmh","User supplies approximate hessian for LMVM solvers","TaoLMVMSetH0",&flg);
138:     if (flg){
139:       info=TaoBLMVMSetH0(tao,TAO_TRUE);CHKERRQ(info);
140:       info=TaoLMVMSetH0(tao,TAO_TRUE);CHKERRQ(info);
141:     }
142:   }
143:   
144:   return(0);
145: }
146:  


149: static char TaoPetscAppXXX[] = "TaoPetscApp";

153: static int TaoAppDestroyPetscAppXXX(void*ctx){
154:   TaoPetscApplication *mctx=(TaoPetscApplication*)ctx;
156:   if (mctx){
157:     delete mctx;
158:   }
159:   return(0);
160: }

164: int TaoAppGetTaoPetscApp(TAO_APPLICATION taoapp,TaoPetscApplication**tpapp){
165:   int ii,info;
166:   MPI_Comm comm;
167:   TaoPetscApplication*ttpapp;

171:   info=TaoAppQueryForObject(taoapp,TaoPetscAppXXX,(void**)&ttpapp); CHKERRQ(info);
172:   if (ttpapp==0){
173:     info=PetscObjectGetComm((PetscObject)taoapp,&comm); CHKERRQ(info);
174:     ttpapp=new TaoPetscApplication(comm);
175:     info=TaoAppAddObject(taoapp,TaoPetscAppXXX,(void*)ttpapp,&ii); CHKERRQ(info);
176:     info=TaoAppSetDestroyRoutine(taoapp,TaoAppDestroyPetscAppXXX, (void*)ttpapp); CHKERRQ(info);
177:   }
178:   ttpapp->papp=taoapp;
179:   *tpapp=ttpapp;
180:   return(0);
181: }

185: /*@C
186:   TaoGetKSP - Gets the linear solver used by the optimization solver.
187:   Application writers should use TaoAppGetKSP if they need direct access
188:   to the PETSc KSP object.
189:   
190:    Input Parameters:
191: .  tao - the TAO solver

193:    Output Parameters:
194: .  ksp - the KSP linear solver used in the optimization solver

196:    Level: developer

198: .keywords: Application

200: .seealso: TaoAppGetKSP()

202: @*/
203: int TaoGetKSP(TAO_SOLVER tao, KSP *ksp)
204: {
205:   int info;
206:   TaoLinearSolver *tsolver;
207:   
210:   if (ksp){
211:     *ksp=0;
212:     info = TaoGetLinearSolver(tao,&tsolver); CHKERRQ(info);
213:     if (tsolver){
214:       TaoLinearSolverPetsc *psolver = dynamic_cast <TaoLinearSolverPetsc *> (tsolver);
215:       *ksp=psolver->GetKSP();
216:     }
217:   } 
218:   return(0);
219: }


224: int TaoAppGetTaoSolver(TAO_APPLICATION taoapp, TAO_SOLVER *tao){
225:   int info;
226:   TaoPetscApplication* taopetscapp;

230:   info = TaoAppGetTaoPetscApp(taoapp,&taopetscapp); CHKERRQ(info);
231:   if (tao){ *tao=taopetscapp->tao; }
232:   return(0);
233: }


238: /*@C
239:   TaoGetVariableBoundVecs - Get the vectors with the
240:   lower and upper bounds in current solver.

242:    Input Parameters:
243: .  tao - the TAO solver

245:    Output Parameter:
246: +  XL - the lower bounds
247: -  XU - the upper bounds

249:    Level: intermediate

251:    Note: 
252:    These vectors should not be destroyed.

254: .keywords: Application, bounds

256: .seealso: TaoAppGetGradientVec(), TaoAppGetSolutionVec(), TaoAppSetVariableBoundsRoutine()
257: @*/
258: int TaoGetVariableBoundVecs(TAO_SOLVER tao, Vec *XL, Vec *XU){
259:   int info;
260:   TaoVec *XXLL,*XXUU;
263:   info = TaoGetVariableBounds(tao,&XXLL,&XXUU); CHKERRQ(info);
264:   if (XL){
265:     *XL=0;
266:     if (XXLL){
267:       TaoVecPetsc *pl = dynamic_cast <TaoVecPetsc *> (XXLL);
268:       *XL = pl->GetVec();
269:     }
270:   }
271:   if (XU){
272:     *XU=0;
273:     if (XXUU){
274:       TaoVecPetsc *pu = dynamic_cast <TaoVecPetsc *> (XXUU);
275:       *XU = pu->GetVec();
276:     }
277:   }
278:   return(0);
279: }

283: /*@C
284:   TaoAppGetGradientVec - Get the vector with the
285:   gradient in the current application.

287:    Input Parameters:
288: .  tao - the solver

290:    Output Parameter:
291: .  G - the gradient vector

293:    Level: intermediate

295:    Note: 
296:    This vector should not be destroyed.

298: .keywords: Application, gradient

300: .seealso:  TaoAppGetSolutionVec()
301: @*/
302: int TaoAppGetGradientVec(TAO_SOLVER tao, Vec *G){
303:   int info;
304:   TaoVec* GG;
307:   info = TaoGetGradient(tao,&GG); CHKERRQ(info);
308:   if (G&&GG) {
309:     TaoVecPetsc *pg = dynamic_cast <TaoVecPetsc *> (GG);
310:     *G = pg->GetVec();
311:   }
312:   return(0);
313: }

317: /*@C
318:   TaoAppGetGessianMat - Get the vector with the
319:   gradient in the current application.

321:    Input Parameters:
322: .  tao - the solver

324:    Output Parameter:
325: .  H - the hessian matrix
326: .  Hpre - the hessian preconditioner

328:    Level: intermediate

330:    Note: 
331:    These matrices should not be destroyed.

333: .keywords: Application, hessian

335: .seealso:  TaoAppGetSolutionVec(), TaoAppGetGradientVec()
336: @*/
337: int TaoAppGetHessianMat(TAO_SOLVER tao, Mat *H, Mat *Hpre){
338:   int info;
339:   TaoMat* TM=0;
340:   MatStructure flag;
343:   info = TaoGetHessian(tao,&TM); CHKERRQ(info);
344:   if (H && TM) {
345:     TaoMatPetsc *tmp = dynamic_cast <TaoMatPetsc *> (TM);
346:     info = tmp->GetMatrix(H, Hpre, &flag); CHKERRQ(info);
347:   }
348:   return(0);
349: }


354: /*@
355:   TaoCopyDualsOfVariableBounds - Copy the current dual variables
356:   corresponding the lower and upper bounds on the variables.
357:   
358:    Input Parameters:
359: .  tao - the solver

361:    Output Parameter:
362: +  DXL - the lower bounds
363: -  DXU - the upper bounds

365:    Level: intermediate

367:    Note:  
368:    Existing vectors of commensurate distribution to the
369:    variable bounds should be passed into this routine.

371:    Note: 
372:    These variables may have to be computed.  It may not be efficient
373:    to call this routine in a Monitor.

375:    Note: 
376:    These variables can be interpreted as the sensitivity of
377:    the objective value with respect to the bounds on the variables.

379: .keywords: Application, bounds, duals

381: .seealso: TaoAppGetGradientVec(), TaoAppGetSolutionVec(), TaoAppSetVariableBoundsRoutine()
382: @*/
383: int TaoCopyDualsOfVariableBounds(TAO_SOLVER tao, Vec DXL, Vec DXU){
384:   int info;
385:   TaoVecPetsc *ddxl,*ddxu;
388:   info = TaoWrapPetscVec(DXL,&ddxl); CHKERRQ(info);
389:   info = TaoWrapPetscVec(DXU,&ddxu); CHKERRQ(info);
390:   info = TaoGetDualVariables(tao, ddxl, ddxu); CHKERRQ(info);
391:   info = TaoVecDestroy(ddxl); CHKERRQ(info);
392:   info = TaoVecDestroy(ddxu); CHKERRQ(info);
393:   return(0);
394: }



400: /*@C
401:    TaoSetInequality Constraints - Set inequality constraints for OOQP

403:    Collective on TAO_SOLVER

405:    Input Parameters:
406: +  taoapp - the TAO_APPLICATION context
407: .  ll - vector to store lower bounds
408: .  uu - vector to store upper bounds
409: -  AA - matrix to store linear inequality constraints

411:    Level: intermediate

413: .keywords: TAO_SOLVER, inequalities

415: @*/
416: int TaoSetInequalityConstraints(TAO_APPLICATION taoapp, Vec ll, Mat A, Vec uu){

418:   int info;
419:   TaoPetscApplication* taopetscapp;
421:   info = TaoAppGetTaoPetscApp(taoapp,&taopetscapp);
422:   info = TaoVecDestroy(taopetscapp->taofvll); CHKERRQ(info); taopetscapp->taofvll=0;
423:   info = TaoWrapPetscVec(ll,&taopetscapp->taofvll); CHKERRQ(info);
424:   info = TaoVecDestroy(taopetscapp->taofvuu); CHKERRQ(info); taopetscapp->taofvuu=0;
425:   info = TaoWrapPetscVec(uu,&taopetscapp->taofvuu); CHKERRQ(info);
426:   info = TaoMatDestroy(taopetscapp->taoAA); CHKERRQ(info); taopetscapp->taoAA=0;
427:   info = TaoWrapPetscMat(A,&taopetscapp->taoAA); CHKERRQ(info);
428:   return(0);
429: }