Actual source code: dmsnes.c
1: #include <petsc/private/snesimpl.h>
2: #include <petsc/private/dmimpl.h>
4: static PetscErrorCode DMSNESUnsetFunctionContext_DMSNES(DMSNES sdm)
5: {
6: PetscFunctionBegin;
7: PetscCall(PetscObjectCompose((PetscObject)sdm, "function ctx", NULL));
8: sdm->functionctxcontainer = NULL;
9: PetscFunctionReturn(PETSC_SUCCESS);
10: }
12: static PetscErrorCode DMSNESUnsetJacobianContext_DMSNES(DMSNES sdm)
13: {
14: PetscFunctionBegin;
15: PetscCall(PetscObjectCompose((PetscObject)sdm, "jacobian ctx", NULL));
16: sdm->jacobianctxcontainer = NULL;
17: PetscFunctionReturn(PETSC_SUCCESS);
18: }
20: static PetscErrorCode DMSNESDestroy(DMSNES *kdm)
21: {
22: PetscFunctionBegin;
23: if (!*kdm) PetscFunctionReturn(PETSC_SUCCESS);
25: if (--((PetscObject)*kdm)->refct > 0) {
26: *kdm = NULL;
27: PetscFunctionReturn(PETSC_SUCCESS);
28: }
29: PetscCall(DMSNESUnsetFunctionContext_DMSNES(*kdm));
30: PetscCall(DMSNESUnsetJacobianContext_DMSNES(*kdm));
31: PetscTryTypeMethod(*kdm, destroy);
32: PetscCall(PetscHeaderDestroy(kdm));
33: PetscFunctionReturn(PETSC_SUCCESS);
34: }
36: PetscErrorCode DMSNESLoad(DMSNES kdm, PetscViewer viewer)
37: {
38: PetscFunctionBegin;
39: PetscCall(PetscViewerBinaryRead(viewer, &kdm->ops->computefunction, 1, NULL, PETSC_FUNCTION));
40: PetscCall(PetscViewerBinaryRead(viewer, &kdm->ops->computejacobian, 1, NULL, PETSC_FUNCTION));
41: PetscFunctionReturn(PETSC_SUCCESS);
42: }
44: PetscErrorCode DMSNESView(DMSNES kdm, PetscViewer viewer)
45: {
46: PetscBool isascii, isbinary;
48: PetscFunctionBegin;
49: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii));
50: PetscCall(PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary));
51: if (isascii) {
52: #if defined(PETSC_SERIALIZE_FUNCTIONS)
53: const char *fname;
55: PetscCall(PetscFPTFind(kdm->ops->computefunction, &fname));
56: if (fname) PetscCall(PetscViewerASCIIPrintf(viewer, "Function used by SNES: %s\n", fname));
57: PetscCall(PetscFPTFind(kdm->ops->computejacobian, &fname));
58: if (fname) PetscCall(PetscViewerASCIIPrintf(viewer, "Jacobian function used by SNES: %s\n", fname));
59: #endif
60: } else if (isbinary) {
61: struct {
62: SNESFunctionFn *func;
63: } funcstruct;
64: struct {
65: SNESJacobianFn *jac;
66: } jacstruct;
67: funcstruct.func = kdm->ops->computefunction;
68: jacstruct.jac = kdm->ops->computejacobian;
69: PetscCall(PetscViewerBinaryWrite(viewer, &funcstruct, 1, PETSC_FUNCTION));
70: PetscCall(PetscViewerBinaryWrite(viewer, &jacstruct, 1, PETSC_FUNCTION));
71: }
72: PetscFunctionReturn(PETSC_SUCCESS);
73: }
75: static PetscErrorCode DMSNESCreate(MPI_Comm comm, DMSNES *kdm)
76: {
77: PetscFunctionBegin;
78: PetscCall(SNESInitializePackage());
79: PetscCall(PetscHeaderCreate(*kdm, DMSNES_CLASSID, "DMSNES", "DMSNES", "DMSNES", comm, DMSNESDestroy, DMSNESView));
80: PetscFunctionReturn(PETSC_SUCCESS);
81: }
83: /* Attaches the DMSNES to the coarse level.
84: * Under what conditions should we copy versus duplicate?
85: */
86: static PetscErrorCode DMCoarsenHook_DMSNES(DM dm, DM dmc, PetscCtx ctx)
87: {
88: PetscFunctionBegin;
89: PetscCall(DMCopyDMSNES(dm, dmc));
90: PetscFunctionReturn(PETSC_SUCCESS);
91: }
93: /* This could restrict auxiliary information to the coarse level.
94: */
95: static PetscErrorCode DMRestrictHook_DMSNES(DM dm, Mat Restrict, Vec rscale, Mat Inject, DM dmc, PetscCtx ctx)
96: {
97: PetscFunctionBegin;
98: PetscFunctionReturn(PETSC_SUCCESS);
99: }
101: /* Attaches the DMSNES to the subdomain. */
102: static PetscErrorCode DMSubDomainHook_DMSNES(DM dm, DM subdm, PetscCtx ctx)
103: {
104: PetscFunctionBegin;
105: PetscCall(DMCopyDMSNES(dm, subdm));
106: PetscFunctionReturn(PETSC_SUCCESS);
107: }
109: /* This could restrict auxiliary information to the coarse level.
110: */
111: static PetscErrorCode DMSubDomainRestrictHook_DMSNES(DM dm, VecScatter gscat, VecScatter lscat, DM subdm, PetscCtx ctx)
112: {
113: PetscFunctionBegin;
114: PetscFunctionReturn(PETSC_SUCCESS);
115: }
117: static PetscErrorCode DMRefineHook_DMSNES(DM dm, DM dmf, PetscCtx ctx)
118: {
119: PetscFunctionBegin;
120: PetscCall(DMCopyDMSNES(dm, dmf));
121: PetscFunctionReturn(PETSC_SUCCESS);
122: }
124: /* This could restrict auxiliary information to the coarse level.
125: */
126: static PetscErrorCode DMInterpolateHook_DMSNES(DM dm, Mat Interp, DM dmf, PetscCtx ctx)
127: {
128: PetscFunctionBegin;
129: PetscFunctionReturn(PETSC_SUCCESS);
130: }
132: /*
133: DMSNESCopy - copies the information in a `DMSNES` to another `DMSNES`
135: Not Collective
137: Input Parameters:
138: + kdm - Original `DMSNES`
139: - nkdm - `DMSNES` to receive the data, should have been created with `DMSNESCreate()`
141: Level: developer
143: .seealso: [](ch_snes), `DMSNES`, `DMSNESCreate()`, `DMSNESDestroy()`
144: */
145: static PetscErrorCode DMSNESCopy(DMSNES kdm, DMSNES nkdm)
146: {
147: PetscFunctionBegin;
150: nkdm->ops->computefunction = kdm->ops->computefunction;
151: nkdm->ops->computejacobian = kdm->ops->computejacobian;
152: nkdm->ops->computegs = kdm->ops->computegs;
153: nkdm->ops->computeobjective = kdm->ops->computeobjective;
154: nkdm->ops->computepjacobian = kdm->ops->computepjacobian;
155: nkdm->ops->computepfunction = kdm->ops->computepfunction;
156: nkdm->ops->destroy = kdm->ops->destroy;
157: nkdm->ops->duplicate = kdm->ops->duplicate;
159: nkdm->gsctx = kdm->gsctx;
160: nkdm->pctx = kdm->pctx;
161: nkdm->objectivectx = kdm->objectivectx;
162: nkdm->originaldm = kdm->originaldm;
163: nkdm->functionctxcontainer = kdm->functionctxcontainer;
164: nkdm->jacobianctxcontainer = kdm->jacobianctxcontainer;
165: if (nkdm->functionctxcontainer) PetscCall(PetscObjectCompose((PetscObject)nkdm, "function ctx", (PetscObject)nkdm->functionctxcontainer));
166: if (nkdm->jacobianctxcontainer) PetscCall(PetscObjectCompose((PetscObject)nkdm, "jacobian ctx", (PetscObject)nkdm->jacobianctxcontainer));
168: /*
169: nkdm->fortran_func_pointers[0] = kdm->fortran_func_pointers[0];
170: nkdm->fortran_func_pointers[1] = kdm->fortran_func_pointers[1];
171: nkdm->fortran_func_pointers[2] = kdm->fortran_func_pointers[2];
172: */
174: /* implementation specific copy hooks */
175: PetscTryTypeMethod(kdm, duplicate, nkdm);
176: PetscFunctionReturn(PETSC_SUCCESS);
177: }
179: /*@C
180: DMGetDMSNES - get read-only private `DMSNES` context from a `DM`
182: Not Collective
184: Input Parameter:
185: . dm - `DM` to be used with `SNES`
187: Output Parameter:
188: . snesdm - private `DMSNES` context
190: Level: developer
192: Note:
193: Use `DMGetDMSNESWrite()` if write access is needed. The DMSNESSetXXX API should be used wherever possible.
195: .seealso: [](ch_snes), `DMSNES`, `DMGetDMSNESWrite()`
196: @*/
197: PetscErrorCode DMGetDMSNES(DM dm, DMSNES *snesdm)
198: {
199: PetscFunctionBegin;
201: *snesdm = (DMSNES)dm->dmsnes;
202: if (!*snesdm) {
203: PetscCall(PetscInfo(dm, "Creating new DMSNES\n"));
204: PetscCall(DMSNESCreate(PetscObjectComm((PetscObject)dm), snesdm));
206: dm->dmsnes = (PetscObject)*snesdm;
207: (*snesdm)->originaldm = dm;
208: PetscCall(DMCoarsenHookAdd(dm, DMCoarsenHook_DMSNES, DMRestrictHook_DMSNES, NULL));
209: PetscCall(DMRefineHookAdd(dm, DMRefineHook_DMSNES, DMInterpolateHook_DMSNES, NULL));
210: PetscCall(DMSubDomainHookAdd(dm, DMSubDomainHook_DMSNES, DMSubDomainRestrictHook_DMSNES, NULL));
211: }
212: PetscFunctionReturn(PETSC_SUCCESS);
213: }
215: /*@C
216: DMGetDMSNESWrite - get write access to private `DMSNES` context from a `DM`
218: Not Collective
220: Input Parameter:
221: . dm - `DM` to be used with `SNES`
223: Output Parameter:
224: . snesdm - private `DMSNES` context
226: Level: developer
228: .seealso: [](ch_snes), `DMSNES`, `DMGetDMSNES()`
229: @*/
230: PetscErrorCode DMGetDMSNESWrite(DM dm, DMSNES *snesdm)
231: {
232: DMSNES sdm;
234: PetscFunctionBegin;
236: PetscCall(DMGetDMSNES(dm, &sdm));
237: PetscCheck(sdm->originaldm, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "DMSNES has a NULL originaldm");
238: if (sdm->originaldm != dm) { /* Copy on write */
239: DMSNES oldsdm = sdm;
240: PetscCall(PetscInfo(dm, "Copying DMSNES due to write\n"));
241: PetscCall(DMSNESCreate(PetscObjectComm((PetscObject)dm), &sdm));
242: PetscCall(DMSNESCopy(oldsdm, sdm));
243: PetscCall(DMSNESDestroy((DMSNES *)&dm->dmsnes));
244: dm->dmsnes = (PetscObject)sdm;
245: sdm->originaldm = dm;
246: }
247: *snesdm = sdm;
248: PetscFunctionReturn(PETSC_SUCCESS);
249: }
251: /*@
252: DMCopyDMSNES - copies a `DMSNES` context to a new `DM`
254: Logically Collective
256: Input Parameters:
257: + dmsrc - `DM` to obtain context from
258: - dmdest - `DM` to add context to
260: Level: developer
262: Note:
263: The context is copied by reference. This function does not ensure that a context exists.
265: .seealso: [](ch_snes), `DMSNES`, `DMGetDMSNES()`, `SNESSetDM()`
266: @*/
267: PetscErrorCode DMCopyDMSNES(DM dmsrc, DM dmdest)
268: {
269: PetscFunctionBegin;
272: if (!dmdest->dmsnes) PetscCall(DMSNESCreate(PetscObjectComm((PetscObject)dmdest), (DMSNES *)&dmdest->dmsnes));
273: PetscCall(DMSNESCopy((DMSNES)dmsrc->dmsnes, (DMSNES)dmdest->dmsnes));
274: PetscCall(DMCoarsenHookAdd(dmdest, DMCoarsenHook_DMSNES, NULL, NULL));
275: PetscCall(DMRefineHookAdd(dmdest, DMRefineHook_DMSNES, NULL, NULL));
276: PetscCall(DMSubDomainHookAdd(dmdest, DMSubDomainHook_DMSNES, DMSubDomainRestrictHook_DMSNES, NULL));
277: PetscFunctionReturn(PETSC_SUCCESS);
278: }
280: /*@C
281: DMSNESSetFunction - set `SNES` residual evaluation function
283: Not Collective
285: Input Parameters:
286: + dm - DM to be used with `SNES`
287: . f - residual evaluation function; see `SNESFunctionFn` for calling sequence
288: - ctx - context for residual evaluation
290: Level: developer
292: Note:
293: `SNESSetFunction()` is normally used, but it calls this function internally because the user context is actually
294: associated with the `DM`. This makes the interface consistent regardless of whether the user interacts with a `DM` or
295: not.
296: If both `f` and `ctx` are `NULL`, it removes the callback.
298: Developer Note:
299: If `DM` took a more central role at some later date, this could become the primary method of setting the residual.
301: .seealso: [](ch_snes), `DMSNES`, `DMSNESSetContext()`, `SNESSetFunction()`, `DMSNESSetJacobian()`, `SNESFunctionFn`
302: @*/
303: PetscErrorCode DMSNESSetFunction(DM dm, SNESFunctionFn *f, PetscCtx ctx)
304: {
305: DMSNES sdm;
307: PetscFunctionBegin;
309: PetscCall(DMGetDMSNESWrite(dm, &sdm));
310: if (f) sdm->ops->computefunction = f;
311: if (ctx) {
312: PetscContainer ctxcontainer;
313: PetscCall(PetscContainerCreate(PetscObjectComm((PetscObject)sdm), &ctxcontainer));
314: PetscCall(PetscContainerSetPointer(ctxcontainer, ctx));
315: PetscCall(PetscObjectCompose((PetscObject)sdm, "function ctx", (PetscObject)ctxcontainer));
316: sdm->functionctxcontainer = ctxcontainer;
317: PetscCall(PetscContainerDestroy(&ctxcontainer));
318: }
319: if (!f && !ctx) {
320: sdm->ops->computefunction = NULL;
321: PetscCall(DMSNESUnsetFunctionContext_DMSNES(sdm));
322: }
323: PetscFunctionReturn(PETSC_SUCCESS);
324: }
326: /*@C
327: DMSNESSetFunctionContextDestroy - set `SNES` residual evaluation context destroy function
329: Not Collective
331: Input Parameters:
332: + dm - `DM` to be used with `SNES`
333: - f - residual evaluation context destroy function, see `PetscCtxDestroyFn` for its calling sequence
335: Level: developer
337: .seealso: [](ch_snes), `DMSNES`, `DMSNESSetFunction()`, `SNESSetFunction()`, `PetscCtxDestroyFn`
338: @*/
339: PetscErrorCode DMSNESSetFunctionContextDestroy(DM dm, PetscCtxDestroyFn *f)
340: {
341: DMSNES sdm;
343: PetscFunctionBegin;
345: PetscCall(DMGetDMSNESWrite(dm, &sdm));
346: if (sdm->functionctxcontainer) PetscCall(PetscContainerSetCtxDestroy(sdm->functionctxcontainer, f));
347: PetscFunctionReturn(PETSC_SUCCESS);
348: }
350: PetscErrorCode DMSNESUnsetFunctionContext_Internal(DM dm)
351: {
352: DMSNES sdm;
354: PetscFunctionBegin;
356: PetscCall(DMGetDMSNESWrite(dm, &sdm));
357: PetscCall(DMSNESUnsetFunctionContext_DMSNES(sdm));
358: PetscFunctionReturn(PETSC_SUCCESS);
359: }
361: /*@C
362: DMSNESSetMFFunction - set `SNES` residual evaluation function used in applying the matrix-free Jacobian with `-snes_mf_operator`
364: Logically Collective
366: Input Parameters:
367: + dm - `DM` to be used with `SNES`
368: . func - residual evaluation function; see `SNESFunctionFn` for calling sequence
369: - ctx - optional function context
371: Level: developer
373: Note:
374: If not provided then the function provided with `SNESSetFunction()` is used.
375: If both `func` and `ctx` are `NULL`, it removes the callback.
377: .seealso: [](ch_snes), `DMSNES`, `DMSNESSetContext()`, `SNESSetFunction()`, `DMSNESSetJacobian()`, `DMSNESSetFunction()`, `SNESFunctionFn`
378: @*/
379: PetscErrorCode DMSNESSetMFFunction(DM dm, SNESFunctionFn *func, PetscCtx ctx)
380: {
381: DMSNES sdm;
383: PetscFunctionBegin;
385: PetscCall(DMGetDMSNESWrite(dm, &sdm));
386: if (func) sdm->ops->computemffunction = func;
387: if (ctx) sdm->mffunctionctx = ctx;
388: if (!func && !ctx) {
389: sdm->ops->computemffunction = NULL;
390: sdm->mffunctionctx = NULL;
391: }
392: PetscFunctionReturn(PETSC_SUCCESS);
393: }
395: /*@C
396: DMSNESGetFunction - get `SNES` residual evaluation function from a `DMSNES` object
398: Not Collective
400: Input Parameter:
401: . dm - `DM` to be used with `SNES`
403: Output Parameters:
404: + f - residual evaluation function; see `SNESFunctionFn` for calling sequence
405: - ctx - context for residual evaluation
407: Level: developer
409: Note:
410: `SNESGetFunction()` is normally used, but it calls this function internally because the user context is actually
411: associated with the `DM`.
413: .seealso: [](ch_snes), `DMSNES`, `DMSNESSetContext()`, `DMSNESSetFunction()`, `SNESSetFunction()`, `SNESFunctionFn`
414: @*/
415: PetscErrorCode DMSNESGetFunction(DM dm, SNESFunctionFn **f, PetscCtxRt ctx)
416: {
417: DMSNES sdm;
419: PetscFunctionBegin;
421: PetscCall(DMGetDMSNES(dm, &sdm));
422: if (f) *f = sdm->ops->computefunction;
423: if (ctx) {
424: if (sdm->functionctxcontainer) PetscCall(PetscContainerGetPointer(sdm->functionctxcontainer, ctx));
425: else *(void **)ctx = NULL;
426: }
427: PetscFunctionReturn(PETSC_SUCCESS);
428: }
430: /*@C
431: DMSNESSetObjective - Sets the objective function minimized by some of the `SNES` linesearch methods into a `DMSNES` object, used instead of the 2-norm of the residual
433: Not Collective
435: Input Parameters:
436: + dm - `DM` to be used with `SNES`
437: . obj - objective evaluation routine; see `SNESObjectiveFn` for the calling sequence
438: - ctx - [optional] user-defined context for private data for the objective evaluation routine (may be `NULL`)
440: Level: developer
442: Note:
443: If both `obj` and `ctx` are `NULL`, it removes the callback.
445: .seealso: [](ch_snes), `DMSNES`, `DMSNESSetContext()`, `SNESGetObjective()`, `DMSNESSetFunction()`, `SNESObjectiveFn`
446: @*/
447: PetscErrorCode DMSNESSetObjective(DM dm, SNESObjectiveFn *obj, PetscCtx ctx)
448: {
449: DMSNES sdm;
451: PetscFunctionBegin;
453: PetscCall(DMGetDMSNESWrite(dm, &sdm));
454: if (obj) sdm->ops->computeobjective = obj;
455: if (ctx) sdm->objectivectx = ctx;
456: if (!obj && !ctx) {
457: sdm->ops->computeobjective = NULL;
458: sdm->objectivectx = NULL;
459: }
460: PetscFunctionReturn(PETSC_SUCCESS);
461: }
463: /*@C
464: DMSNESGetObjective - Returns the objective function set with `DMSNESSetObjective()`
466: Not Collective
468: Input Parameter:
469: . dm - `DM` to be used with `SNES`
471: Output Parameters:
472: + obj - objective evaluation routine (or `NULL`); see `SNESObjectiveFn` for the calling sequence
473: - ctx - the function context (or `NULL`)
475: Level: developer
477: Note:
478: `SNESGetFunction()` is normally used, but it calls this function internally because the user context is actually
479: associated with the `DM`.
481: .seealso: [](ch_snes), `DMSNES`, `DMSNESSetContext()`, `DMSNESSetObjective()`, `SNESSetFunction()`, `SNESObjectiveFn`
482: @*/
483: PetscErrorCode DMSNESGetObjective(DM dm, SNESObjectiveFn **obj, PetscCtxRt ctx)
484: {
485: DMSNES sdm;
487: PetscFunctionBegin;
489: PetscCall(DMGetDMSNES(dm, &sdm));
490: if (obj) *obj = sdm->ops->computeobjective;
491: if (ctx) *(void **)ctx = sdm->objectivectx;
492: PetscFunctionReturn(PETSC_SUCCESS);
493: }
495: /*@C
496: DMSNESSetNGS - set `SNES` Gauss-Seidel relaxation function into a `DMSNES` object
498: Not Collective
500: Input Parameters:
501: + dm - `DM` to be used with `SNES`
502: . f - relaxation function, see `SNESGSFunction`
503: - ctx - context for residual evaluation
505: Level: developer
507: Note:
508: `SNESSetNGS()` is normally used, but it calls this function internally because the user context is actually
509: associated with the `DM`. This makes the interface consistent regardless of whether the user interacts with a `DM` or
510: not.
511: If both `f` and `ctx` are `NULL`, it removes the callback.
513: Developer Note:
514: If `DM` took a more central role at some later date, this could become the primary method of supplying the smoother
516: .seealso: [](ch_snes), `DMSNES`, `DMSNESSetContext()`, `SNESSetFunction()`, `DMSNESSetJacobian()`, `DMSNESSetFunction()`, `SNESGSFunction`
517: @*/
518: PetscErrorCode DMSNESSetNGS(DM dm, PetscErrorCode (*f)(SNES, Vec, Vec, void *), PetscCtx ctx)
519: {
520: DMSNES sdm;
522: PetscFunctionBegin;
524: PetscCall(DMGetDMSNESWrite(dm, &sdm));
525: if (f) sdm->ops->computegs = f;
526: if (ctx) sdm->gsctx = ctx;
527: if (!f && !ctx) {
528: sdm->ops->computegs = NULL;
529: sdm->gsctx = NULL;
530: }
531: PetscFunctionReturn(PETSC_SUCCESS);
532: }
534: /*@C
535: DMSNESGetNGS - get `SNES` Gauss-Seidel relaxation function from a `DMSNES` object
537: Not Collective
539: Input Parameter:
540: . dm - `DM` to be used with `SNES`
542: Output Parameters:
543: + f - relaxation function which performs Gauss-Seidel sweeps, see `SNESSetNGS()`
544: - ctx - context for residual evaluation
546: Level: developer
548: Note:
549: `SNESGetNGS()` is normally used, but it calls this function internally because the user context is actually
550: associated with the `DM`.
552: Developer Note:
553: This makes the interface consistent regardless of whether the user interacts with a `DM` or
554: not. If `DM` took a more central role at some later date, this could become the primary method of setting the residual.
556: .seealso: [](ch_snes), `DMSNES`, `DMSNESSetContext()`, `SNESGetNGS()`, `DMSNESGetJacobian()`, `DMSNESGetFunction()`
557: @*/
558: PetscErrorCode DMSNESGetNGS(DM dm, PetscErrorCode (**f)(SNES, Vec, Vec, void *), PetscCtxRt ctx)
559: {
560: DMSNES sdm;
562: PetscFunctionBegin;
564: PetscCall(DMGetDMSNES(dm, &sdm));
565: if (f) *f = sdm->ops->computegs;
566: if (ctx) *(void **)ctx = sdm->gsctx;
567: PetscFunctionReturn(PETSC_SUCCESS);
568: }
570: /*@C
571: DMSNESSetJacobian - set `SNES` Jacobian evaluation function into a `DMSNES` object
573: Not Collective
575: Input Parameters:
576: + dm - `DM` to be used with `SNES`
577: . J - Jacobian evaluation function, see `SNESJacobianFn`
578: - ctx - context for Jacobian evaluation
580: Level: developer
582: Note:
583: `SNESSetJacobian()` is normally used, but it calls this function internally because the user context is actually
584: associated with the `DM`.
585: If both `J` and `ctx` are `NULL`, it removes the callback.
587: Developer Note:
588: This makes the interface consistent regardless of whether the user interacts with a `DM` or
589: not. If `DM` took a more central role at some later date, this could become the primary method of setting the Jacobian.
591: .seealso: [](ch_snes), `DMSNES`, `DMSNESSetContext()`, `SNESSetFunction()`, `DMSNESGetJacobian()`, `SNESSetJacobian()`, `SNESJacobianFn`
592: @*/
593: PetscErrorCode DMSNESSetJacobian(DM dm, SNESJacobianFn *J, PetscCtx ctx)
594: {
595: DMSNES sdm;
597: PetscFunctionBegin;
599: PetscCall(DMGetDMSNESWrite(dm, &sdm));
600: if (J) sdm->ops->computejacobian = J;
601: if (ctx) {
602: PetscContainer ctxcontainer;
603: PetscCall(PetscContainerCreate(PetscObjectComm((PetscObject)sdm), &ctxcontainer));
604: PetscCall(PetscContainerSetPointer(ctxcontainer, ctx));
605: PetscCall(PetscObjectCompose((PetscObject)sdm, "jacobian ctx", (PetscObject)ctxcontainer));
606: sdm->jacobianctxcontainer = ctxcontainer;
607: PetscCall(PetscContainerDestroy(&ctxcontainer));
608: }
609: if (!J && !ctx) {
610: sdm->ops->computejacobian = NULL;
611: PetscCall(DMSNESUnsetJacobianContext_DMSNES(sdm));
612: }
613: PetscFunctionReturn(PETSC_SUCCESS);
614: }
616: /*@C
617: DMSNESSetJacobianContextDestroy - set `SNES` Jacobian evaluation context destroy function into a `DMSNES` object
619: Not Collective
621: Input Parameters:
622: + dm - `DM` to be used with `SNES`
623: - f - Jacobian evaluation context destroy function, see `PetscCtxDestroyFn` for its calling sequence
625: Level: developer
627: .seealso: [](ch_snes), `DMSNES`, `DMSNESSetJacobian()`
628: @*/
629: PetscErrorCode DMSNESSetJacobianContextDestroy(DM dm, PetscCtxDestroyFn *f)
630: {
631: DMSNES sdm;
633: PetscFunctionBegin;
635: PetscCall(DMGetDMSNESWrite(dm, &sdm));
636: if (sdm->jacobianctxcontainer) PetscCall(PetscContainerSetCtxDestroy(sdm->jacobianctxcontainer, f));
637: PetscFunctionReturn(PETSC_SUCCESS);
638: }
640: PetscErrorCode DMSNESUnsetJacobianContext_Internal(DM dm)
641: {
642: DMSNES sdm;
644: PetscFunctionBegin;
646: PetscCall(DMGetDMSNESWrite(dm, &sdm));
647: PetscCall(DMSNESUnsetJacobianContext_DMSNES(sdm));
648: PetscFunctionReturn(PETSC_SUCCESS);
649: }
651: /*@C
652: DMSNESGetJacobian - get `SNES` Jacobian evaluation function from a `DMSNES` object
654: Not Collective
656: Input Parameter:
657: . dm - `DM` to be used with `SNES`
659: Output Parameters:
660: + J - Jacobian evaluation function; for all calling sequence see `SNESJacobianFn`
661: - ctx - context for residual evaluation
663: Level: developer
665: Note:
666: `SNESGetJacobian()` is normally used, but it calls this function internally because the user context is actually
667: associated with the `DM`. This makes the interface consistent regardless of whether the user interacts with a `DM` or
668: not.
670: Developer Note:
671: If `DM` took a more central role at some later date, this could become the primary method of setting the Jacobian.
673: .seealso: [](ch_snes), `DMSNES`, `DMSNESSetContext()`, `SNESSetFunction()`, `DMSNESSetJacobian()`, `SNESJacobianFn`
674: @*/
675: PetscErrorCode DMSNESGetJacobian(DM dm, SNESJacobianFn **J, PetscCtxRt ctx)
676: {
677: DMSNES sdm;
679: PetscFunctionBegin;
681: PetscCall(DMGetDMSNES(dm, &sdm));
682: if (J) *J = sdm->ops->computejacobian;
683: if (ctx) {
684: if (sdm->jacobianctxcontainer) PetscCall(PetscContainerGetPointer(sdm->jacobianctxcontainer, ctx));
685: else *(void **)ctx = NULL;
686: }
687: PetscFunctionReturn(PETSC_SUCCESS);
688: }
690: /*@C
691: DMSNESSetPicard - set SNES Picard iteration matrix and RHS evaluation functions into a `DMSNES` object
693: Not Collective
695: Input Parameters:
696: + dm - `DM` to be used with `SNES`
697: . b - RHS evaluation function; see `SNESFunctionFn` for calling sequence
698: . J - Picard matrix evaluation function; see `SNESJacobianFn` for calling sequence
699: - ctx - context for residual and matrix evaluation
701: Level: developer
703: .seealso: [](ch_snes), `DMSNES`, `SNESSetPicard()`, `DMSNESSetFunction()`, `DMSNESSetJacobian()`, `SNESFunctionFn`, `SNESJacobianFn`
704: @*/
705: PetscErrorCode DMSNESSetPicard(DM dm, SNESFunctionFn *b, SNESJacobianFn *J, PetscCtx ctx)
706: {
707: DMSNES sdm;
709: PetscFunctionBegin;
711: PetscCall(DMGetDMSNES(dm, &sdm));
712: if (b) sdm->ops->computepfunction = b;
713: if (J) sdm->ops->computepjacobian = J;
714: if (ctx) sdm->pctx = ctx;
715: PetscFunctionReturn(PETSC_SUCCESS);
716: }
718: /*@C
719: DMSNESGetPicard - get `SNES` Picard iteration evaluation functions from a `DMSNES` object
721: Not Collective
723: Input Parameter:
724: . dm - `DM` to be used with `SNES`
726: Output Parameters:
727: + b - RHS evaluation function; see `SNESFunctionFn` for calling sequence
728: . J - Jacobian evaluation function; see `SNESJacobianFn` for calling sequence
729: - ctx - context for residual and matrix evaluation
731: Level: developer
733: .seealso: [](ch_snes), `DMSNES`, `DMSNESSetContext()`, `SNESSetFunction()`, `DMSNESSetJacobian()`, `SNESFunctionFn`, `SNESJacobianFn`
734: @*/
735: PetscErrorCode DMSNESGetPicard(DM dm, SNESFunctionFn **b, SNESJacobianFn **J, PetscCtxRt ctx)
736: {
737: DMSNES sdm;
739: PetscFunctionBegin;
741: PetscCall(DMGetDMSNES(dm, &sdm));
742: if (b) *b = sdm->ops->computepfunction;
743: if (J) *J = sdm->ops->computepjacobian;
744: if (ctx) *(void **)ctx = sdm->pctx;
745: PetscFunctionReturn(PETSC_SUCCESS);
746: }