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: }