Actual source code: taotermcallbacks.c
1: #include <petsc/private/taoimpl.h>
3: typedef struct _n_TaoTerm_Callbacks TaoTerm_Callbacks;
5: struct _n_TaoTerm_Callbacks {
6: Tao tao;
7: PetscErrorCode (*objective)(Tao, Vec, PetscReal *, PetscCtx);
8: PetscErrorCode (*gradient)(Tao, Vec, Vec, PetscCtx);
9: PetscErrorCode (*objectiveandgradient)(Tao, Vec, PetscReal *, Vec, PetscCtx);
10: PetscErrorCode (*hessian)(Tao, Vec, Mat, Mat, PetscCtx);
11: PetscCtx obj_ctx;
12: PetscCtx grad_ctx;
13: PetscCtx objgrad_ctx;
14: PetscCtx hess_ctx;
15: char *obj_name;
16: char *grad_name;
17: char *objgrad_name;
18: char *hess_name;
19: char *set_obj_name;
20: char *set_grad_name;
21: char *set_objgrad_name;
22: char *set_hess_name;
23: };
25: #define PetscCheckTaoTermCallbacksValid(term, tt, params) \
26: do { \
27: PetscCheck((params) == NULL, PetscObjectComm((PetscObject)(term)), PETSC_ERR_ARG_INCOMP, "TAOTERMCALLBACKS does not accept a vector of parameters"); \
28: PetscCheck((tt)->tao != NULL, PetscObjectComm((PetscObject)(term)), PETSC_ERR_ARG_WRONGSTATE, "TAOTERMCALLBACKS does not have an outer Tao"); \
29: } while (0)
31: static PetscErrorCode TaoTermDestroy_Callbacks(TaoTerm term)
32: {
33: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
35: PetscFunctionBegin;
36: // tt->tao is a weak reference, we do not destroy it
37: PetscCall(PetscFree(tt->obj_name));
38: PetscCall(PetscFree(tt->grad_name));
39: PetscCall(PetscFree(tt->objgrad_name));
40: PetscCall(PetscFree(tt->hess_name));
41: PetscCall(PetscFree(tt->set_obj_name));
42: PetscCall(PetscFree(tt->set_grad_name));
43: PetscCall(PetscFree(tt->set_objgrad_name));
44: PetscCall(PetscFree(tt->set_hess_name));
45: PetscCall(PetscFree(tt));
46: term->data = NULL;
47: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksSetObjective_C", NULL));
48: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksGetObjective_C", NULL));
49: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksSetGradient_C", NULL));
50: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksGetGradient_C", NULL));
51: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksSetObjectiveAndGradient_C", NULL));
52: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksGetObjectiveAndGradient_C", NULL));
53: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksSetHessian_C", NULL));
54: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksGetHessian_C", NULL));
55: PetscFunctionReturn(PETSC_SUCCESS);
56: }
58: static PetscErrorCode TaoTermComputeObjective_Callbacks(TaoTerm term, Vec x, Vec params, PetscReal *value)
59: {
60: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
62: PetscFunctionBegin;
63: PetscCheckTaoTermCallbacksValid(term, tt, params);
64: if (tt->objective) {
65: PetscCallBack(tt->obj_name, (*tt->objective)(tt->tao, x, value, tt->obj_ctx));
66: } else if (tt->objectiveandgradient) {
67: Vec dummy;
69: PetscCall(PetscInfo(term, "%s: Duplicating variable vector in order to call func/grad routine\n", ((PetscObject)term)->prefix));
70: PetscCall(VecDuplicate(x, &dummy));
71: PetscCallBack(tt->objgrad_name, (*tt->objectiveandgradient)(tt->tao, x, value, dummy, tt->objgrad_ctx));
72: PetscCall(VecDestroy(&dummy));
73: } else SETERRQ(PetscObjectComm((PetscObject)term), PETSC_ERR_ARG_WRONGSTATE, "Objective routine not set: call %s", tt->set_obj_name);
74: PetscFunctionReturn(PETSC_SUCCESS);
75: }
77: static PetscErrorCode TaoTermComputeGradient_Callbacks(TaoTerm term, Vec x, Vec params, Vec g)
78: {
79: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
81: PetscFunctionBegin;
82: PetscCheckTaoTermCallbacksValid(term, tt, params);
83: if (tt->gradient) {
84: PetscCallBack(tt->grad_name, (*tt->gradient)(tt->tao, x, g, tt->grad_ctx));
85: } else if (tt->objectiveandgradient) {
86: PetscReal dummy;
88: PetscCallBack(tt->objgrad_name, (*tt->objectiveandgradient)(tt->tao, x, &dummy, g, tt->objgrad_ctx));
89: } else SETERRQ(PetscObjectComm((PetscObject)term), PETSC_ERR_ARG_WRONGSTATE, "Gradient routine not set: call %s", tt->set_grad_name);
90: PetscFunctionReturn(PETSC_SUCCESS);
91: }
93: static PetscErrorCode TaoTermComputeObjectiveAndGradient_Callbacks(TaoTerm term, Vec x, Vec params, PetscReal *value, Vec g)
94: {
95: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
97: PetscFunctionBegin;
98: PetscCheckTaoTermCallbacksValid(term, tt, params);
99: if (tt->objectiveandgradient) {
100: PetscCallBack(tt->objgrad_name, (*tt->objectiveandgradient)(tt->tao, x, value, g, tt->objgrad_ctx));
101: } else if (tt->objective && tt->gradient) {
102: PetscCallBack(tt->obj_name, (*tt->objective)(tt->tao, x, value, tt->obj_ctx));
103: PetscCallBack(tt->grad_name, (*tt->gradient)(tt->tao, x, g, tt->grad_ctx));
104: } else SETERRQ(PetscObjectComm((PetscObject)term), PETSC_ERR_ARG_WRONGSTATE, "Objective/gradient routine not set: call %s", tt->set_objgrad_name);
105: PetscFunctionReturn(PETSC_SUCCESS);
106: }
108: static PetscErrorCode TaoTermComputeHessian_Callbacks(TaoTerm term, Vec x, Vec params, Mat H, Mat Hpre)
109: {
110: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
112: PetscFunctionBegin;
113: PetscCheckTaoTermCallbacksValid(term, tt, params);
114: PetscCheck(tt->hessian, PetscObjectComm((PetscObject)term), PETSC_ERR_ARG_WRONGSTATE, "Hessian routine not set: call %s", tt->set_hess_name);
115: PetscCallBack(tt->hess_name, (*tt->hessian)(tt->tao, x, H, Hpre, tt->hess_ctx));
116: PetscFunctionReturn(PETSC_SUCCESS);
117: }
119: static PetscErrorCode TaoTermView_Callbacks(TaoTerm term, PetscViewer viewer)
120: {
121: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
122: PetscBool iascii;
124: PetscFunctionBegin;
125: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
126: if (iascii) {
127: if (tt->tao == NULL) {
128: PetscCall(PetscViewerASCIIPrintf(viewer, "not attached to a Tao\n"));
129: } else {
130: const char *name = "[name omitted]";
131: const char *prefix;
133: if (!PetscCIEnabled) PetscCall(PetscObjectGetName((PetscObject)tt->tao, &name));
134: else if (tt->tao->hdr.name) name = tt->tao->hdr.name;
135: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)tt->tao, &prefix));
136: if (prefix) PetscCall(PetscViewerASCIIPrintf(viewer, "attached to Tao %s (%s)\n", name, prefix));
137: else PetscCall(PetscViewerASCIIPrintf(viewer, "attached to Tao %s\n", name));
138: }
139: }
140: PetscFunctionReturn(PETSC_SUCCESS);
141: }
143: PETSC_INTERN PetscErrorCode TaoTermCallbacksSetObjective(TaoTerm term, PetscErrorCode (*tao_obj)(Tao, Vec, PetscReal *, PetscCtx), PetscCtx ctx)
144: {
145: PetscFunctionBegin;
147: PetscTryMethod(term, "TaoTermCallbacksSetObjective_C", (TaoTerm, PetscErrorCode (*)(Tao, Vec, PetscReal *, PetscCtx), PetscCtx), (term, tao_obj, ctx));
148: PetscFunctionReturn(PETSC_SUCCESS);
149: }
151: static PetscErrorCode TaoTermCallbacksSetObjective_Callbacks(TaoTerm term, PetscErrorCode (*tao_obj)(Tao, Vec, PetscReal *, PetscCtx), PetscCtx ctx)
152: {
153: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
155: PetscFunctionBegin;
156: tt->objective = tao_obj;
157: tt->obj_ctx = ctx;
158: PetscFunctionReturn(PETSC_SUCCESS);
159: }
161: PETSC_INTERN PetscErrorCode TaoTermCallbacksGetObjective(TaoTerm term, PetscErrorCode (**tao_obj)(Tao, Vec, PetscReal *, PetscCtx), PetscCtxRt ctx)
162: {
163: PetscFunctionBegin;
165: if (tao_obj) *tao_obj = NULL;
166: if (ctx) *(void **)ctx = NULL;
167: PetscTryMethod(term, "TaoTermCallbacksGetObjective_C", (TaoTerm, PetscErrorCode (**)(Tao, Vec, PetscReal *, PetscCtx), PetscCtxRt), (term, tao_obj, ctx));
168: PetscFunctionReturn(PETSC_SUCCESS);
169: }
171: static PetscErrorCode TaoTermCallbacksGetObjective_Callbacks(TaoTerm term, PetscErrorCode (**tao_obj)(Tao, Vec, PetscReal *, PetscCtx), PetscCtxRt ctx)
172: {
173: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
175: PetscFunctionBegin;
176: if (tao_obj) *tao_obj = tt->objective;
177: if (ctx) *(void **)ctx = tt->obj_ctx;
178: PetscFunctionReturn(PETSC_SUCCESS);
179: }
181: PETSC_INTERN PetscErrorCode TaoTermCallbacksSetGradient(TaoTerm term, PetscErrorCode (*tao_grad)(Tao, Vec, Vec, PetscCtx), PetscCtx ctx)
182: {
183: PetscFunctionBegin;
185: PetscTryMethod(term, "TaoTermCallbacksSetGradient_C", (TaoTerm, PetscErrorCode (*)(Tao, Vec, Vec, PetscCtx), PetscCtx), (term, tao_grad, ctx));
186: PetscFunctionReturn(PETSC_SUCCESS);
187: }
189: static PetscErrorCode TaoTermCallbacksSetGradient_Callbacks(TaoTerm term, PetscErrorCode (*tao_grad)(Tao, Vec, Vec, PetscCtx), PetscCtx ctx)
190: {
191: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
193: PetscFunctionBegin;
194: tt->gradient = tao_grad;
195: tt->grad_ctx = ctx;
196: PetscFunctionReturn(PETSC_SUCCESS);
197: }
199: PETSC_INTERN PetscErrorCode TaoTermCallbacksGetGradient(TaoTerm term, PetscErrorCode (**tao_grad)(Tao, Vec, Vec, PetscCtx), PetscCtxRt ctx)
200: {
201: PetscFunctionBegin;
203: if (tao_grad) *tao_grad = NULL;
204: if (ctx) *(void **)ctx = NULL;
205: PetscTryMethod(term, "TaoTermCallbacksGetGradient_C", (TaoTerm, PetscErrorCode (**)(Tao, Vec, Vec, PetscCtx), PetscCtxRt), (term, tao_grad, ctx));
206: PetscFunctionReturn(PETSC_SUCCESS);
207: }
209: static PetscErrorCode TaoTermCallbacksGetGradient_Callbacks(TaoTerm term, PetscErrorCode (**tao_grad)(Tao, Vec, Vec, PetscCtx), PetscCtxRt ctx)
210: {
211: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
213: PetscFunctionBegin;
214: if (tao_grad) *tao_grad = tt->gradient;
215: if (ctx) *(void **)ctx = tt->grad_ctx;
216: PetscFunctionReturn(PETSC_SUCCESS);
217: }
219: PETSC_INTERN PetscErrorCode TaoTermCallbacksSetObjectiveAndGradient(TaoTerm term, PetscErrorCode (*tao_objgrad)(Tao, Vec, PetscReal *, Vec, PetscCtx), PetscCtx ctx)
220: {
221: PetscFunctionBegin;
223: PetscTryMethod(term, "TaoTermCallbacksSetObjectiveAndGradient_C", (TaoTerm, PetscErrorCode (*)(Tao, Vec, PetscReal *, Vec, PetscCtx), PetscCtx), (term, tao_objgrad, ctx));
224: PetscFunctionReturn(PETSC_SUCCESS);
225: }
227: static PetscErrorCode TaoTermCallbacksSetObjectiveAndGradient_Callbacks(TaoTerm term, PetscErrorCode (*tao_objgrad)(Tao, Vec, PetscReal *, Vec, PetscCtx), PetscCtx ctx)
228: {
229: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
231: PetscFunctionBegin;
232: tt->objectiveandgradient = tao_objgrad;
233: tt->objgrad_ctx = ctx;
234: PetscFunctionReturn(PETSC_SUCCESS);
235: }
237: PETSC_INTERN PetscErrorCode TaoTermCallbacksGetObjectiveAndGradient(TaoTerm term, PetscErrorCode (**tao_objgrad)(Tao, Vec, PetscReal *, Vec, PetscCtx), PetscCtxRt ctx)
238: {
239: PetscFunctionBegin;
241: if (tao_objgrad) *tao_objgrad = NULL;
242: if (ctx) *(void **)ctx = NULL;
243: PetscTryMethod(term, "TaoTermCallbacksGetObjectiveAndGradient_C", (TaoTerm, PetscErrorCode (**)(Tao, Vec, PetscReal *, Vec, PetscCtx), PetscCtxRt), (term, tao_objgrad, ctx));
244: PetscFunctionReturn(PETSC_SUCCESS);
245: }
247: static PetscErrorCode TaoTermCallbacksGetObjectiveAndGradient_Callbacks(TaoTerm term, PetscErrorCode (**tao_objgrad)(Tao, Vec, PetscReal *, Vec, PetscCtx), PetscCtxRt ctx)
248: {
249: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
251: PetscFunctionBegin;
252: if (tao_objgrad) *tao_objgrad = tt->objectiveandgradient;
253: if (ctx) *(void **)ctx = tt->objgrad_ctx;
254: PetscFunctionReturn(PETSC_SUCCESS);
255: }
257: PETSC_INTERN PetscErrorCode TaoTermCallbacksSetHessian(TaoTerm term, PetscErrorCode (*tao_hess)(Tao, Vec, Mat, Mat, PetscCtx), PetscCtx ctx)
258: {
259: PetscFunctionBegin;
261: PetscTryMethod(term, "TaoTermCallbacksSetHessian_C", (TaoTerm, PetscErrorCode (*)(Tao, Vec, Mat, Mat, PetscCtx), PetscCtx), (term, tao_hess, ctx));
262: PetscFunctionReturn(PETSC_SUCCESS);
263: }
265: static PetscErrorCode TaoTermCallbacksSetHessian_Callbacks(TaoTerm term, PetscErrorCode (*tao_hess)(Tao, Vec, Mat, Mat, PetscCtx), PetscCtx ctx)
266: {
267: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
269: PetscFunctionBegin;
270: tt->hessian = tao_hess;
271: tt->hess_ctx = ctx;
272: PetscFunctionReturn(PETSC_SUCCESS);
273: }
275: PETSC_INTERN PetscErrorCode TaoTermCallbacksGetHessian(TaoTerm term, PetscErrorCode (**tao_hess)(Tao, Vec, Mat, Mat, PetscCtx), PetscCtxRt ctx)
276: {
277: PetscFunctionBegin;
279: if (tao_hess) *tao_hess = NULL;
280: if (ctx) *(void **)ctx = NULL;
281: PetscTryMethod(term, "TaoTermCallbacksGetHessian_C", (TaoTerm, PetscErrorCode (**)(Tao, Vec, Mat, Mat, PetscCtx), PetscCtxRt), (term, tao_hess, ctx));
282: PetscFunctionReturn(PETSC_SUCCESS);
283: }
285: static PetscErrorCode TaoTermCallbacksGetHessian_Callbacks(TaoTerm term, PetscErrorCode (**tao_hess)(Tao, Vec, Mat, Mat, PetscCtx), PetscCtxRt ctx)
286: {
287: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
289: PetscFunctionBegin;
290: if (tao_hess) *tao_hess = tt->hessian;
291: if (ctx) *(void **)ctx = tt->hess_ctx;
292: PetscFunctionReturn(PETSC_SUCCESS);
293: }
295: static PetscErrorCode TaoTermIsObjectiveDefined_Callbacks(TaoTerm term, PetscBool *flg)
296: {
297: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
299: PetscFunctionBegin;
300: *flg = (tt->objective != NULL) ? PETSC_TRUE : PETSC_FALSE;
301: PetscFunctionReturn(PETSC_SUCCESS);
302: }
304: static PetscErrorCode TaoTermIsGradientDefined_Callbacks(TaoTerm term, PetscBool *flg)
305: {
306: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
308: PetscFunctionBegin;
309: *flg = (tt->gradient != NULL) ? PETSC_TRUE : PETSC_FALSE;
310: PetscFunctionReturn(PETSC_SUCCESS);
311: }
313: static PetscErrorCode TaoTermIsObjectiveAndGradientDefined_Callbacks(TaoTerm term, PetscBool *flg)
314: {
315: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
317: PetscFunctionBegin;
318: *flg = (tt->objectiveandgradient != NULL) ? PETSC_TRUE : PETSC_FALSE;
319: PetscFunctionReturn(PETSC_SUCCESS);
320: }
322: static PetscErrorCode TaoTermIsHessianDefined_Callbacks(TaoTerm term, PetscBool *flg)
323: {
324: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)term->data;
326: PetscFunctionBegin;
327: *flg = (tt->hessian != NULL) ? PETSC_TRUE : PETSC_FALSE;
328: PetscFunctionReturn(PETSC_SUCCESS);
329: }
331: static PetscErrorCode TaoTermCreateHessianMatrices_Callbacks(TaoTerm term, Mat *H, Mat *Hpre)
332: {
333: PetscFunctionBegin;
334: if (H) *H = NULL;
335: if (Hpre) *Hpre = NULL;
336: PetscFunctionReturn(PETSC_SUCCESS);
337: }
339: static PetscErrorCode TaoTermCreate_Callbacks_Internal(TaoTerm term, const char obj[], const char set_obj[], const char grad[], const char set_grad[], const char objgrad[], const char set_objgrad[], const char hess[], const char set_hess[])
340: {
341: TaoTerm_Callbacks *tt;
342: char buf[256];
343: size_t len = PETSC_STATIC_ARRAY_LENGTH(buf);
345: PetscFunctionBegin;
346: term->parameters_mode = TAOTERM_PARAMETERS_NONE;
348: PetscCall(PetscNew(&tt));
349: term->data = (void *)tt;
351: term->ops->destroy = TaoTermDestroy_Callbacks;
352: term->ops->objective = TaoTermComputeObjective_Callbacks;
353: term->ops->gradient = TaoTermComputeGradient_Callbacks;
354: term->ops->objectiveandgradient = TaoTermComputeObjectiveAndGradient_Callbacks;
355: term->ops->hessian = TaoTermComputeHessian_Callbacks;
356: term->ops->createhessianmatrices = TaoTermCreateHessianMatrices_Callbacks;
357: term->ops->view = TaoTermView_Callbacks;
358: term->ops->isobjectivedefined = TaoTermIsObjectiveDefined_Callbacks;
359: term->ops->isgradientdefined = TaoTermIsGradientDefined_Callbacks;
360: term->ops->isobjectiveandgradientdefined = TaoTermIsObjectiveAndGradientDefined_Callbacks;
361: term->ops->ishessiandefined = TaoTermIsHessianDefined_Callbacks;
362: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksSetObjective_C", TaoTermCallbacksSetObjective_Callbacks));
363: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksGetObjective_C", TaoTermCallbacksGetObjective_Callbacks));
364: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksSetGradient_C", TaoTermCallbacksSetGradient_Callbacks));
365: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksGetGradient_C", TaoTermCallbacksGetGradient_Callbacks));
366: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksSetObjectiveAndGradient_C", TaoTermCallbacksSetObjectiveAndGradient_Callbacks));
367: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksGetObjectiveAndGradient_C", TaoTermCallbacksGetObjectiveAndGradient_Callbacks));
368: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksSetHessian_C", TaoTermCallbacksSetHessian_Callbacks));
369: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermCallbacksGetHessian_C", TaoTermCallbacksGetHessian_Callbacks));
371: PetscCall(PetscSNPrintf(buf, len, "%s callback", obj ? obj : "unknown objective"));
372: PetscCall(PetscStrallocpy(buf, &tt->obj_name));
373: PetscCall(PetscStrallocpy(set_obj, &tt->set_obj_name));
375: PetscCall(PetscSNPrintf(buf, len, "%s callback", grad ? grad : "unknown gradient"));
376: PetscCall(PetscStrallocpy(buf, &tt->grad_name));
377: PetscCall(PetscStrallocpy(set_grad, &tt->set_grad_name));
379: PetscCall(PetscSNPrintf(buf, len, "%s callback", objgrad ? objgrad : "unknown objective/gradient"));
380: PetscCall(PetscStrallocpy(buf, &tt->objgrad_name));
381: PetscCall(PetscStrallocpy(set_objgrad, &tt->set_objgrad_name));
383: PetscCall(PetscSNPrintf(buf, len, "%s callback", hess ? hess : "unknown hessian"));
384: PetscCall(PetscStrallocpy(buf, &tt->hess_name));
385: PetscCall(PetscStrallocpy(set_hess, &tt->set_hess_name));
386: PetscFunctionReturn(PETSC_SUCCESS);
387: }
389: /*MC
390: TAOTERMCALLBACKS - A `TaoTerm` implementation that accesses the callback functions that have
391: been provided with in `TaoSetObjective()`, `TaoSetGradient()`,
392: `TaoSetObjectiveAndGradient()`, and `TaoSetHessian()`.
394: Level: developer
396: Notes:
397: If you are interested in creating your own term, you should not use this. Use
398: `TAOTERMSHELL` or create your own implementation of `TaoTerm` with
399: `TaoTermRegister()`.
401: A `TAOTERMCALLBACKS` is always `TAOTERM_PARAMETERS_NONE`, so the `params`
402: argument of `TaoTerm` evaluation routines should always be `NULL`.
404: A `TAOTERMCALLBACKS` cannot create Hessian matrices; the user needs to pass the
405: Hessian matrices used in algorithms in `TaoSetHessian()`.
407: Developer Notes:
408: Internally each `Tao` has a `TaoTerm` of type `TAOTERMCALLBACKS` that is updated
409: by the `Tao` callback routines (`TaoSetObjective()`, `TaoSetGradient()`,
410: `TaoSetObjectiveAndGradient()`, and `TaoSetHessian()`).
412: The routines that get the user-defined `Tao` callback functions
413: (`TaoGetObjective()`, `TaoGetObjectiveAndGradient()`, `TaoGetGradient()`,
414: `TaoGetHessian()`) will always return those original callbacks, even if the
415: objective function has been changed by `TaoAddTerm()`,
416: so PETSc/TAO should not assume that those callbacks are valid in any library code.
418: A `TAOTERMCALLBACKS` has a weak-reference to the `Tao` that created it,
419: which may not be the `Tao` currently using it because the term could have been shared
420: using `TaoGetTerm()` and `TaoAddTerm()`.
422: .seealso: [](sec_tao_term),
423: `TaoTerm`,
424: `TaoTermType`
425: M*/
426: PETSC_INTERN PetscErrorCode TaoTermCreate_Callbacks(TaoTerm term)
427: {
428: PetscFunctionBegin;
429: // clang-format off
430: PetscCall(TaoTermCreate_Callbacks_Internal(term,
431: "TaoComputeObjective()", "TaoSetObjective()",
432: "TaoComputeGradient()", "TaoSetGradient()",
433: "TaoComputeObjectiveAndGradient()", "TaoSetObjectiveAndGradient()",
434: "TaoComputeHessian()", "TaoSetHessian()"));
435: // clang-format on
436: PetscFunctionReturn(PETSC_SUCCESS);
437: }
439: PETSC_INTERN PetscErrorCode TaoTermCreateCallbacks(Tao tao, TaoTerm *term)
440: {
441: PetscFunctionBegin;
442: PetscCall(TaoTermCreate(PetscObjectComm((PetscObject)tao), term));
443: PetscCall(TaoTermSetType(*term, TAOTERMCALLBACKS));
444: {
445: TaoTerm_Callbacks *tt = (TaoTerm_Callbacks *)((*term)->data);
447: tt->tao = tao; // weak reference, do not increment reference count
448: }
449: PetscFunctionReturn(PETSC_SUCCESS);
450: }