Actual source code: taoapp_util.c

  1: #include "tao_general.h"
  2: #include "taoapp_petsc.h"  /*I  "tao.h"  I*/
  3: #include "tao.h"

  5: static int Tao_Solve=0;

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

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

 15:    Collective on TAO_APPLICATION

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

 21:    Level: beginner

 23: .keywords: solve

 25: .seealso: TaoAppGetSolutionVec()

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

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

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

 51:    Collective on TAO_SOLVER

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

 57:    Level: intermediate

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

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

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

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

 88:    Collective on TAO_APPLICATION

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

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

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

102: .keywords:  options

104: .seealso: TaoSolveApplication()

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


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

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

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

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

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

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


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

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

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

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

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

195:    Output Parameters:
196: .  ksp - the KSP linear solver used in the optimization solver

198:    Level: developer

200: .keywords: Application

202: .seealso: TaoAppGetKSP()

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


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

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


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

244:    Input Parameters:
245: .  tao - the TAO solver

247:    Output Parameter:
248: +  XL - the lower bounds
249: -  XU - the upper bounds

251:    Level: intermediate

253:    Note: 
254:    These vectors should not be destroyed.

256: .keywords: Application, bounds

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

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

289:    Input Parameters:
290: .  tao - the solver

292:    Output Parameter:
293: .  G - the gradient vector

295:    Level: intermediate

297:    Note: 
298:    This vector should not be destroyed.

300: .keywords: Application, gradient

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

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

323:    Input Parameters:
324: .  tao - the solver

326:    Output Parameter:
327: .  H - the hessian matrix
328: .  Hpre - the hessian preconditioner

330:    Level: intermediate

332:    Note: 
333:    These matrices should not be destroyed.

335: .keywords: Application, hessian

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


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

363:    Output Parameter:
364: +  DXL - the lower bounds
365: -  DXU - the upper bounds

367:    Level: intermediate

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

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

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

381: .keywords: Application, bounds, duals

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



402: /*@C
403:    TaoSetInequality Constraints - Set inequality constraints for OOQP

405:    Collective on TAO_SOLVER

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

413:    Level: intermediate

415: .keywords: TAO_SOLVER, inequalities

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

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