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->ksptmp=0;
 27:   this->tao=0;
 28:   this->papp=0;

 30:   return(0);
 31: }


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

 40:   info = TaoVecDestroy(this->taox); CHKERRQ(info);
 41:   info = TaoMatDestroy(this->taoh); CHKERRQ(info);
 42:   info = TaoVecDestroy(this->taofv); CHKERRQ(info);
 43:   info = TaoMatDestroy(this->taoj); CHKERRQ(info);
 44:   info = TaoVecDestroy(this->taofvll); CHKERRQ(info);
 45:   info = TaoVecDestroy(this->taofvuu); CHKERRQ(info);
 46:   info = TaoMatDestroy(this->taoAA); CHKERRQ(info);
 47:   if (this->ksptmp){ info=KSPDestroy(this->ksptmp); CHKERRQ(info);}
 48:   return(0);
 49: }

 53: TaoPetscApplication::~TaoPetscApplication(){
 54:   this->TakeDown();
 55:   return;
 56: }

 60: int TaoPetscApplication::GetVariableVector(TaoVec **xx){
 61:   int info;
 62:   Vec V;
 64:   info=TaoAppGetSolutionVec(this->papp,&V); CHKERRQ(info);
 65:   if (this->taox==0){
 66:     info = TaoWrapPetscVec(V,&this->taox); CHKERRQ(info);
 67:   }
 68:   info = this->taox->SetVec(V); CHKERRQ(info);
 69:   *xx= this->taox;
 70:   return(0);
 71: }



 77: int TaoPetscApplication::EvaluateObjectiveFunction(TaoVec *xx, double *ff){
 78:   int     info;
 79:   Vec X = (Vec) xx->VecObject;
 81:   info=TaoAppComputeObjective(this->papp,X,ff);CHKERRQ(info);
 82:   return(0);
 83: }


 88: int TaoPetscApplication::EvaluateGradient(TaoVec *xx, TaoVec *gg){
 89:   int     info;
 90:   Vec X = (Vec) xx->VecObject, G = (Vec) gg->VecObject;
 92:   info=TaoAppComputeGradient(this->papp,X,G);CHKERRQ(info);
 93:   return(0);
 94: }



100: int TaoPetscApplication::EvaluateObjectiveAndGradient(TaoVec *xx, double *ff, TaoVec *gg){
101:   int     info;
102:   Vec X = (Vec) xx->VecObject, G = (Vec) gg->VecObject;
104:   info=TaoAppComputeObjectiveAndGradient(this->papp,X,ff,G); CHKERRQ(info);
105:   return(0);
106: }



112: int TaoPetscApplication::GetHessianMatrix(TaoMat **HH){
113:   int info;
114:   Mat H,HP;
115:   MatStructure flag;

118:   info=TaoAppGetHessianMat(this->papp,&H,&HP); CHKERRQ(info);
119:   if (this->taoh==0){
120:     info = TaoWrapPetscMat(H,&this->taoh); CHKERRQ(info);
121:   }
122:   info=this->taoh->GetMatrix(0,0,&flag);CHKERRQ(info);
123:   info=this->taoh->SetMatrix(H,HP,flag);CHKERRQ(info);
124:   *HH=this->taoh;
125:   return(0);
126: }


131: int TaoPetscApplication::EvaluateHessian(TaoVec *xx, TaoMat *HH){
132:   int     info;
133:   Vec X = (Vec) xx->VecObject;
134:   TaoMatPetsc *MM=(TaoMatPetsc*)HH;
135:   Mat H, HPre;
136:   MatStructure pflag;

139:   //  info = TaoAppGetGradientVec(this,&this->FDGrad);CHKERRQ(info);
140:   info = MM->GetMatrix(&H,&HPre,&pflag);CHKERRQ(info);
141:   info = TaoAppComputeHessian(this->papp, X, &H, &HPre, &pflag);
142:   info = MM->SetMatrix(H,HPre,pflag);CHKERRQ(info);
143:   return(0);
144: }


149: int TaoPetscApplication::HessianSolve(TaoVec *vv, TaoVec *ww, TaoTruth *success){
150:   int     info;
151:   Vec V = (Vec) vv->VecObject;
152:   Vec W = (Vec) ww->VecObject;
153:   PetscTruth success2;
155:   info=TaoAppHessianSolve(this->papp,V,W,&success2);CHKERRQ(info);
156:   if (success2==PETSC_TRUE){*success=TAO_TRUE;} else { *success=TAO_FALSE;}
157:   return(0);
158: }

162: int TaoPetscApplication::EvaluateVariableBounds(TaoVec *xxll, TaoVec *xxuu){
163:   int info;
164:   Vec XL=(Vec)xxll->VecObject,XU=(Vec)xxuu->VecObject;
166:   info=TaoAppComputeVariableBounds(this->papp,XL,XU); CHKERRQ(info);
167:   return(0);
168: }



174: int TaoPetscApplication::GetConstraintVector(TaoVec **rr){
175:   int info;
176:   Vec R;
178:   info=TaoAppGetFunctionVec(this->papp,&R); CHKERRQ(info);
179:   if (this->taofv==0){
180:     info = TaoWrapPetscVec(R,&this->taofv); CHKERRQ(info);
181:   }
182:   info = this->taofv->SetVec(R); CHKERRQ(info);
183:   *rr= this->taofv;
184:   return(0);
185: }


190: int TaoPetscApplication::EvaluateConstraints(TaoVec *xx, TaoVec *RR){
191:   int     info;
192:   Vec X = (Vec) xx->VecObject, R = (Vec) RR->VecObject;

195:   info=TaoAppComputeFunction(this->papp,X,R);CHKERRQ(info);
196:   return(0);
197: }


202: int TaoPetscApplication::GetJacobianMatrix(TaoMat **JJ){
203:   int info;
204:   Mat J,JP;
205:   MatStructure flag;

208:   info=TaoAppGetJacobianMat(this->papp,&J,&JP); CHKERRQ(info);
209:   if (this->taoj==0){
210:     info = TaoWrapPetscMat(J,&this->taoj); CHKERRQ(info);
211:   }
212:   info=this->taoj->GetMatrix(0,0,&flag);CHKERRQ(info);
213:   info=this->taoj->SetMatrix(J,JP,flag);CHKERRQ(info);
214:   *JJ=this->taoj;
215:   return(0);
216: }

220: int TaoPetscApplication::EvaluateJacobian(TaoVec *xx, TaoMat *JJ){

222:   int     info;
223:   Vec X = (Vec) xx->VecObject;
224:   Mat J,JPre;
225:   MatStructure pflag;
226:   TaoMatPetsc *JJJ=(TaoMatPetsc*)JJ;

229:   info = JJJ->GetMatrix(&J,&JPre,&pflag);CHKERRQ(info);
230:   info = TaoAppComputeJacobian(this->papp,X,&J,&JPre,&pflag);CHKERRQ(info);
231:   info = JJJ->SetMatrix(J,JPre,pflag);CHKERRQ(info);
232:   return(0);
233: }

237: int TaoPetscApplication::GetInequalityConstraints(TaoVec **LL, TaoMat **JJ, TaoVec**UU){
239:   *LL=this->taofvll;
240:   *JJ=this->taoAA;
241:   *UU=this->taofvuu;
242:   return(0);
243: }

245: int MatCreateADA(Mat,Vec,Vec,Mat*);

249: int TaoPetscApplication::CreateATDAMatrix(TaoMat *A, TaoVec* D,TaoVec*X, TaoMat**AADAA){
250:   int info;
251:   Vec D1 = (Vec) D->VecObject;
252:   Vec D2 = (Vec) X->VecObject;
253:   Mat ADA;
254:   Mat pm=(Mat)A->MatObject;
255:   TaoMatPetsc *TT;
257:   info = MatCreateADA(pm,D1,D2, &ADA); CHKERRQ(info);
258:   info = TaoWrapPetscMat(ADA, &TT); CHKERRQ(info);
259:   info = PetscObjectDereference((PetscObject)ADA);
260:   *AADAA=TT;
261:   return(0);
262: }



268: int TaoPetscApplication::GetLinearSolver(TaoMat *H, int stype, TaoLinearSolver **SS){
269:   MPI_Comm comm2;
270:   int info;
271:   MatType mtype;
272:   PC pc;
273:   KSP newksp;
274:   TaoLinearSolver *S;
275:   PetscTruth flg1,flg2;
276:   Mat pm=(Mat)H->MatObject;

279:   //  info = ((TaoMatPetsc*)H)->CreateLinearSolver(flag,tksp); CHKERRQ(info);
280:   //  info = TaoGetKSP(this->papp,&newksp); CHKERRQ(info);
281:   newksp=0;
282:   if (this->ksptmp){
283:     info = TaoWrapKSP( newksp, (TaoLinearSolverPetsc**)&S ); CHKERRQ(info);
284:     PetscObjectDereference((PetscObject)(ksptmp)); CHKERRQ(info);
285:     *SS=S;
286:     this->ksptmp=0;
287:     return(0);
288:   }

290:   if (pm){
291:     info = PetscObjectGetComm((PetscObject)pm,&comm2); CHKERRQ(info);
292:     info = MatGetType(pm,&mtype); CHKERRQ(info);
293:     info = PetscStrncmp(mtype,MATSEQDENSE,10,&flg1); CHKERRQ(info);
294:     info = PetscStrncmp(mtype,MATMPIDENSE,10,&flg2); CHKERRQ(info);
295:   } else {
296:     comm2 = this->comm;
297:   }
298:   info = KSPCreate(comm,&newksp); CHKERRQ(info);
299:   info = KSPGetPC(newksp,&pc); CHKERRQ(info);
300:   if (flg1==PETSC_TRUE || flg2==PETSC_TRUE){
301:     info=PCSetType(pc,PCJACOBI); CHKERRQ(info);
302:   } else {
303:     info=PCSetType(pc,PCBJACOBI); CHKERRQ(info);
304:   }

306:   if (stype==300){
307:     info=KSPSetType(newksp,KSPCG); CHKERRQ(info);
308:     info=KSPSetFromOptions(newksp); CHKERRQ(info);
309:   } else if (stype==310){
310:     info=KSPSetFromOptions(newksp); CHKERRQ(info);
311:     info=KSPSetType(newksp,KSPGMRES); CHKERRQ(info);
312:   } else if (stype==220){
313:     info=KSPGetPC(newksp,&pc); CHKERRQ(info);
314:     info=PCSetType(pc,PCNONE); CHKERRQ(info);
315:     info=PCSetFromOptions(pc); CHKERRQ(info);
316:     info=KSPSetType(newksp,KSPQCG); CHKERRQ(info);
317:   } else if (stype==110){
318:     info=KSPSetFromOptions(newksp); CHKERRQ(info);
319:     info=KSPSetType(newksp,KSPGMRES); CHKERRQ(info);
320:   } else {
321:     info=KSPSetType(newksp,KSPGMRES); CHKERRQ(info);
322:     info=KSPSetFromOptions(newksp); CHKERRQ(info);
323:   }

325:   info = TaoWrapKSP( newksp, (TaoLinearSolverPetsc**)&S ); CHKERRQ(info);
326:   //  info=KSPDestroy(newksp);CHKERRQ(info); newksp=0;
327:   PetscObjectDereference((PetscObject)(newksp)); CHKERRQ(info);
328:   *SS=S;
329:   //  info = TaoAppSetKSP(this->papp,newksp); CHKERRQ(info);

331:   return(0);
332: }



338: int TaoPetscApplication::InitializeVariables(TaoVec *xx){
340:   return(0);
341: }


346: int TaoPetscApplication::Monitor(){
347:   int info;
349:   info=TaoAppMonitor(this->papp); CHKERRQ(info);
350:   return(0);
351: }

355: int TaoPetscApplication::Monitor2(TaoVec *xx, TaoVec *gl, TaoVec *dx, TaoTruth *flag){
356:   int     info;
357:   PetscTruth flag2;
358:   Vec GL = (Vec) gl->VecObject;
360:   *flag=TAO_FALSE;
361:   info=TaoAppCheckConvergence(this->papp, GL,&flag2); CHKERRQ(info);
362:   if (flag2==PETSC_TRUE) *flag=TAO_TRUE;
363:   return(0);
364: }