1 #include <petsc/private/taoimpl.h> /*I "petsctao.h" I*/
2
3 /*@
4 TaoSetInitialVector - Sets the initial guess for the solve
5
6 Logically collective on Tao
7
8 Input Parameters:
9 + tao - the Tao context
10 - x0 - the initial guess
11
12 Level: beginner
13 .seealso: TaoCreate(), TaoSolve()
14 @*/
15
TaoSetInitialVector(Tao tao,Vec x0)16 PetscErrorCode TaoSetInitialVector(Tao tao, Vec x0)
17 {
18 PetscErrorCode ierr;
19
20 PetscFunctionBegin;
21 PetscValidHeaderSpecific(tao,TAO_CLASSID,1);
22 if (x0) {
23 PetscValidHeaderSpecific(x0,VEC_CLASSID,2);
24 PetscObjectReference((PetscObject)x0);
25 }
26 ierr = VecDestroy(&tao->solution);CHKERRQ(ierr);
27 tao->solution = x0;
28 PetscFunctionReturn(0);
29 }
30
TaoTestGradient(Tao tao,Vec x,Vec g1)31 PetscErrorCode TaoTestGradient(Tao tao,Vec x,Vec g1)
32 {
33 Vec g2,g3;
34 PetscBool complete_print = PETSC_FALSE,test = PETSC_FALSE;
35 PetscReal hcnorm,fdnorm,hcmax,fdmax,diffmax,diffnorm;
36 PetscScalar dot;
37 MPI_Comm comm;
38 PetscViewer viewer,mviewer;
39 PetscViewerFormat format;
40 PetscInt tabs;
41 static PetscBool directionsprinted = PETSC_FALSE;
42 PetscErrorCode ierr;
43
44 PetscFunctionBegin;
45 ierr = PetscObjectOptionsBegin((PetscObject)tao);CHKERRQ(ierr);
46 ierr = PetscOptionsName("-tao_test_gradient","Compare hand-coded and finite difference Gradients","None",&test);CHKERRQ(ierr);
47 ierr = PetscOptionsViewer("-tao_test_gradient_view","View difference between hand-coded and finite difference Gradients element entries","None",&mviewer,&format,&complete_print);CHKERRQ(ierr);
48 ierr = PetscOptionsEnd();CHKERRQ(ierr);
49 if (!test) {
50 if (complete_print) {
51 ierr = PetscViewerDestroy(&mviewer);CHKERRQ(ierr);
52 }
53 PetscFunctionReturn(0);
54 }
55
56 ierr = PetscObjectGetComm((PetscObject)tao,&comm);CHKERRQ(ierr);
57 ierr = PetscViewerASCIIGetStdout(comm,&viewer);CHKERRQ(ierr);
58 ierr = PetscViewerASCIIGetTab(viewer, &tabs);CHKERRQ(ierr);
59 ierr = PetscViewerASCIISetTab(viewer, ((PetscObject)tao)->tablevel);CHKERRQ(ierr);
60 ierr = PetscViewerASCIIPrintf(viewer," ---------- Testing Gradient -------------\n");CHKERRQ(ierr);
61 if (!complete_print && !directionsprinted) {
62 ierr = PetscViewerASCIIPrintf(viewer," Run with -tao_test_gradient_view and optionally -tao_test_gradient <threshold> to show difference\n");CHKERRQ(ierr);
63 ierr = PetscViewerASCIIPrintf(viewer," of hand-coded and finite difference gradient entries greater than <threshold>.\n");CHKERRQ(ierr);
64 }
65 if (!directionsprinted) {
66 ierr = PetscViewerASCIIPrintf(viewer," Testing hand-coded Gradient, if (for double precision runs) ||G - Gfd||/||G|| is\n");CHKERRQ(ierr);
67 ierr = PetscViewerASCIIPrintf(viewer," O(1.e-8), the hand-coded Gradient is probably correct.\n");CHKERRQ(ierr);
68 directionsprinted = PETSC_TRUE;
69 }
70 if (complete_print) {
71 ierr = PetscViewerPushFormat(mviewer,format);CHKERRQ(ierr);
72 }
73
74 ierr = VecDuplicate(x,&g2);CHKERRQ(ierr);
75 ierr = VecDuplicate(x,&g3);CHKERRQ(ierr);
76
77 /* Compute finite difference gradient, assume the gradient is already computed by TaoComputeGradient() and put into g1 */
78 ierr = TaoDefaultComputeGradient(tao,x,g2,NULL);CHKERRQ(ierr);
79
80 ierr = VecNorm(g2,NORM_2,&fdnorm);CHKERRQ(ierr);
81 ierr = VecNorm(g1,NORM_2,&hcnorm);CHKERRQ(ierr);
82 ierr = VecNorm(g2,NORM_INFINITY,&fdmax);CHKERRQ(ierr);
83 ierr = VecNorm(g1,NORM_INFINITY,&hcmax);CHKERRQ(ierr);
84 ierr = VecDot(g1,g2,&dot);CHKERRQ(ierr);
85 ierr = VecCopy(g1,g3);CHKERRQ(ierr);
86 ierr = VecAXPY(g3,-1.0,g2);CHKERRQ(ierr);
87 ierr = VecNorm(g3,NORM_2,&diffnorm);CHKERRQ(ierr);
88 ierr = VecNorm(g3,NORM_INFINITY,&diffmax);CHKERRQ(ierr);
89 ierr = PetscViewerASCIIPrintf(viewer," ||Gfd|| %g, ||G|| = %g, angle cosine = (Gfd'G)/||Gfd||||G|| = %g\n", (double)fdnorm, (double)hcnorm, (double)(PetscRealPart(dot)/(fdnorm*hcnorm)));CHKERRQ(ierr);
90 ierr = PetscViewerASCIIPrintf(viewer," 2-norm ||G - Gfd||/||G|| = %g, ||G - Gfd|| = %g\n",(double)(diffnorm/PetscMax(hcnorm,fdnorm)),(double)diffnorm);CHKERRQ(ierr);
91 ierr = PetscViewerASCIIPrintf(viewer," max-norm ||G - Gfd||/||G|| = %g, ||G - Gfd|| = %g\n",(double)(diffmax/PetscMax(hcmax,fdmax)),(double)diffmax);CHKERRQ(ierr);
92
93 if (complete_print) {
94 ierr = PetscViewerASCIIPrintf(viewer," Hand-coded gradient ----------\n");CHKERRQ(ierr);
95 ierr = VecView(g1,mviewer);CHKERRQ(ierr);
96 ierr = PetscViewerASCIIPrintf(viewer," Finite difference gradient ----------\n");CHKERRQ(ierr);
97 ierr = VecView(g2,mviewer);CHKERRQ(ierr);
98 ierr = PetscViewerASCIIPrintf(viewer," Hand-coded minus finite-difference gradient ----------\n");CHKERRQ(ierr);
99 ierr = VecView(g3,mviewer);CHKERRQ(ierr);
100 }
101 ierr = VecDestroy(&g2);CHKERRQ(ierr);
102 ierr = VecDestroy(&g3);CHKERRQ(ierr);
103
104 if (complete_print) {
105 ierr = PetscViewerPopFormat(mviewer);CHKERRQ(ierr);
106 ierr = PetscViewerDestroy(&mviewer);CHKERRQ(ierr);
107 }
108 ierr = PetscViewerASCIISetTab(viewer,tabs);CHKERRQ(ierr);
109 PetscFunctionReturn(0);
110 }
111
112 /*@
113 TaoComputeGradient - Computes the gradient of the objective function
114
115 Collective on Tao
116
117 Input Parameters:
118 + tao - the Tao context
119 - X - input vector
120
121 Output Parameter:
122 . G - gradient vector
123
124 Options Database Keys:
125 + -tao_test_gradient - compare the user provided gradient with one compute via finite differences to check for errors
126 - -tao_test_gradient_view - display the user provided gradient, the finite difference gradient and the difference between them to help users detect the location of errors in the user provided gradient
127
128 Notes:
129 TaoComputeGradient() is typically used within minimization implementations,
130 so most users would not generally call this routine themselves.
131
132 Level: advanced
133
134 .seealso: TaoComputeObjective(), TaoComputeObjectiveAndGradient(), TaoSetGradientRoutine()
135 @*/
TaoComputeGradient(Tao tao,Vec X,Vec G)136 PetscErrorCode TaoComputeGradient(Tao tao, Vec X, Vec G)
137 {
138 PetscErrorCode ierr;
139 PetscReal dummy;
140
141 PetscFunctionBegin;
142 PetscValidHeaderSpecific(tao,TAO_CLASSID,1);
143 PetscValidHeaderSpecific(X,VEC_CLASSID,2);
144 PetscValidHeaderSpecific(G,VEC_CLASSID,2);
145 PetscCheckSameComm(tao,1,X,2);
146 PetscCheckSameComm(tao,1,G,3);
147 ierr = VecLockReadPush(X);CHKERRQ(ierr);
148 if (tao->ops->computegradient) {
149 ierr = PetscLogEventBegin(TAO_GradientEval,tao,X,G,NULL);CHKERRQ(ierr);
150 PetscStackPush("Tao user gradient evaluation routine");
151 ierr = (*tao->ops->computegradient)(tao,X,G,tao->user_gradP);CHKERRQ(ierr);
152 PetscStackPop;
153 ierr = PetscLogEventEnd(TAO_GradientEval,tao,X,G,NULL);CHKERRQ(ierr);
154 tao->ngrads++;
155 } else if (tao->ops->computeobjectiveandgradient) {
156 ierr = PetscLogEventBegin(TAO_ObjGradEval,tao,X,G,NULL);CHKERRQ(ierr);
157 PetscStackPush("Tao user objective/gradient evaluation routine");
158 ierr = (*tao->ops->computeobjectiveandgradient)(tao,X,&dummy,G,tao->user_objgradP);CHKERRQ(ierr);
159 PetscStackPop;
160 ierr = PetscLogEventEnd(TAO_ObjGradEval,tao,X,G,NULL);CHKERRQ(ierr);
161 tao->nfuncgrads++;
162 } else SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_ARG_WRONGSTATE,"TaoSetGradientRoutine() has not been called");
163 ierr = VecLockReadPop(X);CHKERRQ(ierr);
164
165 ierr = TaoTestGradient(tao,X,G);CHKERRQ(ierr);
166 PetscFunctionReturn(0);
167 }
168
169 /*@
170 TaoComputeObjective - Computes the objective function value at a given point
171
172 Collective on Tao
173
174 Input Parameters:
175 + tao - the Tao context
176 - X - input vector
177
178 Output Parameter:
179 . f - Objective value at X
180
181 Notes:
182 TaoComputeObjective() is typically used within minimization implementations,
183 so most users would not generally call this routine themselves.
184
185 Level: advanced
186
187 .seealso: TaoComputeGradient(), TaoComputeObjectiveAndGradient(), TaoSetObjectiveRoutine()
188 @*/
TaoComputeObjective(Tao tao,Vec X,PetscReal * f)189 PetscErrorCode TaoComputeObjective(Tao tao, Vec X, PetscReal *f)
190 {
191 PetscErrorCode ierr;
192 Vec temp;
193
194 PetscFunctionBegin;
195 PetscValidHeaderSpecific(tao,TAO_CLASSID,1);
196 PetscValidHeaderSpecific(X,VEC_CLASSID,2);
197 PetscCheckSameComm(tao,1,X,2);
198 ierr = VecLockReadPush(X);CHKERRQ(ierr);
199 if (tao->ops->computeobjective) {
200 ierr = PetscLogEventBegin(TAO_ObjectiveEval,tao,X,NULL,NULL);CHKERRQ(ierr);
201 PetscStackPush("Tao user objective evaluation routine");
202 ierr = (*tao->ops->computeobjective)(tao,X,f,tao->user_objP);CHKERRQ(ierr);
203 PetscStackPop;
204 ierr = PetscLogEventEnd(TAO_ObjectiveEval,tao,X,NULL,NULL);CHKERRQ(ierr);
205 tao->nfuncs++;
206 } else if (tao->ops->computeobjectiveandgradient) {
207 ierr = PetscInfo(tao,"Duplicating variable vector in order to call func/grad routine\n");CHKERRQ(ierr);
208 ierr = VecDuplicate(X,&temp);CHKERRQ(ierr);
209 ierr = PetscLogEventBegin(TAO_ObjGradEval,tao,X,NULL,NULL);CHKERRQ(ierr);
210 PetscStackPush("Tao user objective/gradient evaluation routine");
211 ierr = (*tao->ops->computeobjectiveandgradient)(tao,X,f,temp,tao->user_objgradP);CHKERRQ(ierr);
212 PetscStackPop;
213 ierr = PetscLogEventEnd(TAO_ObjGradEval,tao,X,NULL,NULL);CHKERRQ(ierr);
214 ierr = VecDestroy(&temp);CHKERRQ(ierr);
215 tao->nfuncgrads++;
216 } else SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_ARG_WRONGSTATE,"TaoSetObjectiveRoutine() has not been called");
217 ierr = PetscInfo1(tao,"TAO Function evaluation: %20.19e\n",(double)(*f));CHKERRQ(ierr);
218 ierr = VecLockReadPop(X);CHKERRQ(ierr);
219 PetscFunctionReturn(0);
220 }
221
222 /*@
223 TaoComputeObjectiveAndGradient - Computes the objective function value at a given point
224
225 Collective on Tao
226
227 Input Parameters:
228 + tao - the Tao context
229 - X - input vector
230
231 Output Parameter:
232 + f - Objective value at X
233 - g - Gradient vector at X
234
235 Notes:
236 TaoComputeObjectiveAndGradient() is typically used within minimization implementations,
237 so most users would not generally call this routine themselves.
238
239 Level: advanced
240
241 .seealso: TaoComputeGradient(), TaoComputeObjectiveAndGradient(), TaoSetObjectiveRoutine()
242 @*/
TaoComputeObjectiveAndGradient(Tao tao,Vec X,PetscReal * f,Vec G)243 PetscErrorCode TaoComputeObjectiveAndGradient(Tao tao, Vec X, PetscReal *f, Vec G)
244 {
245 PetscErrorCode ierr;
246
247 PetscFunctionBegin;
248 PetscValidHeaderSpecific(tao,TAO_CLASSID,1);
249 PetscValidHeaderSpecific(X,VEC_CLASSID,2);
250 PetscValidHeaderSpecific(G,VEC_CLASSID,4);
251 PetscCheckSameComm(tao,1,X,2);
252 PetscCheckSameComm(tao,1,G,4);
253 ierr = VecLockReadPush(X);CHKERRQ(ierr);
254 if (tao->ops->computeobjectiveandgradient) {
255 ierr = PetscLogEventBegin(TAO_ObjGradEval,tao,X,G,NULL);CHKERRQ(ierr);
256 if (tao->ops->computegradient == TaoDefaultComputeGradient) {
257 ierr = TaoComputeObjective(tao,X,f);CHKERRQ(ierr);
258 ierr = TaoDefaultComputeGradient(tao,X,G,NULL);CHKERRQ(ierr);
259 } else {
260 PetscStackPush("Tao user objective/gradient evaluation routine");
261 ierr = (*tao->ops->computeobjectiveandgradient)(tao,X,f,G,tao->user_objgradP);CHKERRQ(ierr);
262 PetscStackPop;
263 }
264 ierr = PetscLogEventEnd(TAO_ObjGradEval,tao,X,G,NULL);CHKERRQ(ierr);
265 tao->nfuncgrads++;
266 } else if (tao->ops->computeobjective && tao->ops->computegradient) {
267 ierr = PetscLogEventBegin(TAO_ObjectiveEval,tao,X,NULL,NULL);CHKERRQ(ierr);
268 PetscStackPush("Tao user objective evaluation routine");
269 ierr = (*tao->ops->computeobjective)(tao,X,f,tao->user_objP);CHKERRQ(ierr);
270 PetscStackPop;
271 ierr = PetscLogEventEnd(TAO_ObjectiveEval,tao,X,NULL,NULL);CHKERRQ(ierr);
272 tao->nfuncs++;
273 ierr = PetscLogEventBegin(TAO_GradientEval,tao,X,G,NULL);CHKERRQ(ierr);
274 PetscStackPush("Tao user gradient evaluation routine");
275 ierr = (*tao->ops->computegradient)(tao,X,G,tao->user_gradP);CHKERRQ(ierr);
276 PetscStackPop;
277 ierr = PetscLogEventEnd(TAO_GradientEval,tao,X,G,NULL);CHKERRQ(ierr);
278 tao->ngrads++;
279 } else SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_ARG_WRONGSTATE,"TaoSetObjectiveRoutine() or TaoSetGradientRoutine() not set");
280 ierr = PetscInfo1(tao,"TAO Function evaluation: %20.19e\n",(double)(*f));CHKERRQ(ierr);
281 ierr = VecLockReadPop(X);CHKERRQ(ierr);
282
283 ierr = TaoTestGradient(tao,X,G);CHKERRQ(ierr);
284 PetscFunctionReturn(0);
285 }
286
287 /*@C
288 TaoSetObjectiveRoutine - Sets the function evaluation routine for minimization
289
290 Logically collective on Tao
291
292 Input Parameter:
293 + tao - the Tao context
294 . func - the objective function
295 - ctx - [optional] user-defined context for private data for the function evaluation
296 routine (may be NULL)
297
298 Calling sequence of func:
299 $ func (Tao tao, Vec x, PetscReal *f, void *ctx);
300
301 + x - input vector
302 . f - function value
303 - ctx - [optional] user-defined function context
304
305 Level: beginner
306
307 .seealso: TaoSetGradientRoutine(), TaoSetHessianRoutine() TaoSetObjectiveAndGradientRoutine()
308 @*/
TaoSetObjectiveRoutine(Tao tao,PetscErrorCode (* func)(Tao,Vec,PetscReal *,void *),void * ctx)309 PetscErrorCode TaoSetObjectiveRoutine(Tao tao, PetscErrorCode (*func)(Tao, Vec, PetscReal*,void*),void *ctx)
310 {
311 PetscFunctionBegin;
312 PetscValidHeaderSpecific(tao,TAO_CLASSID,1);
313 tao->user_objP = ctx;
314 tao->ops->computeobjective = func;
315 PetscFunctionReturn(0);
316 }
317
318 /*@C
319 TaoSetResidualRoutine - Sets the residual evaluation routine for least-square applications
320
321 Logically collective on Tao
322
323 Input Parameter:
324 + tao - the Tao context
325 . func - the residual evaluation routine
326 - ctx - [optional] user-defined context for private data for the function evaluation
327 routine (may be NULL)
328
329 Calling sequence of func:
330 $ func (Tao tao, Vec x, Vec f, void *ctx);
331
332 + x - input vector
333 . f - function value vector
334 - ctx - [optional] user-defined function context
335
336 Level: beginner
337
338 .seealso: TaoSetObjectiveRoutine(), TaoSetJacobianRoutine()
339 @*/
TaoSetResidualRoutine(Tao tao,Vec res,PetscErrorCode (* func)(Tao,Vec,Vec,void *),void * ctx)340 PetscErrorCode TaoSetResidualRoutine(Tao tao, Vec res, PetscErrorCode (*func)(Tao, Vec, Vec, void*),void *ctx)
341 {
342 PetscErrorCode ierr;
343
344 PetscFunctionBegin;
345 PetscValidHeaderSpecific(tao,TAO_CLASSID,1);
346 PetscValidHeaderSpecific(res,VEC_CLASSID,2);
347 ierr = PetscObjectReference((PetscObject)res);CHKERRQ(ierr);
348 if (tao->ls_res) {
349 ierr = VecDestroy(&tao->ls_res);CHKERRQ(ierr);
350 }
351 tao->ls_res = res;
352 tao->user_lsresP = ctx;
353 tao->ops->computeresidual = func;
354
355 PetscFunctionReturn(0);
356 }
357
358 /*@
359 TaoSetResidualWeights - Give weights for the residual values. A vector can be used if only diagonal terms are used, otherwise a matrix can be give. If this function is not used, or if sigma_v and sigma_w are both NULL, then the default identity matrix will be used for weights.
360
361 Collective on Tao
362
363 Input Parameters:
364 + tao - the Tao context
365 . sigma_v - vector of weights (diagonal terms only)
366 . n - the number of weights (if using off-diagonal)
367 . rows - index list of rows for sigma_w
368 . cols - index list of columns for sigma_w
369 - vals - array of weights
370
371
372
373 Note: Either sigma_v or sigma_w (or both) should be NULL
374
375 Level: intermediate
376
377 .seealso: TaoSetResidualRoutine()
378 @*/
TaoSetResidualWeights(Tao tao,Vec sigma_v,PetscInt n,PetscInt * rows,PetscInt * cols,PetscReal * vals)379 PetscErrorCode TaoSetResidualWeights(Tao tao, Vec sigma_v, PetscInt n, PetscInt *rows, PetscInt *cols, PetscReal *vals)
380 {
381 PetscErrorCode ierr;
382 PetscInt i;
383 PetscFunctionBegin;
384 PetscValidHeaderSpecific(tao,TAO_CLASSID,1);
385 if (sigma_v) {
386 PetscValidHeaderSpecific(sigma_v,VEC_CLASSID,2);
387 ierr = PetscObjectReference((PetscObject)sigma_v);CHKERRQ(ierr);
388 }
389 if (tao->res_weights_v) {
390 ierr = VecDestroy(&tao->res_weights_v);CHKERRQ(ierr);
391 }
392 tao->res_weights_v=sigma_v;
393 if (vals) {
394 if (tao->res_weights_n) {
395 ierr = PetscFree(tao->res_weights_rows);CHKERRQ(ierr);
396 ierr = PetscFree(tao->res_weights_cols);CHKERRQ(ierr);
397 ierr = PetscFree(tao->res_weights_w);CHKERRQ(ierr);
398 }
399 ierr = PetscMalloc1(n,&tao->res_weights_rows);CHKERRQ(ierr);
400 ierr = PetscMalloc1(n,&tao->res_weights_cols);CHKERRQ(ierr);
401 ierr = PetscMalloc1(n,&tao->res_weights_w);CHKERRQ(ierr);
402 tao->res_weights_n=n;
403 for (i=0;i<n;i++) {
404 tao->res_weights_rows[i]=rows[i];
405 tao->res_weights_cols[i]=cols[i];
406 tao->res_weights_w[i]=vals[i];
407 }
408 } else {
409 tao->res_weights_n=0;
410 tao->res_weights_rows=NULL;
411 tao->res_weights_cols=NULL;
412 }
413 PetscFunctionReturn(0);
414 }
415
416 /*@
417 TaoComputeResidual - Computes a least-squares residual vector at a given point
418
419 Collective on Tao
420
421 Input Parameters:
422 + tao - the Tao context
423 - X - input vector
424
425 Output Parameter:
426 . f - Objective vector at X
427
428 Notes:
429 TaoComputeResidual() is typically used within minimization implementations,
430 so most users would not generally call this routine themselves.
431
432 Level: advanced
433
434 .seealso: TaoSetResidualRoutine()
435 @*/
TaoComputeResidual(Tao tao,Vec X,Vec F)436 PetscErrorCode TaoComputeResidual(Tao tao, Vec X, Vec F)
437 {
438 PetscErrorCode ierr;
439
440 PetscFunctionBegin;
441 PetscValidHeaderSpecific(tao,TAO_CLASSID,1);
442 PetscValidHeaderSpecific(X,VEC_CLASSID,2);
443 PetscValidHeaderSpecific(F,VEC_CLASSID,3);
444 PetscCheckSameComm(tao,1,X,2);
445 PetscCheckSameComm(tao,1,F,3);
446 if (tao->ops->computeresidual) {
447 ierr = PetscLogEventBegin(TAO_ObjectiveEval,tao,X,NULL,NULL);CHKERRQ(ierr);
448 PetscStackPush("Tao user least-squares residual evaluation routine");
449 ierr = (*tao->ops->computeresidual)(tao,X,F,tao->user_lsresP);CHKERRQ(ierr);
450 PetscStackPop;
451 ierr = PetscLogEventEnd(TAO_ObjectiveEval,tao,X,NULL,NULL);CHKERRQ(ierr);
452 tao->nfuncs++;
453 } else SETERRQ(PetscObjectComm((PetscObject)tao),PETSC_ERR_ARG_WRONGSTATE,"TaoSetResidualRoutine() has not been called");
454 ierr = PetscInfo(tao,"TAO least-squares residual evaluation.\n");CHKERRQ(ierr);
455 PetscFunctionReturn(0);
456 }
457
458 /*@C
459 TaoSetGradientRoutine - Sets the gradient evaluation routine for minimization
460
461 Logically collective on Tao
462
463 Input Parameter:
464 + tao - the Tao context
465 . func - the gradient function
466 - ctx - [optional] user-defined context for private data for the gradient evaluation
467 routine (may be NULL)
468
469 Calling sequence of func:
470 $ func (Tao tao, Vec x, Vec g, void *ctx);
471
472 + x - input vector
473 . g - gradient value (output)
474 - ctx - [optional] user-defined function context
475
476 Level: beginner
477
478 .seealso: TaoSetObjectiveRoutine(), TaoSetHessianRoutine() TaoSetObjectiveAndGradientRoutine()
479 @*/
TaoSetGradientRoutine(Tao tao,PetscErrorCode (* func)(Tao,Vec,Vec,void *),void * ctx)480 PetscErrorCode TaoSetGradientRoutine(Tao tao, PetscErrorCode (*func)(Tao, Vec, Vec, void*),void *ctx)
481 {
482 PetscFunctionBegin;
483 PetscValidHeaderSpecific(tao,TAO_CLASSID,1);
484 tao->user_gradP = ctx;
485 tao->ops->computegradient = func;
486 PetscFunctionReturn(0);
487 }
488
489 /*@C
490 TaoSetObjectiveAndGradientRoutine - Sets a combined objective function and gradient evaluation routine for minimization
491
492 Logically collective on Tao
493
494 Input Parameter:
495 + tao - the Tao context
496 . func - the gradient function
497 - ctx - [optional] user-defined context for private data for the gradient evaluation
498 routine (may be NULL)
499
500 Calling sequence of func:
501 $ func (Tao tao, Vec x, PetscReal *f, Vec g, void *ctx);
502
503 + x - input vector
504 . f - objective value (output)
505 . g - gradient value (output)
506 - ctx - [optional] user-defined function context
507
508 Level: beginner
509
510 .seealso: TaoSetObjectiveRoutine(), TaoSetHessianRoutine() TaoSetObjectiveAndGradientRoutine()
511 @*/
TaoSetObjectiveAndGradientRoutine(Tao tao,PetscErrorCode (* func)(Tao,Vec,PetscReal *,Vec,void *),void * ctx)512 PetscErrorCode TaoSetObjectiveAndGradientRoutine(Tao tao, PetscErrorCode (*func)(Tao, Vec, PetscReal *, Vec, void*), void *ctx)
513 {
514 PetscFunctionBegin;
515 PetscValidHeaderSpecific(tao,TAO_CLASSID,1);
516 tao->user_objgradP = ctx;
517 tao->ops->computeobjectiveandgradient = func;
518 PetscFunctionReturn(0);
519 }
520
521 /*@
522 TaoIsObjectiveDefined -- Checks to see if the user has
523 declared an objective-only routine. Useful for determining when
524 it is appropriate to call TaoComputeObjective() or
525 TaoComputeObjectiveAndGradient()
526
527 Collective on Tao
528
529 Input Parameter:
530 + tao - the Tao context
531 - ctx - PETSC_TRUE if objective function routine is set by user,
532 PETSC_FALSE otherwise
533 Level: developer
534
535 .seealso: TaoSetObjectiveRoutine(), TaoIsGradientDefined(), TaoIsObjectiveAndGradientDefined()
536 @*/
TaoIsObjectiveDefined(Tao tao,PetscBool * flg)537 PetscErrorCode TaoIsObjectiveDefined(Tao tao, PetscBool *flg)
538 {
539 PetscFunctionBegin;
540 PetscValidHeaderSpecific(tao,TAO_CLASSID,1);
541 if (tao->ops->computeobjective == NULL) *flg = PETSC_FALSE;
542 else *flg = PETSC_TRUE;
543 PetscFunctionReturn(0);
544 }
545
546 /*@
547 TaoIsGradientDefined -- Checks to see if the user has
548 declared an objective-only routine. Useful for determining when
549 it is appropriate to call TaoComputeGradient() or
550 TaoComputeGradientAndGradient()
551
552 Not Collective
553
554 Input Parameter:
555 + tao - the Tao context
556 - ctx - PETSC_TRUE if gradient routine is set by user, PETSC_FALSE otherwise
557 Level: developer
558
559 .seealso: TaoSetGradientRoutine(), TaoIsObjectiveDefined(), TaoIsObjectiveAndGradientDefined()
560 @*/
TaoIsGradientDefined(Tao tao,PetscBool * flg)561 PetscErrorCode TaoIsGradientDefined(Tao tao, PetscBool *flg)
562 {
563 PetscFunctionBegin;
564 PetscValidHeaderSpecific(tao,TAO_CLASSID,1);
565 if (tao->ops->computegradient == NULL) *flg = PETSC_FALSE;
566 else *flg = PETSC_TRUE;
567 PetscFunctionReturn(0);
568 }
569
570 /*@
571 TaoIsObjectiveAndGradientDefined -- Checks to see if the user has
572 declared a joint objective/gradient routine. Useful for determining when
573 it is appropriate to call TaoComputeObjective() or
574 TaoComputeObjectiveAndGradient()
575
576 Not Collective
577
578 Input Parameter:
579 + tao - the Tao context
580 - ctx - PETSC_TRUE if objective/gradient routine is set by user, PETSC_FALSE otherwise
581 Level: developer
582
583 .seealso: TaoSetObjectiveAndGradientRoutine(), TaoIsObjectiveDefined(), TaoIsGradientDefined()
584 @*/
TaoIsObjectiveAndGradientDefined(Tao tao,PetscBool * flg)585 PetscErrorCode TaoIsObjectiveAndGradientDefined(Tao tao, PetscBool *flg)
586 {
587 PetscFunctionBegin;
588 PetscValidHeaderSpecific(tao,TAO_CLASSID,1);
589 if (tao->ops->computeobjectiveandgradient == NULL) *flg = PETSC_FALSE;
590 else *flg = PETSC_TRUE;
591 PetscFunctionReturn(0);
592 }
593