Actual source code: taoapp_petsc.c

  1: #include "taoapp_petsc.h"  /*I  "tao.h"  I*/
  2: #include "taoapp.h"
  3: #include "src/petsctao/include/taopetsc.h"


 10: TaoPetscApplication::TaoPetscApplication(MPI_Comm mpicomm){
 11:    this->comm=mpicomm;
 12:    this->Setup();
 13:    return;
 14: }

 18: int TaoPetscApplication::Setup(){


 22:   this->taox=0;this->taoh=0;
 23:   this->taofv=0; this->taoj=0;
 24:   this->taofvll=0; this->taofvuu=0; this->taoAA=0;

 26:   this->tao=0;
 27:   this->papp=0;

 29:   this->ksptmp = 0;
 30:   return(0);
 31: }


 36: int TaoPetscApplication::TakeDown()
 37: {
 38:   int info;

 41:   info = TaoVecDestroy(taox); CHKERRQ(info);
 42:   info = TaoMatDestroy(taoh); CHKERRQ(info);
 43:   info = TaoVecDestroy(taofv); CHKERRQ(info);
 44:   info = TaoMatDestroy(taoj); CHKERRQ(info);
 45:   info = TaoVecDestroy(taofvll); CHKERRQ(info);
 46:   info = TaoVecDestroy(taofvuu); CHKERRQ(info);
 47:   info = TaoMatDestroy(taoAA); CHKERRQ(info);
 48:   if (ksptmp) {
 49:     delete ksptmp;
 50:     ksptmp=0;
 51:   }
 52:   return(0);
 53: }

 57: TaoPetscApplication::~TaoPetscApplication()
 58: {
 59:   this->TakeDown();
 60:   return;
 61: }

 65: int TaoPetscApplication::GetComm(MPI_Comm *c) {
 67:   *c = comm;
 68:   return(0);
 69: }

 73: int TaoPetscApplication::GetVariableVector(TaoVec **xx){
 74:   int info;
 75:   Vec V;
 77:   info=TaoAppGetSolutionVec(this->papp,&V); CHKERRQ(info);
 78:   if (this->taox==0){
 79:     info = TaoWrapPetscVec(V,&this->taox); CHKERRQ(info);
 80:   }
 81:   info = this->taox->SetVec(V); CHKERRQ(info);
 82:   *xx= this->taox;
 83:   return(0);
 84: }



 90: int TaoPetscApplication::EvaluateObjectiveFunction(TaoVec *tx, double *ff)
 91: {
 92:   TaoVecPetsc *px = dynamic_cast <TaoVecPetsc *> (tx);

 94:   int     info;
 96:   info=TaoAppComputeObjective(papp,px->GetVec(),ff);CHKERRQ(info);
 97:   return(0);
 98: }


103: int TaoPetscApplication::EvaluateGradient(TaoVec *tx, TaoVec *tg)
104: {
105:   TaoVecPetsc *px = dynamic_cast <TaoVecPetsc *> (tx);
106:   TaoVecPetsc *pg = dynamic_cast <TaoVecPetsc *> (tg);
107:   int     info;
109:   info=TaoAppComputeGradient(papp,px->GetVec(),pg->GetVec());CHKERRQ(info);
110:   return(0);
111: }



117: int TaoPetscApplication::EvaluateObjectiveAndGradient(TaoVec *tx, double *ff, 
118:                                                       TaoVec *tg)
119: {
120:   TaoVecPetsc *px = dynamic_cast <TaoVecPetsc *> (tx);
121:   TaoVecPetsc *pg = dynamic_cast <TaoVecPetsc *> (tg);
122:   int     info;
124:   info=TaoAppComputeObjectiveAndGradient(papp,px->GetVec(),ff,pg->GetVec()); CHKERRQ(info);
125:   return(0);
126: }



132: int TaoPetscApplication::GetHessianMatrix(TaoMat **HH){
133:   int info;
134:   Mat H,HP;
135:   MatStructure flag;

138:   info=TaoAppGetHessianMat(papp,&H,&HP); CHKERRQ(info);
139:   if (this->taoh==0){
140:     info = TaoWrapPetscMat(H,&this->taoh); CHKERRQ(info);
141:   }
142:   info=this->taoh->GetMatrix(0,0,&flag);CHKERRQ(info);
143:   info=this->taoh->SetMatrix(H,HP,flag);CHKERRQ(info);
144:   *HH=this->taoh;
145:   return(0);
146: }


151: int TaoPetscApplication::EvaluateHessian(TaoVec *tx, TaoMat *th)
152: {
153:   TaoVecPetsc *px = dynamic_cast <TaoVecPetsc *> (tx);
154:   TaoMatPetsc *ph = dynamic_cast <TaoMatPetsc *> (th);

156:   int     info;
157:   Mat H, HPre;
158:   MatStructure pflag;

161:   //  info = TaoAppGetGradientVec(this,&this->FDGrad);CHKERRQ(info);
162:   info = ph->GetMatrix(&H,&HPre,&pflag); CHKERRQ(info);
163:   info = TaoAppComputeHessian(this->papp, px->GetVec(), &H, &HPre, &pflag); CHKERRQ(info);
164:   info = ph->SetMatrix(H,HPre,pflag); CHKERRQ(info);
165:   return(0);
166: }

170: int TaoPetscApplication::HessianSolve(TaoVec *tv, TaoVec *tw, TaoTruth *success)
171: {
172:   TaoVecPetsc *pv = dynamic_cast <TaoVecPetsc *> (tv);
173:   TaoVecPetsc *pw = dynamic_cast <TaoVecPetsc *> (tw);

175:   int     info;
176:   PetscTruth success2;
178:   info=TaoAppHessianSolve(papp,pv->GetVec(),pw->GetVec(),&success2);CHKERRQ(info);
179:   if (success2==PETSC_TRUE){*success=TAO_TRUE;} else { *success=TAO_FALSE;}
180:   return(0);
181: }

185: int TaoPetscApplication::EvaluateVariableBounds(TaoVec *tl, TaoVec *tu)
186: {
187:   TaoVecPetsc *pl = dynamic_cast <TaoVecPetsc *> (tl);
188:   TaoVecPetsc *pu = dynamic_cast <TaoVecPetsc *> (tu);
189:   int info;
191:   info=TaoAppComputeVariableBounds(this->papp,pl->GetVec(),pu->GetVec()); CHKERRQ(info);
192:   return(0);
193: }



199: int TaoPetscApplication::GetConstraintVector(TaoVec **rr){
200:   int info;
201:   Vec R;
203:   info=TaoAppGetFunctionVec(this->papp,&R); CHKERRQ(info);
204:   if (this->taofv==0){
205:     info = TaoWrapPetscVec(R,&this->taofv); CHKERRQ(info);
206:   }
207:   info = this->taofv->SetVec(R); CHKERRQ(info);
208:   *rr= this->taofv;
209:   return(0);
210: }


215: int TaoPetscApplication::EvaluateConstraints(TaoVec *tx, TaoVec *tr)
216: {
217:   TaoVecPetsc *px = dynamic_cast <TaoVecPetsc *> (tx);
218:   TaoVecPetsc *pr = dynamic_cast <TaoVecPetsc *> (tr);
219:   int     info;

222:   info=TaoAppComputeFunction(papp,px->GetVec(),pr->GetVec());CHKERRQ(info);
223:   return(0);
224: }


229: int TaoPetscApplication::GetJacobianMatrix(TaoMat **JJ){
230:   int info;
231:   Mat J,JP;
232:   MatStructure flag;

235:   info=TaoAppGetJacobianMat(this->papp,&J,&JP); CHKERRQ(info);
236:   if (this->taoj==0){
237:     info = TaoWrapPetscMat(J,&this->taoj); CHKERRQ(info);
238:   }
239:   info=this->taoj->GetMatrix(0,0,&flag);CHKERRQ(info);
240:   info=this->taoj->SetMatrix(J,JP,flag);CHKERRQ(info);
241:   *JJ=this->taoj;
242:   return(0);
243: }

247: int TaoPetscApplication::EvaluateJacobian(TaoVec *tx, TaoMat *tj)
248: {
249:   TaoVecPetsc *px = dynamic_cast <TaoVecPetsc *> (tx);
250:   TaoMatPetsc *pj = dynamic_cast <TaoMatPetsc *> (tj);

252:   int     info;
253:   Mat J,JPre;
254:   MatStructure pflag;

257:   info = pj->GetMatrix(&J,&JPre,&pflag);CHKERRQ(info);
258:   info = TaoAppComputeJacobian(papp,px->GetVec(),&J,&JPre,&pflag);CHKERRQ(info);
259:   info = pj->SetMatrix(J,JPre,pflag);CHKERRQ(info);
260:   return(0);
261: }

265: int TaoPetscApplication::GetInequalityConstraints(TaoVec **LL, TaoMat **JJ, TaoVec**UU){
267:   *LL=this->taofvll;
268:   *JJ=this->taoAA;
269:   *UU=this->taofvuu;
270:   return(0);
271: }

275: int TaoPetscApplication::GetLinearSolver(TaoMat *tH, int stype, TaoLinearSolver **SS)
276: {
277:   TaoMatPetsc *pH = dynamic_cast <TaoMatPetsc *> (tH);

279:   Mat pm, pmpre;
280:   MatStructure pflag;

282:   MPI_Comm comm2;
283:   int info;
284:   MatType mtype;
285:   PC pc;
286:   KSP newksp;
287:   TaoLinearSolverPetsc *S;
288:   PetscTruth flg1,flg2;


292:   info = pH->GetMatrix(&pm,&pmpre,&pflag); CHKERRQ(info);

294:   info = TaoAppGetKSP(papp,&newksp); CHKERRQ(info);
295:   if (this->ksptmp==0){
296:     info = TaoWrapKSP( newksp, &S ); CHKERRQ(info);
297:     this->ksptmp=S;
298:   }

300:   info=this->ksptmp->SetKSP(newksp);CHKERRQ(info);
301:   *SS=this->ksptmp;

303:   if (pm) {
304:     info = PetscObjectGetComm((PetscObject)pm,&comm2); CHKERRQ(info);
305:     info = MatGetType(pm,&mtype); CHKERRQ(info);
306:     info = PetscStrncmp(mtype,MATSEQDENSE,10,&flg1); CHKERRQ(info);
307:     info = PetscStrncmp(mtype,MATMPIDENSE,10,&flg2); CHKERRQ(info);
308:   } 
309:   else {
310:     comm2 = this->comm;
311:   }

313:   info = KSPGetPC(newksp,&pc); CHKERRQ(info);
314:   info=KSPSetFromOptions(newksp); CHKERRQ(info);
315:   if (stype==220){
316:     info=PCSetFromOptions(pc); CHKERRQ(info);
317:     info=PCSetType(pc,PCJACOBI); CHKERRQ(info);
318:     info=KSPSetType(newksp,KSPSTCG); CHKERRQ(info);
319:   }

321:   return(0);
322: }

326: int TaoPetscApplication::InitializeVariables(TaoVec *xx){
328:   return(0);
329: }


334: int TaoPetscApplication::Monitor(){
335:   int info;
337:   info=TaoAppMonitor(this->papp); CHKERRQ(info);
338:   return(0);
339: }

343: int TaoPetscApplication::Monitor2(TaoVec *xx, TaoVec *tg, TaoVec *dx, TaoTruth *flag){
344:   int     info;
345:   PetscTruth flag2;

348:   if (tg) {
349:     TaoVecPetsc *pg = dynamic_cast <TaoVecPetsc *> (tg);
350:     *flag=TAO_FALSE;
351:     info=TaoAppCheckConvergence(papp, pg->GetVec(), &flag2); CHKERRQ(info);
352:     if (flag2==PETSC_TRUE) *flag=TAO_TRUE;
353:   } 
354:   return(0);
355: }