Actual source code: taoappobject.c

  1: #include "taoappobject.h"     /*I "tao_solver.h"  I*/
  2: #include "tao_general.h"
  3: #include "taovec.h"


  8: /*@C
  9:    EvaluateObjectiveFunction - Evaluate the objective function
 10:    at the point x.

 12:    Collective on TAO_SOLVER

 14:    Input Parameters:
 15: .  xx - the point at which to evaluate the obective function

 17:    Output Parameters:
 18: .  ff - objective function value

 20: .seealso TaoSetPetscGradient(), EvaluateObjectiveFunction(), EvaluateObjectiveAndGradientFunction()

 22:    Level: intermediate

 24: .keywords: application, context

 26: @*/
 27: int TaoApplication::EvaluateObjectiveFunction(TaoVec *xx, double *ff){
 28:   TaoFunctionBegin;
 29:   SETERRQ(56,"Operation not defined");
 30:   //  TaoFunctionReturn(0);
 31: }

 35: /*@C
 36:    EvaluateGradient - Evaluate the gradient of the
 37:    objective function at the point x.

 39:    Collective on TAO_SOLVER

 41:    Input Parameters:
 42: .  xx - the point at which to evaluate the gradient

 44:    Output Parameters:
 45: .  gg - the gradient vector

 47: .seealso TaoSetPetscGradient(), EvaluateObjectiveFunction(), EvaluateObjectiveAndGradientFunction()

 49:    Level: intermediate

 51: .keywords: application, context

 53: @*/
 54: int TaoApplication::EvaluateGradient(TaoVec *xx, TaoVec *gg){
 55:   TaoFunctionBegin;
 56:   SETERRQ(56,"Operation not defined");
 57:   //  TaoFunctionReturn(0);
 58: }


 63: /*@C
 64:    EvaluateObjectiveAndGradient - Evaluate the objective function
 65:    and its gradient at the point x.

 67:    Collective on TAO_SOLVER

 69:    Input Parameters:
 70: .  xx - the point at which to evaluate the gradient

 72:    Output Parameters:
 73: +  ff - the objective value
 74: -  gg - the gradient vector

 76: .seealso TaoSetPetscFunctionGradient(), EvaluateGradient(), EvaluateObjectiveFunction()

 78:    Level: intermediate

 80: .keywords: application, gradient

 82: @*/
 83: int TaoApplication::EvaluateObjectiveAndGradient(TaoVec *xx, double *ff, TaoVec *gg){
 84:   int info;
 85:   TaoFunctionBegin;
 86:   info=this->EvaluateObjectiveFunction(xx,ff);CHKERRQ(info);
 87:   info=this->EvaluateGradient(xx,gg);CHKERRQ(info);
 88:   TaoFunctionReturn(0);
 89: }

 93: /*@C
 94:    EvaluateHessian - Evaluate the Hessian of the objective function
 95:    at the point x.

 97:    Collective on TAO_SOLVER

 99:    Input Parameters:
100: .  xx - the point at which to evaluate the gradient

102:    Output Parameters:
103: .  HH - the Hessian matrix

105: .seealso TaoSetPetscHessian()
106:    Level: intermediate

108: .keywords: application, hessian

110: @*/
111: int TaoApplication::EvaluateHessian(TaoVec *xx, TaoMat *HH){
112:   TaoFunctionBegin;
113:   SETERRQ(56,"Operation not defined");
114:   //  TaoFunctionReturn(0);
115: }

119: /*@C
120:    HessianSolve - Solve a linear system involving the Hessian
121: matrix, or apply an inverse Hessian operator to a vector.

123:    Collective on TAO_SOLVER

125:    Input Parameters:
126: .  rhs - the right-hand side

128:    Output Parameters:
129: +  vv - the solution
130: -  success - flag states whether solution was found.

132: .seealso TaoSetMethod(), TaoLMVMSetSize(), EvaluateHessian()

134:    Level: intermediate

136: .keywords: application, Hessian, LMVM

138: @*/
139: int TaoApplication::HessianSolve(TaoVec *rhs, TaoVec*vv, TaoTruth *success){
140:   int info;
141:   TaoFunctionBegin;
142:   info=vv->CopyFrom(rhs);CHKERRQ(info);
143:   *success=TAO_TRUE;
144:   TaoFunctionReturn(0);
145: }


150: /*@C
151:    EvaluateConstraints - Evaluate the constraint functions
152:    at the point x.

154:    Collective on TAO_SOLVER

156:    Input Parameters:
157: .  xx - the point at which to evaluate the gradient

159:    Output Parameters:
160: .  RR - the constraint vector

162: .seealso TaoSetPetscConstraintsFunction()
163:    Level: intermediate

165: .keywords: application, constraints

167: @*/
168: int TaoApplication::EvaluateConstraints(TaoVec *xx, TaoVec *RR){
169:   TaoFunctionBegin;
170:   SETERRQ(56,"Operation not defined");
171:   //  TaoFunctionReturn(0);
172: }

176: /*@C
177:    EvaluateJacobian - Evaluate the Jacobian of the constraint functions
178:    at the point x.

180:    Collective on TAO_SOLVER

182:    Input Parameters:
183: .  xx - the point at which to evaluate the gradient

185:    Output Parameters:
186: .  JJ - the Jacobian matrix

188: .seealso TaoSetPetscJacobian()

190:    Level: intermediate

192: .keywords: application, constraints

194: @*/
195: int TaoApplication::EvaluateJacobian(TaoVec *xx, TaoMat *JJ){
196:   TaoFunctionBegin;
197:   SETERRQ(56,"Operation not defined");
198:   //  TaoFunctionReturn(0);
199: }

203: /*@C
204:    InitializeVariables - Initialize the variables for the optimization
205:    solver.  Set an initial point.


208:    Collective on TAO_SOLVER

210:    Input Parameters:
211: .  xx - the variable vector

213:    Level: intermediate

215: .keywords: application, constraints

217: @*/
218: int TaoApplication::InitializeVariables(TaoVec *xx){
219:   int info;
220:   TaoFunctionBegin;
221:   info=xx->SetToZero();CHKERRQ(info);
222:   TaoFunctionReturn(0);
223: }

227: /*@C
228:    GetVariableVector - Sets the pointer to a vector that will
229:    be used to store the variables.

231:    Collective on TAO_SOLVER

233:    Output Parameters:
234: .  xx - vector to the variable vector

236:    Level: intermediate

238: .keywords: application, variables

240: @*/
241: int TaoApplication::GetVariableVector(TaoVec **xx){
242:   TaoFunctionBegin;
243:   *xx=0;
244:   TaoFunctionReturn(0);
245: }


250: /*@C
251:    EvaluateVariableBounds - Set these vectors equal to the lower
252:    and upper bounds on the variables.

254:    Collective on TAO_SOLVER

256:    Input Parameters:
257: +  xxll - vector vector of lower bounds
258: -  xxuu - vector vector of upper bounds

260:    Level: intermediate

262: .keywords: application, bounds

264: @*/
265: int TaoApplication::EvaluateVariableBounds(TaoVec *xxll, TaoVec *xxuu)
266: {
267:   int info;

269:   TaoFunctionBegin;

271:   if (xxll) {
272:     info = xxll->SetToConstant(-TAO_INFINITY); CHKERRQ(info);
273:   }

275:   if (xxuu) {
276:     info = xxuu->SetToConstant(TAO_INFINITY); CHKERRQ(info);
277:   }
278:   TaoFunctionReturn(0);
279: }


284: /*@C
285:    GetHessianMatrix - Sets the pointer to a matrix that will
286:    be used to store the Hessian of the objective function.

288:    Collective on TAO_SOLVER

290:    Output Parameters:
291: .  HH - vector to the Hessian

293:    Level: intermediate

295: .keywords: application, gradient

297: @*/
298: int TaoApplication::GetHessianMatrix(TaoMat **HH){
299:   TaoFunctionBegin;
300:   *HH=0;
301:   TaoFunctionReturn(0);
302: }

306: /*@C
307:    GetJacobianMatrix - Sets the pointer to a matrix that will
308:    be used to store the Jacobian

310:    Collective on TAO_SOLVER

312:    Output Parameters:
313: .  JJ - matrix to the Jacobian

315:    Level: intermediate

317: .keywords: application, gradient

319: @*/
320: int TaoApplication::GetJacobianMatrix(TaoMat **JJ){
321:   TaoFunctionBegin;
322:   *JJ=0;
323:   TaoFunctionReturn(0);
324: }

326: int TaoApplication::GetConstraintVector(TaoVec **rr){
327:   TaoFunctionBegin;
328:   *rr=0;
329:   TaoFunctionReturn(0);
330: }


335: int TaoApplication::GetInequalityConstraints(TaoVec**ll, TaoMat **AA, TaoVec **uu){
336:   TaoFunctionBegin;
337:   *ll=0;
338:   *AA=0;
339:   *uu=0;
340:   TaoFunctionReturn(0);
341: }

345: /*@C
346:    Monitor - Monitor the current solution

348:    Collective on TAO_SOLVER

350:    Level: intermediate

352:    Note: This routine is called after each iteration

354: .keywords: application, monitor

356: @*/
357: int TaoApplication::Monitor(){
358:   TaoFunctionBegin;
359:   TaoFunctionReturn(0);
360: }

364: /*@C
365:    Monitor2 - 

367:    Collective on TAO_SOLVER

369:    Input Parameters:
370: .  xx - the current solution
371: .  gl - the gradient of the Lagrangian function
372: .  dx - the step direction

374:    Output Parameters:
375: .  *term - TAO_TRUE if the solver should stop iterating and TAO_FALSE if the solver should continue;

377: .seealso TaoSetPetscGradient(),

379:    Level: intermediate

381: .keywords: application, context

383: @*/
384: int TaoApplication::Monitor2(TaoVec *xx, TaoVec *gg, TaoVec *ff, TaoTruth *term){
385:   TaoFunctionBegin;
386:   *term=TAO_FALSE;
387:   TaoFunctionReturn(0);
388: }



394: /*@C
395:    GetLinearSolver - Create and a linear solver used to solve this matrix

397:    Collective on TAO_SOLVER

399:    input Parameters:
400: +  H - the operator to be used to solve
401: -  flag - indicates properties solver must have

403:    Output Parameters:
404: .  tksp - the linear solver

406:    Types of linear solvers are:
407: +  100 - Solver for any square matrix
408: .  110 - GMRES
409: .  200 - Any symmetric matrix
410: .  210 - MINRES
411: .  220 - CG with Trust  Region
412: .  300 - Any symmetric positive definite
413: -  310 - Conjugate Gradient

415:    Level: intermediate

417: .keywords: application, linear solver

419: @*/
420: int TaoApplication::GetLinearSolver(TaoMat *H, int flag, TaoLinearSolver **tksp){
421:   TaoFunctionBegin;
422:   SETERRQ(56,"Operation not defined");
423: }

427: /*@
428:    TaoDestroyApplication - Destroys the TAO application

430:    Collective on TAO_APPLICATION

432:    Input Parameters:
433: .  myapp - pointer to the TaoApplication object (TaoApplication*),


436:    Level: developer

438: .keywords: application

440: @*/
441: int TaoDestroyApplication(TaoApplication *myapp){
442:   TaoFunctionBegin;
443:   delete myapp;
444:   TaoFunctionReturn(0);
445: }