Actual source code: oldroutines.c
1: #include "tao.h" /*I "tao.h" I*/
4: class TaoPetscApplication;
5: int TaoAppGetTaoPetscApp(TAO_APPLICATION, TaoPetscApplication**);
7: typedef struct {
8: /* Function Gradient Evaluation over single element */
9: int (*computef)(TAO_SOLVER, Vec, double*, void*);
10: int (*computeg)(TAO_SOLVER, Vec, Vec, void*);
11: int (*computefg)(TAO_SOLVER, Vec, double*, Vec, void*);
12: int (*computeh)(TAO_SOLVER, Vec, Mat*, Mat*, MatStructure*, void*);
13: int (*computer)(TAO_SOLVER, Vec, Vec, void*);
14: int (*computej)(TAO_SOLVER, Vec, Mat*, void*);
15: TAO_SOLVER tao;
16: void *fctx;
17: void *fgctx;
18: void *gctx;
19: void *hctx;
20: void *rctx;
21: void *jctx;
22: } TaoPETScOLD;
27: static int TaoPetscAppF(TAO_APPLICATION taoapp , Vec X , double *f, void*ctx){
28: int info;
29: TaoPETScOLD* mctx= (TaoPETScOLD*)ctx;
31: info=(*mctx->computef)(mctx->tao,X,f,mctx->fctx); CHKERRQ(info);
32: return(0);
33: }
37: static int TaoPetscAppG(TAO_APPLICATION taoapp , Vec X , Vec G, void*ctx){
38: int info;
39: TaoPETScOLD* mctx= (TaoPETScOLD*)ctx;
41: info=(*mctx->computeg)(mctx->tao,X,G,mctx->gctx); CHKERRQ(info);
42: return(0);
43: }
47: static int TaoPetscAppFG(TAO_APPLICATION taoapp , Vec X , double *f, Vec G, void* ctx){
48: int info;
49: TaoPETScOLD* mctx= (TaoPETScOLD*)ctx;
51: info=(*mctx->computefg)(mctx->tao,X,f,G,mctx->fgctx); CHKERRQ(info);
52: return(0);
53: }
57: static int TaoPetscAppH(TAO_APPLICATION taoapp , Vec X , Mat *H, Mat *HP, MatStructure *flag, void* ctx){
58: int info;
59: TaoPETScOLD* mctx= (TaoPETScOLD*)ctx;
61: info=(*mctx->computeh)(mctx->tao,X,H,HP,flag,mctx->hctx); CHKERRQ(info);
62: return(0);
63: }
67: static int TaoPetscAppR(TAO_APPLICATION taoapp , Vec X , Vec R, void*ctx){
68: int info;
69: TaoPETScOLD* mctx= (TaoPETScOLD*)ctx;
71: info=(*mctx->computer)(mctx->tao,X,R,mctx->rctx); CHKERRQ(info);
72: return(0);
73: }
77: static int TaoPetscAppJ(TAO_APPLICATION taoapp , Vec X , Mat *J, Mat *JP, MatStructure *flag, void* ctx){
78: int info;
79: TaoPETScOLD* mctx= (TaoPETScOLD*)ctx;
81: info=(*mctx->computej)(mctx->tao,X,J,mctx->jctx); CHKERRQ(info);
82: return(0);
83: }
88: /* @C
89: TaoSetPetscFunction - Sets a routine that evaluates the function at
90: the specified point.
92: Collective on TAO_SOLVER
94: Input Parameters:
95: + taoapp - the TAO_APPLICATION context
96: . x - variable vector that stores the solution
97: . func - function evaluation routine
98: - ctx - [optional] user-defined context for private data for the
99: function and gradient evaluation routine (may be TAO_NULL)
101: Note: This function is no longer supported. Use TaoSetObjectiveFunction() and TaoSetVariableVec() instead.
103: Level: intermediate
105: .seealso: TaoAppSetObjectiveRoutine(), TaoAppSetInitialSolutionVec()
106: @ */
107: int TaoSetPetscFunction(TAO_APPLICATION taoapp, Vec X, int (*func)(TAO_SOLVER,Vec,double*,void*),void *ctx){
108: int info;
109: TaoPETScOLD* fgctx= (TaoPETScOLD*)ctx;
111: PetscNew(TaoPETScOLD,&fgctx);
112: fgctx->computef = func;
113: fgctx->fctx=ctx;
114: fgctx->tao = 0;
115: info = TaoAppSetInitialSolutionVec(taoapp,X); CHKERRQ(info);
116: info = TaoAppSetObjectiveRoutine(taoapp,TaoPetscAppF,(void*)fgctx); CHKERRQ(info);
117: info = TaoAppSetDestroyRoutine(taoapp,TaoApplicationFreeMemory, (void*)fgctx); CHKERRQ(info);
118: return(0);
119: }
124: /* @C
125: TaoSetPetscGradient - Sets the gradient evaluation routine and gradient
126: vector for use by the TAO_SOLVER routines.
128: Collective on TAO_SOLVER
130: Input Parameters:
131: + taoapp - the TAO_APPLICATION context
132: . g - vector to store gradient
133: . grad - gradient evaluation routine
134: - ctx - [optional] user-defined function context
136: Note:
137: This routine is no longer supported.
139: .seealso: TaoAppSetGradientRoutine()
141: @ */
142: int TaoSetPetscGradient(TAO_APPLICATION taoapp, Vec g, int (*grad)(TAO_SOLVER,Vec,Vec,void*),void *ctx){
143: int info;
144: TaoPETScOLD* fgctx= (TaoPETScOLD*)ctx;
146: PetscNew(TaoPETScOLD,&fgctx);
147: fgctx->computeg = grad;
148: fgctx->gctx=ctx;
149: fgctx->tao = 0;
150: info = TaoAppSetGradientRoutine(taoapp,TaoPetscAppG,(void*)fgctx); CHKERRQ(info);
151: info = TaoAppSetDestroyRoutine(taoapp,TaoApplicationFreeMemory, (void*)fgctx); CHKERRQ(info);
152: return(0);
153: }
157: /* @C
158: TaoSetPetscFunctionGradient - Sets a routine for function and gradient evaluation.
160: Collective on TAO_SOLVER
162: Input Parameters:
163: + taoapp - the TAO_APPLICATION context
164: . x - vector to store solution
165: . g - vector to store gradient
166: . funcgrad - routine for evaluating the function and gradient
167: - ctx - optional user-defined context for private data for the
168: function and gradient evaluation routine (may be TAO_NULL)
170: Level: intermediate
172: Note:
173: This routine is no longer supported.
175: .seealso: TaoAppSetObjectiveAndGradientRoutine(), TaoSetVariableVec()
177: @ */
178: int TaoSetPetscFunctionGradient(TAO_APPLICATION taoapp, Vec x, Vec g, int (*funcgrad)(TAO_SOLVER,Vec,double*,Vec, void*),void *ctx){
179: int info;
180: TaoPETScOLD* fgctx= (TaoPETScOLD*)ctx;
182: PetscNew(TaoPETScOLD,&fgctx);
183: fgctx->computefg = funcgrad;
184: fgctx->fgctx=ctx;
185: fgctx->tao = 0;
186: info = TaoAppSetObjectiveAndGradientRoutine(taoapp,TaoPetscAppFG,(void*)fgctx); CHKERRQ(info);
187: info = TaoAppSetDestroyRoutine(taoapp,TaoApplicationFreeMemory, (void*)fgctx); CHKERRQ(info);
188: info = TaoAppSetInitialSolutionVec(taoapp,x); CHKERRQ(info);
189: return(0);
190: }
194: /* @C
195: TaoSetPetscHessian - Sets the function to compute the Hessian as well as the
196: location to store the matrix.
198: Collective on TAO_SOLVER and Mat
200: Input Parameters:
201: + taoapp - the TAO_APPLICATION context
202: . H - Hessian matrix
203: . Hpre - preconditioner matrix (usually same as the Hessian)
204: . hess - Hessian evaluation routine
205: - ctx - [optional] user-defined context for private data for the
206: Hessian evaluation routine (may be TAO_NULL)
208: Level: intermediate
210: Note: This routine is no longer supported.
212: .seealso: TaoAppSetHessianRoutine(), TaoAppSetHessianMat()
213: @ */
214: int TaoSetPetscHessian(TAO_APPLICATION taoapp, Mat H, Mat Hpre, int (*hess)(TAO_SOLVER,Vec,Mat*,Mat*,MatStructure*,void*),void *ctx){
215: int info;
216: TaoPETScOLD* fgctx= (TaoPETScOLD*)ctx;
218: PetscNew(TaoPETScOLD,&fgctx);
219: fgctx->computeh = hess;
220: fgctx->hctx=ctx;
221: fgctx->tao = 0;
222: info = TaoAppSetHessianRoutine(taoapp,TaoPetscAppH,(void*)fgctx);CHKERRQ(info);
223: info = TaoAppSetHessianMat(taoapp,H,Hpre);CHKERRQ(info);
224: info = TaoAppSetDestroyRoutine(taoapp,TaoApplicationFreeMemory, (void*)fgctx); CHKERRQ(info);
225: return(0);
226: }
231: /* @C
232: TaoSetPetscConstraintsFunction - Sets the routine that evaluates
233: equality constraint functions.
235: Collective on TAO_APPLICATION
237: Input Parameters:
238: + taoapp - the TAO_APPLICATION context
239: . r - vector to constrainf function values
240: . func - the constraint function evaluation routine
241: - ctx - [optional] user-defined function context
243: Calling sequence of func:
244: $ func (TAO_SOLVER solver,Vec x,Vec r,void *ctx);
246: + solver - the TAO_SOLVER solver context
247: . x - input vector
248: . r - constraint vector
249: - ctx - user-defined function gradient context set from TaoSetPetscFunction()
251: Level: intermediate
253: .keywords: TAO_SOLVER, set, function
255: .seealso: TaoGetGradient(), TaoSetPetscHessian()
257: @ */
258: int TaoSetPetscConstraintsFunction(TAO_APPLICATION taoapp, Vec r, int (*func)(TAO_SOLVER,Vec,Vec,void*),void
259: *ctx){
260: int info;
261: TaoPETScOLD* fgctx= (TaoPETScOLD*)ctx;
263: PetscNew(TaoPETScOLD,&fgctx);
264: fgctx->computer = func;
265: fgctx->rctx=ctx;
266: fgctx->tao = 0;
267: info = TaoAppSetConstraintRoutine(taoapp,TaoPetscAppR,(void*)fgctx); CHKERRQ(info);
268: info = TaoAppSetDestroyRoutine(taoapp,TaoApplicationFreeMemory, (void*)fgctx); CHKERRQ(info);
269: info = TaoAppSetFunctionVec(taoapp,r); CHKERRQ(info);
270: return(0);
271: }
276: /* @C
277: TaoSetPetscJacobian - Sets the function to compute the Jacobian of
278: the equality constraint function as well as the
279: location to store the matrix.
281: Collective on TAO_APPLICATION and Mat
283: Input Parameters:
284: + solver - the TAO_APPLICATION context
285: . J - Jacobian matrix
286: . jac - Jacobian evaluation routine
287: - ctx - [optional] user-defined context for private data for the
288: Hessian evaluation routine (may be TAO_NULL)
290: Calling sequence of func:
291: $ jac (TAO_SOLVER solver,Vec x,Mat *J,void *ctx);
293: + solver - the TAO_SOLVER solver context
294: . x - input vector
295: . J - Jacobian matrix
296: - ctx - [optional] user-defined Hessian context
298: Notes:
300: The function jac() takes Mat * as the matrix arguments rather than Mat.
301: This allows the Jacobian evaluation routine to replace J with a
302: completely new new matrix structure (not just different matrix elements)
303: when appropriate, for instance, if the nonzero structure is changing
304: throughout the global iterations.
306: Level: intermediate
308: .keywords: TAO_SOLVER, Jacobian
310: .seealso: TaoSetPetscFunction(), TaoSetPetscConstraintsFunction()
311: @ */
312: int TaoSetPetscJacobian(TAO_APPLICATION taoapp,Mat J,int (*jac)(TAO_SOLVER,Vec,Mat*,void*),void *ctx){
313: int info;
314: TaoPETScOLD* fgctx= (TaoPETScOLD*)ctx;
316: PetscNew(TaoPETScOLD,&fgctx);
317: fgctx->computej = jac;
318: fgctx->jctx=ctx;
319: fgctx->tao = 0;
320: info = TaoAppSetJacobianRoutine(taoapp,TaoPetscAppJ,(void*)fgctx);CHKERRQ(info);
321: info = TaoAppSetDestroyRoutine(taoapp,TaoApplicationFreeMemory, (void*)fgctx); CHKERRQ(info);
322: info = TaoAppSetJacobianMat(taoapp,J,J); CHKERRQ(info);
323: return(0);
324: }
331: int TaoSetFromOptions(TAO_APPLICATION taoapp){
332: int info;
334: info=TaoAppSetFromOptions(taoapp); CHKERRQ(info);
335: return(0);
336: }
338: class TaoApplication;
343: int TaoSetApplication(TAO_SOLVER tao, TAO_APPLICATION taoapp){
344: int info;
345: TaoPetscApplication* tpapp;
347: info=TaoAppGetTaoPetscApp(taoapp, &tpapp); CHKERRQ(info);
348: info=TaoSetApplication(tao,(TaoApplication*)(tpapp)); CHKERRQ(info);
349: return(0);
350: }
355: int TaoApplicationDestroy(TAO_APPLICATION taoapp){
356: int info;
358: info=TaoAppDestroy(taoapp); CHKERRQ(info);
359: return(0);
360: }