Actual source code: taotermshell.c
1: #include <petsc/private/taoimpl.h>
3: typedef struct _n_TaoTerm_Shell TaoTerm_Shell;
5: struct _n_TaoTerm_Shell {
6: PetscContainer ctxcontainer;
7: PetscBool3 iscomputehessianfdpossible;
8: };
10: /*@C
11: TaoTermShellSetContextDestroy - Set a method to destroy the context resources when a `TAOTERMSHELL` is destroyed
13: Logically collective
15: Input Parameters:
16: + term - a `TaoTerm` of type `TAOTERMSHELL`
17: - destroy - the context destroy function
19: Level: intermediate
21: .seealso: [](sec_tao_term), `TaoTerm`, `TAOTERMSHELL`, `TaoTermShellSetContext()`, `TaoTermShellGetContext()`
22: @*/
23: PetscErrorCode TaoTermShellSetContextDestroy(TaoTerm term, PetscCtxDestroyFn *destroy)
24: {
25: PetscFunctionBegin;
27: PetscTryMethod(term, "TaoTermShellSetContextDestroy_C", (TaoTerm, PetscCtxDestroyFn *), (term, destroy));
28: PetscFunctionReturn(PETSC_SUCCESS);
29: }
31: static PetscErrorCode TaoTermShellSetContextDestroy_Shell(TaoTerm term, PetscCtxDestroyFn *f)
32: {
33: TaoTerm_Shell *shell = (TaoTerm_Shell *)term->data;
35: PetscFunctionBegin;
36: if (shell->ctxcontainer) PetscCall(PetscContainerSetCtxDestroy(shell->ctxcontainer, f));
37: PetscFunctionReturn(PETSC_SUCCESS);
38: }
40: /*@
41: TaoTermShellGetContext - Get the context for a `TAOTERMSHELL`
43: Not collective
45: Input Parameter:
46: . term - a `TaoTerm` of type `TAOTERMSHELL`
48: Output Parameter:
49: . ctx - a context
51: Level: intermediate
53: .seealso: [](sec_tao_term), `TaoTerm`, `TAOTERMSHELL`, `TaoTermShellSetContext()`, `TaoTermShellSetContextDestroy()`
54: @*/
55: PetscErrorCode TaoTermShellGetContext(TaoTerm term, PetscCtxRt ctx)
56: {
57: PetscFunctionBegin;
59: PetscAssertPointer(ctx, 2);
60: PetscUseMethod(term, "TaoTermShellGetContext_C", (TaoTerm, PetscCtxRt), (term, ctx));
61: PetscFunctionReturn(PETSC_SUCCESS);
62: }
64: static PetscErrorCode TaoTermShellGetContext_Shell(TaoTerm term, PetscCtxRt ctx)
65: {
66: TaoTerm_Shell *shell = (TaoTerm_Shell *)term->data;
68: PetscFunctionBegin;
69: if (shell->ctxcontainer) PetscCall(PetscContainerGetPointer(shell->ctxcontainer, (void **)ctx));
70: else *(void **)ctx = NULL;
71: PetscFunctionReturn(PETSC_SUCCESS);
72: }
74: /*@
75: TaoTermShellSetContext - Set a context for a `TAOTERMSHELL`
77: Logically collective
79: Input Parameters:
80: + term - a `TaoTerm` of type `TAOTERMSHELL`
81: - ctx - a context
83: Level: intermediate
85: Note:
86: The context can be accessed in callbacks using `TaoTermShellGetContext()`
88: .seealso: [](sec_tao_term), `TaoTerm`, `TAOTERMSHELL`, `TaoTermShellGetContext()`, `TaoTermShellSetContextDestroy()`
89: @*/
90: PetscErrorCode TaoTermShellSetContext(TaoTerm term, PetscCtx ctx)
91: {
92: PetscFunctionBegin;
94: PetscTryMethod(term, "TaoTermShellSetContext_C", (TaoTerm, PetscCtx), (term, ctx));
95: PetscFunctionReturn(PETSC_SUCCESS);
96: }
98: static PetscErrorCode TaoTermShellSetContext_Shell(TaoTerm term, PetscCtx ctx)
99: {
100: TaoTerm_Shell *shell = (TaoTerm_Shell *)term->data;
102: PetscFunctionBegin;
103: if (ctx) {
104: PetscContainer ctxcontainer;
106: PetscCall(PetscContainerCreate(PetscObjectComm((PetscObject)term), &ctxcontainer));
107: PetscCall(PetscContainerSetPointer(ctxcontainer, ctx));
108: PetscCall(PetscObjectCompose((PetscObject)term, "TaoTermShell ctx", (PetscObject)ctxcontainer));
109: shell->ctxcontainer = ctxcontainer;
110: PetscCall(PetscContainerDestroy(&ctxcontainer));
111: } else {
112: PetscCall(PetscObjectCompose((PetscObject)term, "TaoTermShell ctx", NULL));
113: shell->ctxcontainer = NULL;
114: }
115: PetscFunctionReturn(PETSC_SUCCESS);
116: }
118: static PetscErrorCode TaoTermView_Shell(TaoTerm term, PetscViewer viewer)
119: {
120: TaoTerm_Shell *shell = (TaoTerm_Shell *)term->data;
121: PetscBool iascii;
123: PetscFunctionBegin;
124: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &iascii));
125: if (iascii) {
126: PetscBool any;
127: if (shell->ctxcontainer) PetscCall(PetscViewerASCIIPrintf(viewer, "Context has been set\n"));
128: else PetscCall(PetscViewerASCIIPrintf(viewer, "No context has been set\n"));
130: PetscCall(PetscViewerASCIIPrintf(viewer, "The following methods have been set:"));
131: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
132: any = PETSC_FALSE;
133: if (term->ops->objective) {
134: any = PETSC_TRUE;
135: PetscCall(PetscViewerASCIIPrintf(viewer, " objective,"));
136: }
137: if (term->ops->gradient) {
138: any = PETSC_TRUE;
139: PetscCall(PetscViewerASCIIPrintf(viewer, " gradient,"));
140: }
141: if (term->ops->objectiveandgradient) {
142: any = PETSC_TRUE;
143: PetscCall(PetscViewerASCIIPrintf(viewer, " objectiveandgradient,"));
144: }
145: if (term->ops->hessian) {
146: any = PETSC_TRUE;
147: PetscCall(PetscViewerASCIIPrintf(viewer, " hessian,"));
148: }
149: if (any == PETSC_FALSE) PetscCall(PetscViewerASCIIPrintf(viewer, " (none)"));
150: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
151: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
152: }
153: PetscFunctionReturn(PETSC_SUCCESS);
154: }
156: static PetscErrorCode TaoTermDestroy_Shell(TaoTerm term)
157: {
158: TaoTerm_Shell *shell = (TaoTerm_Shell *)term->data;
160: PetscFunctionBegin;
161: PetscCall(TaoTermShellSetContext_Shell(term, NULL));
162: PetscCall(PetscFree(shell));
163: term->data = NULL;
164: term->ops->objective = NULL;
165: term->ops->gradient = NULL;
166: term->ops->objectiveandgradient = NULL;
167: term->ops->hessian = NULL;
168: term->ops->view = NULL;
170: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetContextDestroy_C", NULL));
171: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetContext_C", NULL));
172: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellGetContext_C", NULL));
173: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetObjective_C", NULL));
174: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetGradient_C", NULL));
175: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetObjectiveAndGradient_C", NULL));
176: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetHessian_C", NULL));
177: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetIsComputeHessianFDPossible_C", NULL));
178: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetView_C", NULL));
179: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetCreateSolutionVec_C", NULL));
180: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetCreateParametersVec_C", NULL));
181: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetCreateHessianMatrices_C", NULL));
182: PetscFunctionReturn(PETSC_SUCCESS);
183: }
185: /*@C
186: TaoTermShellSetObjective - Set the objective function of a `TAOTERMSHELL`
188: Logically collective
190: Input Parameters:
191: + term - a `TaoTerm` of type `TAOTERMSHELL`
192: - objective - a `TaoTermObjectiveFn` function pointer
194: Level: intermediate
196: .seealso: [](sec_tao_term), `TaoTerm`, `TAOTERMSHELL`, `TaoTermShellGetContext()`, `TaoTermShellSetContextDestroy()`,
197: `TaoTermShellSetGradient()`,
198: `TaoTermShellSetObjectiveAndGradient()`,
199: `TaoTermShellSetHessian()`,
200: `TaoTermShellSetView()`,
201: `TaoTermObjectiveFn`
202: @*/
203: PetscErrorCode TaoTermShellSetObjective(TaoTerm term, TaoTermObjectiveFn *objective)
204: {
205: PetscFunctionBegin;
207: PetscTryMethod(term, "TaoTermShellSetObjective_C", (TaoTerm, TaoTermObjectiveFn *), (term, objective));
208: PetscFunctionReturn(PETSC_SUCCESS);
209: }
211: static PetscErrorCode TaoTermShellSetObjective_Shell(TaoTerm term, TaoTermObjectiveFn *objective)
212: {
213: PetscFunctionBegin;
214: term->ops->objective = objective;
215: PetscFunctionReturn(PETSC_SUCCESS);
216: }
218: /*@C
219: TaoTermShellSetGradient - Set the gradient function of a `TAOTERMSHELL`
221: Logically collective
223: Input Parameters:
224: + term - a `TaoTerm` of type `TAOTERMSHELL`
225: - gradient - a `TaoTermGradientFn` function pointer
227: Level: intermediate
229: .seealso: [](sec_tao_term), `TaoTerm`, `TAOTERMSHELL`, `TaoTermShellGetContext()`, `TaoTermShellSetContextDestroy()`,
230: `TaoTermShellSetObjective()`,
231: `TaoTermShellSetObjectiveAndGradient()`,
232: `TaoTermShellSetHessian()`,
233: `TaoTermShellSetView()`,
234: `TaoTermGradientFn`
235: @*/
236: PetscErrorCode TaoTermShellSetGradient(TaoTerm term, TaoTermGradientFn *gradient)
237: {
238: PetscFunctionBegin;
240: PetscTryMethod(term, "TaoTermShellSetGradient_C", (TaoTerm, TaoTermGradientFn *), (term, gradient));
241: PetscFunctionReturn(PETSC_SUCCESS);
242: }
244: static PetscErrorCode TaoTermShellSetGradient_Shell(TaoTerm term, TaoTermGradientFn *gradient)
245: {
246: PetscFunctionBegin;
247: term->ops->gradient = gradient;
248: PetscFunctionReturn(PETSC_SUCCESS);
249: }
251: /*@C
252: TaoTermShellSetObjectiveAndGradient - Set the objective and gradient function of a `TAOTERMSHELL`
254: Logically collective
256: Input Parameters:
257: + term - a `TaoTerm` of type `TAOTERMSHELL`
258: - objandgrad - a `TaoTermObjectiveAndGradientFn` function pointer
260: Level: intermediate
262: .seealso: [](sec_tao_term), `TaoTerm`, `TAOTERMSHELL`, `TaoTermShellGetContext()`, `TaoTermShellSetContextDestroy()`,
263: `TaoTermShellSetObjective()`,
264: `TaoTermShellSetGradient()`,
265: `TaoTermShellSetHessian()`,
266: `TaoTermShellSetView()`,
267: `TaoTermObjectiveAndGradientFn`
268: @*/
269: PetscErrorCode TaoTermShellSetObjectiveAndGradient(TaoTerm term, TaoTermObjectiveAndGradientFn *objandgrad)
270: {
271: PetscFunctionBegin;
273: PetscTryMethod(term, "TaoTermShellSetObjectiveAndGradient_C", (TaoTerm, TaoTermObjectiveAndGradientFn *), (term, objandgrad));
274: PetscFunctionReturn(PETSC_SUCCESS);
275: }
277: static PetscErrorCode TaoTermShellSetObjectiveAndGradient_Shell(TaoTerm term, TaoTermObjectiveAndGradientFn *objandgrad)
278: {
279: PetscFunctionBegin;
280: term->ops->objectiveandgradient = objandgrad;
281: PetscFunctionReturn(PETSC_SUCCESS);
282: }
284: /*@C
285: TaoTermShellSetHessian - Set the Hessian function of a `TAOTERMSHELL`
287: Logically collective
289: Input Parameters:
290: + term - a `TaoTerm` of type `TAOTERMSHELL`
291: - hessian - a `TaoTermHessianFn` function pointer
293: Level: intermediate
295: .seealso: [](sec_tao_term), `TaoTerm`, `TAOTERMSHELL`, `TaoTermShellGetContext()`, `TaoTermShellSetContextDestroy()`,
296: `TaoTermShellSetObjective()`,
297: `TaoTermShellSetGradient()`,
298: `TaoTermShellSetObjectiveAndGradient()`,
299: `TaoTermShellSetView()`,
300: `TaoTermHessianFn`
301: @*/
302: PetscErrorCode TaoTermShellSetHessian(TaoTerm term, TaoTermHessianFn *hessian)
303: {
304: PetscFunctionBegin;
306: PetscTryMethod(term, "TaoTermShellSetHessian_C", (TaoTerm, TaoTermHessianFn *), (term, hessian));
307: PetscFunctionReturn(PETSC_SUCCESS);
308: }
310: static PetscErrorCode TaoTermShellSetHessian_Shell(TaoTerm term, TaoTermHessianFn *hessian)
311: {
312: PetscFunctionBegin;
313: term->ops->hessian = hessian;
314: PetscFunctionReturn(PETSC_SUCCESS);
315: }
317: static PetscErrorCode TaoTermIsComputeHessianFDPossible_Shell(TaoTerm term, PetscBool3 *ispossible)
318: {
319: TaoTerm_Shell *shell = (TaoTerm_Shell *)term->data;
321: PetscFunctionBegin;
322: *ispossible = shell->iscomputehessianfdpossible;
323: PetscFunctionReturn(PETSC_SUCCESS);
324: }
326: /*@C
327: TaoTermShellSetIsComputeHessianFDPossible - Set whether this term can compute Hessian with finite differences for a `TAOTERMSHELL`
329: Logically collective
331: Input Parameters:
332: + term - a `TaoTerm` of type `TAOTERMSHELL`
333: - ispossible - whether Hessian computation with finite differences is possible
335: Level: intermediate
337: .seealso: [](sec_tao_term), `TaoTerm`, `TAOTERMSHELL`, `TaoTermShellGetContext()`, `TaoTermShellSetContextDestroy()`,
338: `TaoTermShellSetObjective()`,
339: `TaoTermShellSetGradient()`,
340: `TaoTermShellSetObjectiveAndGradient()`,
341: `TaoTermShellSetHessian()`,
342: `TaoTermIsComputeHessianFDPossible()`
343: @*/
344: PetscErrorCode TaoTermShellSetIsComputeHessianFDPossible(TaoTerm term, PetscBool3 ispossible)
345: {
346: PetscFunctionBegin;
348: PetscTryMethod(term, "TaoTermShellSetIsComputeHessianFDPossible_C", (TaoTerm, PetscBool3), (term, ispossible));
349: PetscFunctionReturn(PETSC_SUCCESS);
350: }
352: static PetscErrorCode TaoTermShellSetIsComputeHessianFDPossible_Shell(TaoTerm term, PetscBool3 ispossible)
353: {
354: TaoTerm_Shell *shell = (TaoTerm_Shell *)term->data;
356: PetscFunctionBegin;
357: shell->iscomputehessianfdpossible = ispossible;
358: term->ops->iscomputehessianfdpossible = TaoTermIsComputeHessianFDPossible_Shell;
359: PetscFunctionReturn(PETSC_SUCCESS);
360: }
362: /*@C
363: TaoTermShellSetView - Set the view function of a `TAOTERMSHELL`
365: Logically collective
367: Input Parameters:
368: + term - a `TaoTerm` of type `TAOTERMSHELL`
369: - view - a function with the same signature as `TaoTermView()`
371: Calling sequence of `view`:
372: + term - the `TaoTerm`
373: - viewer - a `PetscViewer`
375: Level: intermediate
377: .seealso: [](sec_tao_term), `TaoTerm`, `TAOTERMSHELL`, `TaoTermShellGetContext()`, `TaoTermShellSetContextDestroy()`,
378: `TaoTermShellSetObjective()`,
379: `TaoTermShellSetGradient()`,
380: `TaoTermShellSetObjectiveAndGradient()`,
381: `TaoTermShellSetHessian()`
382: @*/
383: PetscErrorCode TaoTermShellSetView(TaoTerm term, PetscErrorCode (*view)(TaoTerm term, PetscViewer viewer))
384: {
385: PetscFunctionBegin;
387: PetscTryMethod(term, "TaoTermShellSetView_C", (TaoTerm, PetscErrorCode (*)(TaoTerm, PetscViewer)), (term, view));
388: PetscFunctionReturn(PETSC_SUCCESS);
389: }
391: static PetscErrorCode TaoTermShellSetView_Shell(TaoTerm term, PetscErrorCode (*view)(TaoTerm, PetscViewer))
392: {
393: PetscFunctionBegin;
394: term->ops->view = view;
395: PetscFunctionReturn(PETSC_SUCCESS);
396: }
398: /*@C
399: TaoTermShellSetCreateSolutionVec - Set the routine that creates solution vector for a `TaoTerm` of type `TAOTERMSHELL`
401: Logically collective
403: Input Parameters:
404: + term - a `TaoTerm` of type `TAOTERMSHELL`
405: - createsolutionvec - a function with the same signature as `TaoTermCreateSolutionVec()`
407: Calling sequence of `createsolutionvec`:
408: + term - the `TaoTerm`
409: - solution - a solution vector for `term`
411: Level: intermediate
413: .seealso: [](sec_tao_term), `TaoTerm`, `TAOTERMSHELL`, `TaoTermShellGetContext()`, `TaoTermShellSetContextDestroy()`,
414: `TaoTermShellSetCreateHessianMatrices()`
415: @*/
416: PetscErrorCode TaoTermShellSetCreateSolutionVec(TaoTerm term, PetscErrorCode (*createsolutionvec)(TaoTerm term, Vec *solution))
417: {
418: PetscFunctionBegin;
420: PetscTryMethod(term, "TaoTermShellSetCreateSolutionVec_C", (TaoTerm, PetscErrorCode (*)(TaoTerm, Vec *)), (term, createsolutionvec));
421: PetscFunctionReturn(PETSC_SUCCESS);
422: }
424: /*@C
425: TaoTermShellSetCreateParametersVec - Set the routine that creates parameters vector for a `TaoTerm` of type `TAOTERMSHELL`
427: Logically collective
429: Input Parameters:
430: + term - a `TaoTerm` of type `TAOTERMSHELL`
431: - createparametersvec - a function with the same signature as `TaoTermCreateParametersVec()`
433: Calling sequence of `createparametersvec`:
434: + term - the `TaoTerm`
435: - parameters - a parameters vector for `term`
437: Level: intermediate
439: .seealso: [](sec_tao_term), `TaoTerm`, `TAOTERMSHELL`, `TaoTermShellGetContext()`, `TaoTermShellSetContextDestroy()`,
440: `TaoTermShellSetCreateHessianMatrices()`
441: @*/
442: PetscErrorCode TaoTermShellSetCreateParametersVec(TaoTerm term, PetscErrorCode (*createparametersvec)(TaoTerm term, Vec *parameters))
443: {
444: PetscFunctionBegin;
446: PetscTryMethod(term, "TaoTermShellSetCreateParametersVec_C", (TaoTerm, PetscErrorCode (*)(TaoTerm, Vec *)), (term, createparametersvec));
447: PetscFunctionReturn(PETSC_SUCCESS);
448: }
450: static PetscErrorCode TaoTermShellSetCreateSolutionVec_Shell(TaoTerm term, PetscErrorCode (*createsolutionvec)(TaoTerm, Vec *))
451: {
452: PetscFunctionBegin;
453: term->ops->createsolutionvec = createsolutionvec;
454: PetscFunctionReturn(PETSC_SUCCESS);
455: }
457: static PetscErrorCode TaoTermShellSetCreateParametersVec_Shell(TaoTerm term, PetscErrorCode (*createparametersvec)(TaoTerm, Vec *))
458: {
459: PetscFunctionBegin;
460: term->ops->createparametersvec = createparametersvec;
461: PetscFunctionReturn(PETSC_SUCCESS);
462: }
464: /*@C
465: TaoTermShellSetCreateHessianMatrices - Set the routine that creates Hessian matrices for a `TaoTerm` of type `TAOTERMSHELL`
467: Logically collective
469: Input Parameters:
470: + term - a `TaoTerm` of type `TAOTERMSHELL`
471: - createmats - a function with the same signature as `TaoTermCreateHessianMatrices()`
473: Calling sequence of `createmats`:
474: + f - the `TaoTerm`
475: . H - (optional) a matrix of the appropriate type and size for the Hessian of `term`
476: - Hpre - (optional) a matrix of the appropriate type and size for constructing a preconditioner for the Hessian of `term`
478: Level: intermediate
480: .seealso: [](sec_tao_term), `TaoTerm`, `TAOTERMSHELL`, `TaoTermShellGetContext()`, `TaoTermShellSetContextDestroy()`,
481: `TaoTermShellSetCreateSolutionVec()`, `TaoTermShellSetCreateParametersVec()`
482: @*/
483: PetscErrorCode TaoTermShellSetCreateHessianMatrices(TaoTerm term, PetscErrorCode (*createmats)(TaoTerm f, Mat *H, Mat *Hpre))
484: {
485: PetscFunctionBegin;
487: PetscTryMethod(term, "TaoTermShellSetCreateHessianMatrices_C", (TaoTerm, PetscErrorCode (*)(TaoTerm, Mat *, Mat *)), (term, createmats));
488: PetscFunctionReturn(PETSC_SUCCESS);
489: }
491: static PetscErrorCode TaoTermShellSetCreateHessianMatrices_Shell(TaoTerm term, PetscErrorCode (*createhessianmatrices)(TaoTerm, Mat *, Mat *))
492: {
493: PetscFunctionBegin;
494: term->ops->createhessianmatrices = createhessianmatrices;
495: PetscFunctionReturn(PETSC_SUCCESS);
496: }
498: /*MC
499: TAOTERMSHELL - A `TaoTerm` that uses user-provided function callbacks for its operations
501: Level: intermediate
503: .seealso: [](sec_tao_term), `TaoTerm`, `TaoTermShellGetContext()`, `TaoTermShellSetContextDestroy()`, `TaoTermCreateShell()`,
504: `TaoTermShellSetObjective()`,
505: `TaoTermShellSetGradient()`,
506: `TaoTermShellSetObjectiveAndGradient()`,
507: `TaoTermShellSetHessian()`
508: M*/
509: PETSC_INTERN PetscErrorCode TaoTermCreate_Shell(TaoTerm term)
510: {
511: TaoTerm_Shell *shell;
513: PetscFunctionBegin;
514: PetscCall(PetscNew(&shell));
515: term->data = (void *)shell;
517: shell->iscomputehessianfdpossible = PETSC_BOOL3_UNKNOWN;
518: term->ops->iscomputehessianfdpossible = NULL;
519: term->ops->destroy = TaoTermDestroy_Shell;
520: term->ops->view = TaoTermView_Shell;
522: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetContextDestroy_C", TaoTermShellSetContextDestroy_Shell));
523: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetContext_C", TaoTermShellSetContext_Shell));
524: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellGetContext_C", TaoTermShellGetContext_Shell));
525: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetObjective_C", TaoTermShellSetObjective_Shell));
526: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetGradient_C", TaoTermShellSetGradient_Shell));
527: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetObjectiveAndGradient_C", TaoTermShellSetObjectiveAndGradient_Shell));
528: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetHessian_C", TaoTermShellSetHessian_Shell));
529: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetIsComputeHessianFDPossible_C", TaoTermShellSetIsComputeHessianFDPossible_Shell));
530: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetView_C", TaoTermShellSetView_Shell));
531: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetCreateSolutionVec_C", TaoTermShellSetCreateSolutionVec_Shell));
532: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetCreateParametersVec_C", TaoTermShellSetCreateParametersVec_Shell));
533: PetscCall(PetscObjectComposeFunction((PetscObject)term, "TaoTermShellSetCreateHessianMatrices_C", TaoTermShellSetCreateHessianMatrices_Shell));
534: PetscFunctionReturn(PETSC_SUCCESS);
535: }
537: /*@C
538: TaoTermCreateShell - Create a `TaoTerm` of type `TAOTERMSHELL` that is ready to accept user-provided callback operations.
540: Collective
542: Input Parameters:
543: + comm - the MPI communicator for computing the term
544: . ctx - (optional) a context to be used by routines
545: - destroy - (optional) a routine to destroy the context when `term` is destroyed
547: Output Parameter:
548: . term - a `TaoTerm` of type `TAOTERMSHELL`
550: Level: intermediate
552: .seealso: [](sec_tao_term), `TaoTerm`, `TAOTERMSHELL`
553: @*/
554: PetscErrorCode TaoTermCreateShell(MPI_Comm comm, PetscCtx ctx, PetscCtxDestroyFn *destroy, TaoTerm *term)
555: {
556: PetscFunctionBegin;
557: PetscAssertPointer(term, 4);
558: PetscCall(TaoTermCreate(comm, term));
559: PetscCall(TaoTermSetType(*term, TAOTERMSHELL));
560: if (ctx) PetscCall(TaoTermShellSetContext(*term, ctx));
561: if (destroy) PetscCall(TaoTermShellSetContextDestroy(*term, destroy));
562: PetscFunctionReturn(PETSC_SUCCESS);
563: }