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: PetscObjectCompose((PetscObject)sdm, "function ctx", NULL);
7: sdm->functionctxcontainer = NULL;
8: return 0;
9: }
11: static PetscErrorCode DMSNESUnsetJacobianContext_DMSNES(DMSNES sdm)
12: {
13: PetscObjectCompose((PetscObject)sdm, "jacobian ctx", NULL);
14: sdm->jacobianctxcontainer = NULL;
15: return 0;
16: }
18: static PetscErrorCode DMSNESDestroy(DMSNES *kdm)
19: {
20: if (!*kdm) return 0;
22: if (--((PetscObject)(*kdm))->refct > 0) {
23: *kdm = NULL;
24: return 0;
25: }
26: DMSNESUnsetFunctionContext_DMSNES(*kdm);
27: DMSNESUnsetJacobianContext_DMSNES(*kdm);
28: if ((*kdm)->ops->destroy) ((*kdm)->ops->destroy)(*kdm);
29: PetscHeaderDestroy(kdm);
30: return 0;
31: }
33: PetscErrorCode DMSNESLoad(DMSNES kdm, PetscViewer viewer)
34: {
35: PetscViewerBinaryRead(viewer, &kdm->ops->computefunction, 1, NULL, PETSC_FUNCTION);
36: PetscViewerBinaryRead(viewer, &kdm->ops->computejacobian, 1, NULL, PETSC_FUNCTION);
37: return 0;
38: }
40: PetscErrorCode DMSNESView(DMSNES kdm, PetscViewer viewer)
41: {
42: PetscBool isascii, isbinary;
44: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERASCII, &isascii);
45: PetscObjectTypeCompare((PetscObject)viewer, PETSCVIEWERBINARY, &isbinary);
46: if (isascii) {
47: #if defined(PETSC_SERIALIZE_FUNCTIONS) && defined(PETSC_SERIALIZE_FUNCTIONS_VIEW)
48: const char *fname;
50: PetscFPTFind(kdm->ops->computefunction, &fname);
51: if (fname) PetscViewerASCIIPrintf(viewer, "Function used by SNES: %s\n", fname);
52: PetscFPTFind(kdm->ops->computejacobian, &fname);
53: if (fname) PetscViewerASCIIPrintf(viewer, "Jacobian function used by SNES: %s\n", fname);
54: #endif
55: } else if (isbinary) {
56: struct {
57: PetscErrorCode (*func)(SNES, Vec, Vec, void *);
58: } funcstruct;
59: struct {
60: PetscErrorCode (*jac)(SNES, Vec, Mat, Mat, void *);
61: } jacstruct;
62: funcstruct.func = kdm->ops->computefunction;
63: jacstruct.jac = kdm->ops->computejacobian;
64: PetscViewerBinaryWrite(viewer, &funcstruct, 1, PETSC_FUNCTION);
65: PetscViewerBinaryWrite(viewer, &jacstruct, 1, PETSC_FUNCTION);
66: }
67: return 0;
68: }
70: static PetscErrorCode DMSNESCreate(MPI_Comm comm, DMSNES *kdm)
71: {
72: SNESInitializePackage();
73: PetscHeaderCreate(*kdm, DMSNES_CLASSID, "DMSNES", "DMSNES", "DMSNES", comm, DMSNESDestroy, DMSNESView);
74: return 0;
75: }
77: /* Attaches the DMSNES to the coarse level.
78: * Under what conditions should we copy versus duplicate?
79: */
80: static PetscErrorCode DMCoarsenHook_DMSNES(DM dm, DM dmc, void *ctx)
81: {
82: DMCopyDMSNES(dm, dmc);
83: return 0;
84: }
86: /* This could restrict auxiliary information to the coarse level.
87: */
88: static PetscErrorCode DMRestrictHook_DMSNES(DM dm, Mat Restrict, Vec rscale, Mat Inject, DM dmc, void *ctx)
89: {
90: return 0;
91: }
93: /* Attaches the DMSNES to the subdomain. */
94: static PetscErrorCode DMSubDomainHook_DMSNES(DM dm, DM subdm, void *ctx)
95: {
96: DMCopyDMSNES(dm, subdm);
97: return 0;
98: }
100: /* This could restrict auxiliary information to the coarse level.
101: */
102: static PetscErrorCode DMSubDomainRestrictHook_DMSNES(DM dm, VecScatter gscat, VecScatter lscat, DM subdm, void *ctx)
103: {
104: return 0;
105: }
107: static PetscErrorCode DMRefineHook_DMSNES(DM dm, DM dmf, void *ctx)
108: {
109: DMCopyDMSNES(dm, dmf);
110: return 0;
111: }
113: /* This could restrict auxiliary information to the coarse level.
114: */
115: static PetscErrorCode DMInterpolateHook_DMSNES(DM dm, Mat Interp, DM dmf, void *ctx)
116: {
117: return 0;
118: }
120: /*@C
121: DMSNESCopy - copies the information in a `DMSNES` to another `DMSNES`
123: Not Collective
125: Input Parameters:
126: + kdm - Original `DMSNES`
127: - nkdm - `DMSNES` to receive the data, should have been created with `DMSNESCreate()`
129: Level: developer
131: .seealso: `DMSNES`, `DMSNESCreate()`, `DMSNESDestroy()`
132: @*/
133: PetscErrorCode DMSNESCopy(DMSNES kdm, DMSNES nkdm)
134: {
137: nkdm->ops->computefunction = kdm->ops->computefunction;
138: nkdm->ops->computejacobian = kdm->ops->computejacobian;
139: nkdm->ops->computegs = kdm->ops->computegs;
140: nkdm->ops->computeobjective = kdm->ops->computeobjective;
141: nkdm->ops->computepjacobian = kdm->ops->computepjacobian;
142: nkdm->ops->computepfunction = kdm->ops->computepfunction;
143: nkdm->ops->destroy = kdm->ops->destroy;
144: nkdm->ops->duplicate = kdm->ops->duplicate;
146: nkdm->gsctx = kdm->gsctx;
147: nkdm->pctx = kdm->pctx;
148: nkdm->objectivectx = kdm->objectivectx;
149: nkdm->originaldm = kdm->originaldm;
150: nkdm->functionctxcontainer = kdm->functionctxcontainer;
151: nkdm->jacobianctxcontainer = kdm->jacobianctxcontainer;
152: if (nkdm->functionctxcontainer) PetscObjectCompose((PetscObject)nkdm, "function ctx", (PetscObject)nkdm->functionctxcontainer);
153: if (nkdm->jacobianctxcontainer) PetscObjectCompose((PetscObject)nkdm, "jacobian ctx", (PetscObject)nkdm->jacobianctxcontainer);
155: /*
156: nkdm->fortran_func_pointers[0] = kdm->fortran_func_pointers[0];
157: nkdm->fortran_func_pointers[1] = kdm->fortran_func_pointers[1];
158: nkdm->fortran_func_pointers[2] = kdm->fortran_func_pointers[2];
159: */
161: /* implementation specific copy hooks */
162: PetscTryTypeMethod(kdm, duplicate, nkdm);
163: return 0;
164: }
166: /*@C
167: DMGetDMSNES - get read-only private `DMSNES` context from a `DM`
169: Not Collective
171: Input Parameter:
172: . dm - `DM` to be used with `SNES`
174: Output Parameter:
175: . snesdm - private `DMSNES` context
177: Level: developer
179: Note:
180: Use `DMGetDMSNESWrite()` if write access is needed. The DMSNESSetXXX API should be used wherever possible.
182: .seealso: `DMGetDMSNESWrite()`
183: @*/
184: PetscErrorCode DMGetDMSNES(DM dm, DMSNES *snesdm)
185: {
187: *snesdm = (DMSNES)dm->dmsnes;
188: if (!*snesdm) {
189: PetscInfo(dm, "Creating new DMSNES\n");
190: DMSNESCreate(PetscObjectComm((PetscObject)dm), snesdm);
192: dm->dmsnes = (PetscObject)*snesdm;
193: (*snesdm)->originaldm = dm;
194: DMCoarsenHookAdd(dm, DMCoarsenHook_DMSNES, DMRestrictHook_DMSNES, NULL);
195: DMRefineHookAdd(dm, DMRefineHook_DMSNES, DMInterpolateHook_DMSNES, NULL);
196: DMSubDomainHookAdd(dm, DMSubDomainHook_DMSNES, DMSubDomainRestrictHook_DMSNES, NULL);
197: }
198: return 0;
199: }
201: /*@C
202: DMGetDMSNESWrite - get write access to private `DMSNES` context from a `DM`
204: Not Collective
206: Input Parameter:
207: . dm - `DM` to be used with `SNES`
209: Output Parameter:
210: . snesdm - private `DMSNES` context
212: Level: developer
214: .seealso: `DMGetDMSNES()`
215: @*/
216: PetscErrorCode DMGetDMSNESWrite(DM dm, DMSNES *snesdm)
217: {
218: DMSNES sdm;
221: DMGetDMSNES(dm, &sdm);
223: if (sdm->originaldm != dm) { /* Copy on write */
224: DMSNES oldsdm = sdm;
225: PetscInfo(dm, "Copying DMSNES due to write\n");
226: DMSNESCreate(PetscObjectComm((PetscObject)dm), &sdm);
227: DMSNESCopy(oldsdm, sdm);
228: DMSNESDestroy((DMSNES *)&dm->dmsnes);
229: dm->dmsnes = (PetscObject)sdm;
230: sdm->originaldm = dm;
231: }
232: *snesdm = sdm;
233: return 0;
234: }
236: /*@C
237: DMCopyDMSNES - copies a `DMSNES` context to a new `DM`
239: Logically Collective
241: Input Parameters:
242: + dmsrc - `DM` to obtain context from
243: - dmdest - `DM` to add context to
245: Level: developer
247: Note:
248: The context is copied by reference. This function does not ensure that a context exists.
250: .seealso: `DMGetDMSNES()`, `SNESSetDM()`
251: @*/
252: PetscErrorCode DMCopyDMSNES(DM dmsrc, DM dmdest)
253: {
256: if (!dmdest->dmsnes) DMSNESCreate(PetscObjectComm((PetscObject)dmdest), (DMSNES *)&dmdest->dmsnes);
257: DMSNESCopy((DMSNES)dmsrc->dmsnes, (DMSNES)dmdest->dmsnes);
258: DMCoarsenHookAdd(dmdest, DMCoarsenHook_DMSNES, NULL, NULL);
259: DMRefineHookAdd(dmdest, DMRefineHook_DMSNES, NULL, NULL);
260: DMSubDomainHookAdd(dmdest, DMSubDomainHook_DMSNES, DMSubDomainRestrictHook_DMSNES, NULL);
261: return 0;
262: }
264: /*@C
265: DMSNESSetFunction - set `SNES` residual evaluation function
267: Not Collective
269: Input Parameters:
270: + dm - DM to be used with `SNES`
271: . f - residual evaluation function; see `SNESFunction` for details
272: - ctx - context for residual evaluation
274: Level: advanced
276: Note:
277: `SNESSetFunction()` is normally used, but it calls this function internally because the user context is actually
278: associated with the `DM`. This makes the interface consistent regardless of whether the user interacts with a `DM` or
279: not.
281: Developer Note:
282: If `DM` took a more central role at some later date, this could become the primary method of setting the residual.
284: .seealso: `DMSNESSetContext()`, `SNESSetFunction()`, `DMSNESSetJacobian()`, `SNESFunction`
285: @*/
286: PetscErrorCode DMSNESSetFunction(DM dm, PetscErrorCode (*f)(SNES, Vec, Vec, void *), void *ctx)
287: {
288: DMSNES sdm;
291: DMGetDMSNESWrite(dm, &sdm);
292: if (f) sdm->ops->computefunction = f;
293: if (ctx) {
294: PetscContainer ctxcontainer;
295: PetscContainerCreate(PetscObjectComm((PetscObject)sdm), &ctxcontainer);
296: PetscContainerSetPointer(ctxcontainer, ctx);
297: PetscObjectCompose((PetscObject)sdm, "function ctx", (PetscObject)ctxcontainer);
298: sdm->functionctxcontainer = ctxcontainer;
299: PetscContainerDestroy(&ctxcontainer);
300: }
301: return 0;
302: }
304: /*@C
305: DMSNESSetFunctionContextDestroy - set `SNES` residual evaluation context destroy function
307: Not Collective
309: Input Parameters:
310: + dm - `DM` to be used with `SNES`
311: - f - residual evaluation context destroy function
313: Level: advanced
315: .seealso: `DMSNESSetFunction()`, `SNESSetFunction()`
316: @*/
317: PetscErrorCode DMSNESSetFunctionContextDestroy(DM dm, PetscErrorCode (*f)(void *))
318: {
319: DMSNES sdm;
322: DMGetDMSNESWrite(dm, &sdm);
323: if (sdm->functionctxcontainer) PetscContainerSetUserDestroy(sdm->functionctxcontainer, f);
324: return 0;
325: }
327: PetscErrorCode DMSNESUnsetFunctionContext_Internal(DM dm)
328: {
329: DMSNES sdm;
332: DMGetDMSNESWrite(dm, &sdm);
333: DMSNESUnsetFunctionContext_DMSNES(sdm);
334: return 0;
335: }
337: /*@C
338: DMSNESSetMFFunction - set `SNES` residual evaluation function used in applying the matrix-free Jacobian with -snes_mf_operator
340: Logically Collective
342: Input Parameters:
343: + dm - `DM` to be used with `SNES`
344: - f - residual evaluation function; see `SNESFunction` for details
346: Level: advanced
348: .seealso: `DMSNESSetContext()`, `SNESSetFunction()`, `DMSNESSetJacobian()`, `SNESFunction`, `DMSNESSetFunction()`
349: @*/
350: PetscErrorCode DMSNESSetMFFunction(DM dm, PetscErrorCode (*f)(SNES, Vec, Vec, void *), void *ctx)
351: {
352: DMSNES sdm;
355: if (f || ctx) DMGetDMSNESWrite(dm, &sdm);
356: if (f) sdm->ops->computemffunction = f;
357: if (ctx) sdm->mffunctionctx = ctx;
358: return 0;
359: }
361: /*@C
362: DMSNESGetFunction - get `SNES` residual evaluation function
364: Not Collective
366: Input Parameter:
367: . dm - `DM` to be used with `SNES`
369: Output Parameters:
370: + f - residual evaluation function; see `SNESFunction` for details
371: - ctx - context for residual evaluation
373: Level: advanced
375: Note:
376: `SNESGetFunction()` is normally used, but it calls this function internally because the user context is actually
377: associated with the `DM`.
379: .seealso: `DMSNESSetContext()`, `DMSNESSetFunction()`, `SNESSetFunction()`, `SNESFunction`
380: @*/
381: PetscErrorCode DMSNESGetFunction(DM dm, PetscErrorCode (**f)(SNES, Vec, Vec, void *), void **ctx)
382: {
383: DMSNES sdm;
386: DMGetDMSNES(dm, &sdm);
387: if (f) *f = sdm->ops->computefunction;
388: if (ctx) {
389: if (sdm->functionctxcontainer) PetscContainerGetPointer(sdm->functionctxcontainer, ctx);
390: else *ctx = NULL;
391: }
392: return 0;
393: }
395: /*@C
396: DMSNESSetObjective - set `SNES` objective evaluation function
398: Not Collective
400: Input Parameters:
401: + dm - `DM` to be used with `SNES`
402: . obj - objective evaluation function; see `SNESObjectiveFunction` for details
403: - ctx - context for residual evaluation
405: Level: advanced
407: .seealso: `DMSNESSetContext()`, `SNESGetObjective()`, `DMSNESSetFunction()`
408: @*/
409: PetscErrorCode DMSNESSetObjective(DM dm, PetscErrorCode (*obj)(SNES, Vec, PetscReal *, void *), void *ctx)
410: {
411: DMSNES sdm;
414: if (obj || ctx) DMGetDMSNESWrite(dm, &sdm);
415: if (obj) sdm->ops->computeobjective = obj;
416: if (ctx) sdm->objectivectx = ctx;
417: return 0;
418: }
420: /*@C
421: DMSNESGetObjective - get `SNES` objective evaluation function
423: Not Collective
425: Input Parameter:
426: . dm - `DM` to be used with `SNES`
428: Output Parameters:
429: + obj- residual evaluation function; see `SNESObjectiveFunction` for details
430: - ctx - context for residual evaluation
432: Level: advanced
434: Note:
435: `SNESGetFunction()` is normally used, but it calls this function internally because the user context is actually
436: associated with the `DM`.
438: .seealso: `DMSNESSetContext()`, `DMSNESSetObjective()`, `SNESSetFunction()`
439: @*/
440: PetscErrorCode DMSNESGetObjective(DM dm, PetscErrorCode (**obj)(SNES, Vec, PetscReal *, void *), void **ctx)
441: {
442: DMSNES sdm;
445: DMGetDMSNES(dm, &sdm);
446: if (obj) *obj = sdm->ops->computeobjective;
447: if (ctx) *ctx = sdm->objectivectx;
448: return 0;
449: }
451: /*@C
452: DMSNESSetNGS - set `SNES` Gauss-Seidel relaxation function
454: Not Collective
456: Input Parameters:
457: + dm - `DM` to be used with `SNES`
458: . f - relaxation function, see `SNESGSFunction`
459: - ctx - context for residual evaluation
461: Level: advanced
463: Note:
464: `SNESSetNGS()` is normally used, but it calls this function internally because the user context is actually
465: associated with the `DM`. This makes the interface consistent regardless of whether the user interacts with a `DM` or
466: not.
468: Dveloper Note:
469: If `DM `took a more central role at some later date, this could become the primary method of supplying the smoother
471: .seealso: `DMSNESSetContext()`, `SNESSetFunction()`, `DMSNESSetJacobian()`, `DMSNESSetFunction()`, `SNESGSFunction`
472: @*/
473: PetscErrorCode DMSNESSetNGS(DM dm, PetscErrorCode (*f)(SNES, Vec, Vec, void *), void *ctx)
474: {
475: DMSNES sdm;
478: if (f || ctx) DMGetDMSNESWrite(dm, &sdm);
479: if (f) sdm->ops->computegs = f;
480: if (ctx) sdm->gsctx = ctx;
481: return 0;
482: }
484: /*@C
485: DMSNESGetNGS - get `SNES` Gauss-Seidel relaxation function
487: Not Collective
489: Input Parameter:
490: . dm - `DM` to be used with `SNES`
492: Output Parameters:
493: + f - relaxation function which performs Gauss-Seidel sweeps, see `SNESGSFunction`
494: - ctx - context for residual evaluation
496: Level: advanced
498: Note:
499: `SNESGetNGS()` is normally used, but it calls this function internally because the user context is actually
500: associated with the `DM`.
502: Developer Note:
503: This makes the interface consistent regardless of whether the user interacts with a `DM` or
504: not. If `DM` took a more central role at some later date, this could become the primary method of setting the residual.
506: .seealso: `DMSNESSetContext()`, `SNESGetNGS()`, `DMSNESGetJacobian()`, `DMSNESGetFunction()`, `SNESNGSFunction`
507: @*/
508: PetscErrorCode DMSNESGetNGS(DM dm, PetscErrorCode (**f)(SNES, Vec, Vec, void *), void **ctx)
509: {
510: DMSNES sdm;
513: DMGetDMSNES(dm, &sdm);
514: if (f) *f = sdm->ops->computegs;
515: if (ctx) *ctx = sdm->gsctx;
516: return 0;
517: }
519: /*@C
520: DMSNESSetJacobian - set `SNES` Jacobian evaluation function
522: Not Collective
524: Input Parameters:
525: + dm - `DM` to be used with `SNES`
526: . J - Jacobian evaluation function
527: - ctx - context for residual evaluation
529: Level: advanced
531: Note:
532: `SNESSetJacobian()` is normally used, but it calls this function internally because the user context is actually
533: associated with the `DM`.
535: Developer Note:
536: This makes the interface consistent regardless of whether the user interacts with a `DM` or
537: not. If `DM` took a more central role at some later date, this could become the primary method of setting the Jacobian.
539: .seealso: `DMSNESSetContext()`, `SNESSetFunction()`, `DMSNESGetJacobian()`, `SNESSetJacobian()`, `SNESJacobianFunction`
540: @*/
541: PetscErrorCode DMSNESSetJacobian(DM dm, PetscErrorCode (*J)(SNES, Vec, Mat, Mat, void *), void *ctx)
542: {
543: DMSNES sdm;
546: if (J || ctx) DMGetDMSNESWrite(dm, &sdm);
547: if (J) sdm->ops->computejacobian = J;
548: if (ctx) {
549: PetscContainer ctxcontainer;
550: PetscContainerCreate(PetscObjectComm((PetscObject)sdm), &ctxcontainer);
551: PetscContainerSetPointer(ctxcontainer, ctx);
552: PetscObjectCompose((PetscObject)sdm, "jacobian ctx", (PetscObject)ctxcontainer);
553: sdm->jacobianctxcontainer = ctxcontainer;
554: PetscContainerDestroy(&ctxcontainer);
555: }
556: return 0;
557: }
559: /*@C
560: DMSNESSetJacobianContextDestroy - set `SNES` Jacobian evaluation context destroy function
562: Not Collective
564: Input Parameters:
565: + dm - `DM` to be used with `SNES`
566: - f - Jacobian evaluation context destroy function
568: Level: advanced
570: .seealso: `DMSNESSetJacobian()`
571: @*/
572: PetscErrorCode DMSNESSetJacobianContextDestroy(DM dm, PetscErrorCode (*f)(void *))
573: {
574: DMSNES sdm;
577: DMGetDMSNESWrite(dm, &sdm);
578: if (sdm->jacobianctxcontainer) PetscContainerSetUserDestroy(sdm->jacobianctxcontainer, f);
579: return 0;
580: }
582: PetscErrorCode DMSNESUnsetJacobianContext_Internal(DM dm)
583: {
584: DMSNES sdm;
587: DMGetDMSNESWrite(dm, &sdm);
588: DMSNESUnsetJacobianContext_DMSNES(sdm);
589: return 0;
590: }
592: /*@C
593: DMSNESGetJacobian - get `SNES` Jacobian evaluation function
595: Not Collective
597: Input Parameter:
598: . dm - `DM` to be used with `SNES`
600: Output Parameters:
601: + J - Jacobian evaluation function; see `SNESJacobianFunction` for all calling sequence
602: - ctx - context for residual evaluation
604: Level: advanced
606: Note:
607: `SNESGetJacobian()` is normally used, but it calls this function internally because the user context is actually
608: associated with the `DM`. This makes the interface consistent regardless of whether the user interacts with a `DM` or
609: not.
611: Developer Note:
612: If `DM` took a more central role at some later date, this could become the primary method of setting the Jacobian.
614: .seealso: `DMSNESSetContext()`, `SNESSetFunction()`, `DMSNESSetJacobian()`, `SNESJacobianFunction`
615: @*/
616: PetscErrorCode DMSNESGetJacobian(DM dm, PetscErrorCode (**J)(SNES, Vec, Mat, Mat, void *), void **ctx)
617: {
618: DMSNES sdm;
621: DMGetDMSNES(dm, &sdm);
622: if (J) *J = sdm->ops->computejacobian;
623: if (ctx) {
624: if (sdm->jacobianctxcontainer) PetscContainerGetPointer(sdm->jacobianctxcontainer, ctx);
625: else *ctx = NULL;
626: }
627: return 0;
628: }
630: /*@C
631: DMSNESSetPicard - set SNES Picard iteration matrix and RHS evaluation functions.
633: Not Collective
635: Input Parameters:
636: + dm - `DM` to be used with `SNES`
637: . b - RHS evaluation function
638: . J - Picard matrix evaluation function
639: - ctx - context for residual evaluation
641: Level: advanced
643: .seealso: `SNESSetPicard()`, `DMSNESSetFunction()`, `DMSNESSetJacobian()`
644: @*/
645: PetscErrorCode DMSNESSetPicard(DM dm, PetscErrorCode (*b)(SNES, Vec, Vec, void *), PetscErrorCode (*J)(SNES, Vec, Mat, Mat, void *), void *ctx)
646: {
647: DMSNES sdm;
650: DMGetDMSNES(dm, &sdm);
651: if (b) sdm->ops->computepfunction = b;
652: if (J) sdm->ops->computepjacobian = J;
653: if (ctx) sdm->pctx = ctx;
654: return 0;
655: }
657: /*@C
658: DMSNESGetPicard - get `SNES` Picard iteration evaluation functions
660: Not Collective
662: Input Parameter:
663: . dm - `DM` to be used with `SNES`
665: Output Parameters:
666: + b - RHS evaluation function; see `SNESFunction` for details
667: . J - RHS evaluation function; see `SNESJacobianFunction` for detailsa
668: - ctx - context for residual evaluation
670: Level: advanced
672: .seealso: `DMSNESSetContext()`, `SNESSetFunction()`, `DMSNESSetJacobian()`
673: @*/
674: PetscErrorCode DMSNESGetPicard(DM dm, PetscErrorCode (**b)(SNES, Vec, Vec, void *), PetscErrorCode (**J)(SNES, Vec, Mat, Mat, void *), void **ctx)
675: {
676: DMSNES sdm;
679: DMGetDMSNES(dm, &sdm);
680: if (b) *b = sdm->ops->computepfunction;
681: if (J) *J = sdm->ops->computepjacobian;
682: if (ctx) *ctx = sdm->pctx;
683: return 0;
684: }