Actual source code: dmdasnes.c
1: #include <petscdmda.h>
2: #include <petsc/private/dmimpl.h>
3: #include <petsc/private/snesimpl.h>
5: /* This structure holds the user-provided DMDA callbacks */
6: typedef struct {
7: PetscErrorCode (*residuallocal)(DMDALocalInfo *, void *, void *, void *);
8: PetscErrorCode (*jacobianlocal)(DMDALocalInfo *, void *, Mat, Mat, void *);
9: PetscErrorCode (*objectivelocal)(DMDALocalInfo *, void *, PetscReal *, void *);
11: PetscErrorCode (*residuallocalvec)(DMDALocalInfo *, Vec, Vec, void *);
12: PetscErrorCode (*jacobianlocalvec)(DMDALocalInfo *, Vec, Mat, Mat, void *);
13: PetscErrorCode (*objectivelocalvec)(DMDALocalInfo *, Vec, PetscReal *, void *);
14: void *residuallocalctx;
15: void *jacobianlocalctx;
16: void *objectivelocalctx;
17: InsertMode residuallocalimode;
19: /* For Picard iteration defined locally */
20: PetscErrorCode (*rhsplocal)(DMDALocalInfo *, void *, void *, void *);
21: PetscErrorCode (*jacobianplocal)(DMDALocalInfo *, void *, Mat, Mat, void *);
22: void *picardlocalctx;
23: } DMSNES_DA;
25: static PetscErrorCode DMSNESDestroy_DMDA(DMSNES sdm)
26: {
27: PetscFree(sdm->data);
28: return 0;
29: }
31: static PetscErrorCode DMSNESDuplicate_DMDA(DMSNES oldsdm, DMSNES sdm)
32: {
33: PetscNew((DMSNES_DA **)&sdm->data);
34: if (oldsdm->data) PetscMemcpy(sdm->data, oldsdm->data, sizeof(DMSNES_DA));
35: return 0;
36: }
38: static PetscErrorCode DMDASNESGetContext(DM dm, DMSNES sdm, DMSNES_DA **dmdasnes)
39: {
40: *dmdasnes = NULL;
41: if (!sdm->data) {
42: PetscNew((DMSNES_DA **)&sdm->data);
43: sdm->ops->destroy = DMSNESDestroy_DMDA;
44: sdm->ops->duplicate = DMSNESDuplicate_DMDA;
45: }
46: *dmdasnes = (DMSNES_DA *)sdm->data;
47: return 0;
48: }
50: static PetscErrorCode SNESComputeFunction_DMDA(SNES snes, Vec X, Vec F, void *ctx)
51: {
52: DM dm;
53: DMSNES_DA *dmdasnes = (DMSNES_DA *)ctx;
54: DMDALocalInfo info;
55: Vec Xloc;
56: void *x, *f;
62: SNESGetDM(snes, &dm);
63: DMGetLocalVector(dm, &Xloc);
64: DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc);
65: DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc);
66: DMDAGetLocalInfo(dm, &info);
67: switch (dmdasnes->residuallocalimode) {
68: case INSERT_VALUES: {
69: PetscLogEventBegin(SNES_FunctionEval, snes, X, F, 0);
70: if (dmdasnes->residuallocalvec) PetscCallBack("SNES DMDA local callback function", (*dmdasnes->residuallocalvec)(&info, Xloc, F, dmdasnes->residuallocalctx));
71: else {
72: DMDAVecGetArray(dm, Xloc, &x);
73: DMDAVecGetArray(dm, F, &f);
74: PetscCallBack("SNES DMDA local callback function", (*dmdasnes->residuallocal)(&info, x, f, dmdasnes->residuallocalctx));
75: DMDAVecRestoreArray(dm, Xloc, &x);
76: DMDAVecRestoreArray(dm, F, &f);
77: }
78: PetscLogEventEnd(SNES_FunctionEval, snes, X, F, 0);
79: } break;
80: case ADD_VALUES: {
81: Vec Floc;
82: DMGetLocalVector(dm, &Floc);
83: VecZeroEntries(Floc);
84: PetscLogEventBegin(SNES_FunctionEval, snes, X, F, 0);
85: if (dmdasnes->residuallocalvec) PetscCallBack("SNES DMDA local callback function", (*dmdasnes->residuallocalvec)(&info, Xloc, Floc, dmdasnes->residuallocalctx));
86: else {
87: DMDAVecGetArray(dm, Xloc, &x);
88: DMDAVecGetArray(dm, Floc, &f);
89: PetscCallBack("SNES DMDA local callback function", (*dmdasnes->residuallocal)(&info, x, f, dmdasnes->residuallocalctx));
90: DMDAVecRestoreArray(dm, Xloc, &x);
91: DMDAVecRestoreArray(dm, Floc, &f);
92: }
93: PetscLogEventEnd(SNES_FunctionEval, snes, X, F, 0);
94: VecZeroEntries(F);
95: DMLocalToGlobalBegin(dm, Floc, ADD_VALUES, F);
96: DMLocalToGlobalEnd(dm, Floc, ADD_VALUES, F);
97: DMRestoreLocalVector(dm, &Floc);
98: } break;
99: default:
100: SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_INCOMP, "Cannot use imode=%d", (int)dmdasnes->residuallocalimode);
101: }
102: DMRestoreLocalVector(dm, &Xloc);
103: if (snes->domainerror) VecSetInf(F);
104: return 0;
105: }
107: static PetscErrorCode SNESComputeObjective_DMDA(SNES snes, Vec X, PetscReal *ob, void *ctx)
108: {
109: DM dm;
110: DMSNES_DA *dmdasnes = (DMSNES_DA *)ctx;
111: DMDALocalInfo info;
112: Vec Xloc;
113: void *x;
119: SNESGetDM(snes, &dm);
120: DMGetLocalVector(dm, &Xloc);
121: DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc);
122: DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc);
123: DMDAGetLocalInfo(dm, &info);
124: if (dmdasnes->objectivelocalvec) PetscCallBack("SNES DMDA local callback objective", (*dmdasnes->objectivelocalvec)(&info, Xloc, ob, dmdasnes->objectivelocalctx));
125: else {
126: DMDAVecGetArray(dm, Xloc, &x);
127: PetscCallBack("SNES DMDA local callback objective", (*dmdasnes->objectivelocal)(&info, x, ob, dmdasnes->objectivelocalctx));
128: DMDAVecRestoreArray(dm, Xloc, &x);
129: }
130: DMRestoreLocalVector(dm, &Xloc);
131: return 0;
132: }
134: /* Routine is called by example, hence must be labeled PETSC_EXTERN */
135: PETSC_EXTERN PetscErrorCode SNESComputeJacobian_DMDA(SNES snes, Vec X, Mat A, Mat B, void *ctx)
136: {
137: DM dm;
138: DMSNES_DA *dmdasnes = (DMSNES_DA *)ctx;
139: DMDALocalInfo info;
140: Vec Xloc;
141: void *x;
144: SNESGetDM(snes, &dm);
146: if (dmdasnes->jacobianlocal || dmdasnes->jacobianlocalctx) {
147: DMGetLocalVector(dm, &Xloc);
148: DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc);
149: DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc);
150: DMDAGetLocalInfo(dm, &info);
151: if (dmdasnes->jacobianlocalvec) PetscCallBack("SNES DMDA local callback Jacobian", (*dmdasnes->jacobianlocalvec)(&info, Xloc, A, B, dmdasnes->jacobianlocalctx));
152: else {
153: DMDAVecGetArray(dm, Xloc, &x);
154: PetscCallBack("SNES DMDA local callback Jacobian", (*dmdasnes->jacobianlocal)(&info, x, A, B, dmdasnes->jacobianlocalctx));
155: DMDAVecRestoreArray(dm, Xloc, &x);
156: }
157: DMRestoreLocalVector(dm, &Xloc);
158: } else {
159: MatFDColoring fdcoloring;
160: PetscObjectQuery((PetscObject)dm, "DMDASNES_FDCOLORING", (PetscObject *)&fdcoloring);
161: if (!fdcoloring) {
162: ISColoring coloring;
164: DMCreateColoring(dm, dm->coloringtype, &coloring);
165: MatFDColoringCreate(B, coloring, &fdcoloring);
166: switch (dm->coloringtype) {
167: case IS_COLORING_GLOBAL:
168: MatFDColoringSetFunction(fdcoloring, (PetscErrorCode(*)(void))SNESComputeFunction_DMDA, dmdasnes);
169: break;
170: default:
171: SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_SUP, "No support for coloring type '%s'", ISColoringTypes[dm->coloringtype]);
172: }
173: PetscObjectSetOptionsPrefix((PetscObject)fdcoloring, ((PetscObject)dm)->prefix);
174: MatFDColoringSetFromOptions(fdcoloring);
175: MatFDColoringSetUp(B, coloring, fdcoloring);
176: ISColoringDestroy(&coloring);
177: PetscObjectCompose((PetscObject)dm, "DMDASNES_FDCOLORING", (PetscObject)fdcoloring);
178: PetscObjectDereference((PetscObject)fdcoloring);
180: /* The following breaks an ugly reference counting loop that deserves a paragraph. MatFDColoringApply() will call
181: * VecDuplicate() with the state Vec and store inside the MatFDColoring. This Vec will duplicate the Vec, but the
182: * MatFDColoring is composed with the DM. We dereference the DM here so that the reference count will eventually
183: * drop to 0. Note the code in DMDestroy() that exits early for a negative reference count. That code path will be
184: * taken when the PetscObjectList for the Vec inside MatFDColoring is destroyed.
185: */
186: PetscObjectDereference((PetscObject)dm);
187: }
188: MatFDColoringApply(B, fdcoloring, X, snes);
189: }
190: /* This will be redundant if the user called both, but it's too common to forget. */
191: if (A != B) {
192: MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
193: MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
194: }
195: return 0;
196: }
198: /*@C
199: DMDASNESSetFunctionLocal - set a local residual evaluation function for use with `DMDA`
201: Logically Collective
203: Input Parameters:
204: + dm - `DM` to associate callback with
205: . imode - `INSERT_VALUES` if local function computes owned part, `ADD_VALUES` if it contributes to ghosted part
206: . func - local residual evaluation
207: - ctx - optional context for local residual evaluation
209: Calling sequence:
210: For PetscErrorCode (*func)(DMDALocalInfo *info,void *x, void *f, void *ctx),
211: + info - `DMDALocalInfo` defining the subdomain to evaluate the residual on
212: . x - dimensional pointer to state at which to evaluate residual (e.g. PetscScalar *x or **x or ***x)
213: . f - dimensional pointer to residual, write the residual here (e.g. PetscScalar *f or **f or ***f)
214: - ctx - optional context passed above
216: Level: beginner
218: .seealso: `DMDA`, `DMDASNESSetJacobianLocal()`, `DMSNESSetFunction()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
219: @*/
220: PetscErrorCode DMDASNESSetFunctionLocal(DM dm, InsertMode imode, PetscErrorCode (*func)(DMDALocalInfo *, void *, void *, void *), void *ctx)
221: {
222: DMSNES sdm;
223: DMSNES_DA *dmdasnes;
226: DMGetDMSNESWrite(dm, &sdm);
227: DMDASNESGetContext(dm, sdm, &dmdasnes);
229: dmdasnes->residuallocalimode = imode;
230: dmdasnes->residuallocal = func;
231: dmdasnes->residuallocalctx = ctx;
233: DMSNESSetFunction(dm, SNESComputeFunction_DMDA, dmdasnes);
234: if (!sdm->ops->computejacobian) { /* Call us for the Jacobian too, can be overridden by the user. */
235: DMSNESSetJacobian(dm, SNESComputeJacobian_DMDA, dmdasnes);
236: }
237: return 0;
238: }
240: /*@C
241: DMDASNESSetFunctionLocalVec - set a local residual evaluation function that operates on a local vector for `DMDA`
243: Logically Collective
245: Input Parameters:
246: + dm - `DM` to associate callback with
247: . imode - `INSERT_VALUES` if local function computes owned part, `ADD_VALUES` if it contributes to ghosted part
248: . func - local residual evaluation
249: - ctx - optional context for local residual evaluation
251: Calling sequence:
252: For PetscErrorCode (*func)(DMDALocalInfo *info,Vec x, Vec f, void *ctx),
253: + info - `DMDALocalInfo` defining the subdomain to evaluate the residual on
254: . x - state vector at which to evaluate residual
255: . f - residual vector
256: - ctx - optional context passed above
258: Level: beginner
260: .seealso: `DMDA`, `DMDASNESSetFunctionLocal()`, `DMDASNESSetJacobianLocalVec()`, `DMSNESSetFunction()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
261: @*/
262: PetscErrorCode DMDASNESSetFunctionLocalVec(DM dm, InsertMode imode, PetscErrorCode (*func)(DMDALocalInfo *, Vec, Vec, void *), void *ctx)
263: {
264: DMSNES sdm;
265: DMSNES_DA *dmdasnes;
268: DMGetDMSNESWrite(dm, &sdm);
269: DMDASNESGetContext(dm, sdm, &dmdasnes);
271: dmdasnes->residuallocalimode = imode;
272: dmdasnes->residuallocalvec = func;
273: dmdasnes->residuallocalctx = ctx;
275: DMSNESSetFunction(dm, SNESComputeFunction_DMDA, dmdasnes);
276: if (!sdm->ops->computejacobian) { /* Call us for the Jacobian too, can be overridden by the user. */
277: DMSNESSetJacobian(dm, SNESComputeJacobian_DMDA, dmdasnes);
278: }
279: return 0;
280: }
282: /*@C
283: DMDASNESSetJacobianLocal - set a local Jacobian evaluation function for use with `DMDA`
285: Logically Collective
287: Input Parameters:
288: + dm - `DM` to associate callback with
289: . func - local Jacobian evaluation
290: - ctx - optional context for local Jacobian evaluation
292: Calling sequence:
293: For PetscErrorCode (*func)(DMDALocalInfo *info,void *x,Mat J,Mat M,void *ctx),
294: + info - `DMDALocalInfo` defining the subdomain to evaluate the Jacobian at
295: . x - dimensional pointer to state at which to evaluate Jacobian (e.g. PetscScalar *x or **x or ***x)
296: . J - Mat object for the Jacobian
297: . M - Mat object for the Jacobian preconditioner matrix
298: - ctx - optional context passed above
300: Level: beginner
302: .seealso: `DMDA`, `DMDASNESSetFunctionLocal()`, `DMSNESSetJacobian()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
303: @*/
304: PetscErrorCode DMDASNESSetJacobianLocal(DM dm, PetscErrorCode (*func)(DMDALocalInfo *, void *, Mat, Mat, void *), void *ctx)
305: {
306: DMSNES sdm;
307: DMSNES_DA *dmdasnes;
310: DMGetDMSNESWrite(dm, &sdm);
311: DMDASNESGetContext(dm, sdm, &dmdasnes);
313: dmdasnes->jacobianlocal = func;
314: dmdasnes->jacobianlocalctx = ctx;
316: DMSNESSetJacobian(dm, SNESComputeJacobian_DMDA, dmdasnes);
317: return 0;
318: }
320: /*@C
321: DMDASNESSetJacobianLocalVec - set a local Jacobian evaluation function that operates on a local vector with `DMDA`
323: Logically Collective
325: Input Parameters:
326: + dm - `DM` to associate callback with
327: . func - local Jacobian evaluation
328: - ctx - optional context for local Jacobian evaluation
330: Calling sequence:
331: For PetscErrorCode (*func)(DMDALocalInfo *info,Vec x,Mat J,Mat M,void *ctx),
332: + info - `DMDALocalInfo` defining the subdomain to evaluate the Jacobian at
333: . x - state vector at which to evaluate Jacobian
334: . J - Mat object for the Jacobian
335: . M - Mat object for the Jacobian preconditioner matrix
336: - ctx - optional context passed above
338: Level: beginner
340: .seealso: `DMDA`, `DMDASNESSetJacobianLocal()`, `DMDASNESSetFunctionLocalVec()`, `DMSNESSetJacobian()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
341: @*/
342: PetscErrorCode DMDASNESSetJacobianLocalVec(DM dm, PetscErrorCode (*func)(DMDALocalInfo *, Vec, Mat, Mat, void *), void *ctx)
343: {
344: DMSNES sdm;
345: DMSNES_DA *dmdasnes;
348: DMGetDMSNESWrite(dm, &sdm);
349: DMDASNESGetContext(dm, sdm, &dmdasnes);
351: dmdasnes->jacobianlocalvec = func;
352: dmdasnes->jacobianlocalctx = ctx;
354: DMSNESSetJacobian(dm, SNESComputeJacobian_DMDA, dmdasnes);
355: return 0;
356: }
358: /*@C
359: DMDASNESSetObjectiveLocal - set a local residual evaluation function to used with a `DMDA`
361: Logically Collective
363: Input Parameters:
364: + dm - `DM` to associate callback with
365: . func - local objective evaluation
366: - ctx - optional context for local residual evaluation
368: Calling sequence for func:
369: + info - `DMDALocalInfo` defining the subdomain to evaluate the residual on
370: . x - dimensional pointer to state at which to evaluate residual
371: . ob - eventual objective value
372: - ctx - optional context passed above
374: Level: beginner
376: .seealso: `DMDA`, `DMSNESSetFunction()`, `DMDASNESSetJacobianLocal()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
377: @*/
378: PetscErrorCode DMDASNESSetObjectiveLocal(DM dm, DMDASNESObjective func, void *ctx)
379: {
380: DMSNES sdm;
381: DMSNES_DA *dmdasnes;
384: DMGetDMSNESWrite(dm, &sdm);
385: DMDASNESGetContext(dm, sdm, &dmdasnes);
387: dmdasnes->objectivelocal = func;
388: dmdasnes->objectivelocalctx = ctx;
390: DMSNESSetObjective(dm, SNESComputeObjective_DMDA, dmdasnes);
391: return 0;
392: }
394: /*@C
395: DMDASNESSetObjectiveLocal - set a local residual evaluation function that operates on a local vector with `DMDA`
397: Logically Collective
399: Input Parameters:
400: + dm - `DM` to associate callback with
401: . func - local objective evaluation
402: - ctx - optional context for local residual evaluation
404: Calling sequence
405: For PetscErrorCode (*func)(DMDALocalInfo *info,Vec x,PetscReal *ob,void *ctx);
406: + info - `DMDALocalInfo` defining the subdomain to evaluate the residual on
407: . x - state vector at which to evaluate residual
408: . ob - eventual objective value
409: - ctx - optional context passed above
411: Level: beginner
413: .seealso: `DMDA`, `DMDASNESSetObjectiveLocal()`, `DMSNESSetFunction()`, `DMDASNESSetJacobianLocalVec()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
414: @*/
415: PetscErrorCode DMDASNESSetObjectiveLocalVec(DM dm, DMDASNESObjectiveVec func, void *ctx)
416: {
417: DMSNES sdm;
418: DMSNES_DA *dmdasnes;
421: DMGetDMSNESWrite(dm, &sdm);
422: DMDASNESGetContext(dm, sdm, &dmdasnes);
424: dmdasnes->objectivelocalvec = func;
425: dmdasnes->objectivelocalctx = ctx;
427: DMSNESSetObjective(dm, SNESComputeObjective_DMDA, dmdasnes);
428: return 0;
429: }
431: static PetscErrorCode SNESComputePicard_DMDA(SNES snes, Vec X, Vec F, void *ctx)
432: {
433: DM dm;
434: DMSNES_DA *dmdasnes = (DMSNES_DA *)ctx;
435: DMDALocalInfo info;
436: Vec Xloc;
437: void *x, *f;
443: SNESGetDM(snes, &dm);
444: DMGetLocalVector(dm, &Xloc);
445: DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc);
446: DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc);
447: DMDAGetLocalInfo(dm, &info);
448: DMDAVecGetArray(dm, Xloc, &x);
449: switch (dmdasnes->residuallocalimode) {
450: case INSERT_VALUES: {
451: DMDAVecGetArray(dm, F, &f);
452: PetscCallBack("SNES Picard DMDA local callback function", (*dmdasnes->rhsplocal)(&info, x, f, dmdasnes->picardlocalctx));
453: DMDAVecRestoreArray(dm, F, &f);
454: } break;
455: case ADD_VALUES: {
456: Vec Floc;
457: DMGetLocalVector(dm, &Floc);
458: VecZeroEntries(Floc);
459: DMDAVecGetArray(dm, Floc, &f);
460: PetscCallBack("SNES Picard DMDA local callback function", (*dmdasnes->rhsplocal)(&info, x, f, dmdasnes->picardlocalctx));
461: DMDAVecRestoreArray(dm, Floc, &f);
462: VecZeroEntries(F);
463: DMLocalToGlobalBegin(dm, Floc, ADD_VALUES, F);
464: DMLocalToGlobalEnd(dm, Floc, ADD_VALUES, F);
465: DMRestoreLocalVector(dm, &Floc);
466: } break;
467: default:
468: SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_INCOMP, "Cannot use imode=%d", (int)dmdasnes->residuallocalimode);
469: }
470: DMDAVecRestoreArray(dm, Xloc, &x);
471: DMRestoreLocalVector(dm, &Xloc);
472: return 0;
473: }
475: static PetscErrorCode SNESComputePicardJacobian_DMDA(SNES snes, Vec X, Mat A, Mat B, void *ctx)
476: {
477: DM dm;
478: DMSNES_DA *dmdasnes = (DMSNES_DA *)ctx;
479: DMDALocalInfo info;
480: Vec Xloc;
481: void *x;
484: SNESGetDM(snes, &dm);
486: DMGetLocalVector(dm, &Xloc);
487: DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc);
488: DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc);
489: DMDAGetLocalInfo(dm, &info);
490: DMDAVecGetArray(dm, Xloc, &x);
491: PetscCallBack("SNES Picard DMDA local callback Jacobian", (*dmdasnes->jacobianplocal)(&info, x, A, B, dmdasnes->picardlocalctx));
492: DMDAVecRestoreArray(dm, Xloc, &x);
493: DMRestoreLocalVector(dm, &Xloc);
494: if (A != B) {
495: MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
496: MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
497: }
498: return 0;
499: }
501: /*@C
502: DMDASNESSetPicardLocal - set a local right hand side and matrix evaluation function for Picard iteration with `DMDA`
504: Logically Collective
506: Input Parameters:
507: + dm - `DM` to associate callback with
508: . imode - `INSERT_VALUES` if local function computes owned part, `ADD_VALUES` if it contributes to ghosted part
509: . func - local residual evaluation
510: - ctx - optional context for local residual evaluation
512: Calling sequence for func:
513: + info - `DMDALocalInfo` defining the subdomain to evaluate the residual on
514: . x - dimensional pointer to state at which to evaluate residual
515: . f - dimensional pointer to residual, write the residual here
516: - ctx - optional context passed above
518: Note:
519: The user must use `SNESSetFunction`(snes,NULL,`SNESPicardComputeFunction`,&user));
520: in their code before calling this routine.
522: Level: beginner
524: .seealso: `DMDA`, `DMSNESSetFunction()`, `DMDASNESSetJacobian()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
525: @*/
526: PetscErrorCode DMDASNESSetPicardLocal(DM dm, InsertMode imode, PetscErrorCode (*func)(DMDALocalInfo *, void *, void *, void *), PetscErrorCode (*jac)(DMDALocalInfo *, void *, Mat, Mat, void *), void *ctx)
527: {
528: DMSNES sdm;
529: DMSNES_DA *dmdasnes;
532: DMGetDMSNESWrite(dm, &sdm);
533: DMDASNESGetContext(dm, sdm, &dmdasnes);
535: dmdasnes->residuallocalimode = imode;
536: dmdasnes->rhsplocal = func;
537: dmdasnes->jacobianplocal = jac;
538: dmdasnes->picardlocalctx = ctx;
540: DMSNESSetPicard(dm, SNESComputePicard_DMDA, SNESComputePicardJacobian_DMDA, dmdasnes);
541: DMSNESSetMFFunction(dm, SNESComputeFunction_DMDA, dmdasnes);
542: return 0;
543: }