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: PetscFunctionBegin;
28: PetscCall(PetscFree(sdm->data));
29: PetscFunctionReturn(PETSC_SUCCESS);
30: }
32: static PetscErrorCode DMSNESDuplicate_DMDA(DMSNES oldsdm, DMSNES sdm)
33: {
34: PetscFunctionBegin;
35: PetscCall(PetscNew((DMSNES_DA **)&sdm->data));
36: if (oldsdm->data) PetscCall(PetscMemcpy(sdm->data, oldsdm->data, sizeof(DMSNES_DA)));
37: PetscFunctionReturn(PETSC_SUCCESS);
38: }
40: static PetscErrorCode DMDASNESGetContext(DM dm, DMSNES sdm, DMSNES_DA **dmdasnes)
41: {
42: PetscFunctionBegin;
43: *dmdasnes = NULL;
44: if (!sdm->data) {
45: PetscCall(PetscNew((DMSNES_DA **)&sdm->data));
46: sdm->ops->destroy = DMSNESDestroy_DMDA;
47: sdm->ops->duplicate = DMSNESDuplicate_DMDA;
48: }
49: *dmdasnes = (DMSNES_DA *)sdm->data;
50: PetscFunctionReturn(PETSC_SUCCESS);
51: }
53: static PetscErrorCode SNESComputeFunction_DMDA(SNES snes, Vec X, Vec F, void *ctx)
54: {
55: DM dm;
56: DMSNES_DA *dmdasnes = (DMSNES_DA *)ctx;
57: DMDALocalInfo info;
58: Vec Xloc;
59: void *x, *f;
61: PetscFunctionBegin;
65: PetscCheck(dmdasnes->residuallocal || dmdasnes->residuallocalvec, PetscObjectComm((PetscObject)snes), PETSC_ERR_PLIB, "Corrupt context");
66: PetscCall(SNESGetDM(snes, &dm));
67: PetscCall(DMGetLocalVector(dm, &Xloc));
68: PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc));
69: PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc));
70: PetscCall(DMDAGetLocalInfo(dm, &info));
71: switch (dmdasnes->residuallocalimode) {
72: case INSERT_VALUES: {
73: PetscCall(PetscLogEventBegin(SNES_FunctionEval, snes, X, F, 0));
74: if (dmdasnes->residuallocalvec) PetscCallBack("SNES DMDA local callback function", (*dmdasnes->residuallocalvec)(&info, Xloc, F, dmdasnes->residuallocalctx));
75: else {
76: PetscCall(DMDAVecGetArray(dm, Xloc, &x));
77: PetscCall(DMDAVecGetArray(dm, F, &f));
78: PetscCallBack("SNES DMDA local callback function", (*dmdasnes->residuallocal)(&info, x, f, dmdasnes->residuallocalctx));
79: PetscCall(DMDAVecRestoreArray(dm, Xloc, &x));
80: PetscCall(DMDAVecRestoreArray(dm, F, &f));
81: }
82: PetscCall(PetscLogEventEnd(SNES_FunctionEval, snes, X, F, 0));
83: } break;
84: case ADD_VALUES: {
85: Vec Floc;
86: PetscCall(DMGetLocalVector(dm, &Floc));
87: PetscCall(VecZeroEntries(Floc));
88: PetscCall(PetscLogEventBegin(SNES_FunctionEval, snes, X, F, 0));
89: if (dmdasnes->residuallocalvec) PetscCallBack("SNES DMDA local callback function", (*dmdasnes->residuallocalvec)(&info, Xloc, Floc, dmdasnes->residuallocalctx));
90: else {
91: PetscCall(DMDAVecGetArray(dm, Xloc, &x));
92: PetscCall(DMDAVecGetArray(dm, Floc, &f));
93: PetscCallBack("SNES DMDA local callback function", (*dmdasnes->residuallocal)(&info, x, f, dmdasnes->residuallocalctx));
94: PetscCall(DMDAVecRestoreArray(dm, Xloc, &x));
95: PetscCall(DMDAVecRestoreArray(dm, Floc, &f));
96: }
97: PetscCall(PetscLogEventEnd(SNES_FunctionEval, snes, X, F, 0));
98: PetscCall(VecZeroEntries(F));
99: PetscCall(DMLocalToGlobalBegin(dm, Floc, ADD_VALUES, F));
100: PetscCall(DMLocalToGlobalEnd(dm, Floc, ADD_VALUES, F));
101: PetscCall(DMRestoreLocalVector(dm, &Floc));
102: } break;
103: default:
104: SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_INCOMP, "Cannot use imode=%d", (int)dmdasnes->residuallocalimode);
105: }
106: PetscCall(DMRestoreLocalVector(dm, &Xloc));
107: if (snes->domainerror) PetscCall(VecSetInf(F));
108: PetscFunctionReturn(PETSC_SUCCESS);
109: }
111: static PetscErrorCode SNESComputeObjective_DMDA(SNES snes, Vec X, PetscReal *ob, void *ctx)
112: {
113: DM dm;
114: DMSNES_DA *dmdasnes = (DMSNES_DA *)ctx;
115: DMDALocalInfo info;
116: Vec Xloc;
117: void *x;
119: PetscFunctionBegin;
123: PetscCheck(dmdasnes->objectivelocal || dmdasnes->objectivelocalvec, PetscObjectComm((PetscObject)snes), PETSC_ERR_PLIB, "Corrupt context");
124: PetscCall(SNESGetDM(snes, &dm));
125: PetscCall(DMGetLocalVector(dm, &Xloc));
126: PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc));
127: PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc));
128: PetscCall(DMDAGetLocalInfo(dm, &info));
129: if (dmdasnes->objectivelocalvec) PetscCallBack("SNES DMDA local callback objective", (*dmdasnes->objectivelocalvec)(&info, Xloc, ob, dmdasnes->objectivelocalctx));
130: else {
131: PetscCall(DMDAVecGetArray(dm, Xloc, &x));
132: PetscCallBack("SNES DMDA local callback objective", (*dmdasnes->objectivelocal)(&info, x, ob, dmdasnes->objectivelocalctx));
133: PetscCall(DMDAVecRestoreArray(dm, Xloc, &x));
134: }
135: PetscCall(DMRestoreLocalVector(dm, &Xloc));
136: PetscFunctionReturn(PETSC_SUCCESS);
137: }
139: /* Routine is called by example, hence must be labeled PETSC_EXTERN */
140: PETSC_EXTERN PetscErrorCode SNESComputeJacobian_DMDA(SNES snes, Vec X, Mat A, Mat B, void *ctx)
141: {
142: DM dm;
143: DMSNES_DA *dmdasnes = (DMSNES_DA *)ctx;
144: DMDALocalInfo info;
145: Vec Xloc;
146: void *x;
148: PetscFunctionBegin;
149: PetscCheck(dmdasnes->residuallocal || dmdasnes->residuallocalvec, PetscObjectComm((PetscObject)snes), PETSC_ERR_PLIB, "Corrupt context");
150: PetscCall(SNESGetDM(snes, &dm));
152: if (dmdasnes->jacobianlocal || dmdasnes->jacobianlocalctx) {
153: PetscCall(DMGetLocalVector(dm, &Xloc));
154: PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc));
155: PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc));
156: PetscCall(DMDAGetLocalInfo(dm, &info));
157: if (dmdasnes->jacobianlocalvec) PetscCallBack("SNES DMDA local callback Jacobian", (*dmdasnes->jacobianlocalvec)(&info, Xloc, A, B, dmdasnes->jacobianlocalctx));
158: else {
159: PetscCall(DMDAVecGetArray(dm, Xloc, &x));
160: PetscCallBack("SNES DMDA local callback Jacobian", (*dmdasnes->jacobianlocal)(&info, x, A, B, dmdasnes->jacobianlocalctx));
161: PetscCall(DMDAVecRestoreArray(dm, Xloc, &x));
162: }
163: PetscCall(DMRestoreLocalVector(dm, &Xloc));
164: } else {
165: MatFDColoring fdcoloring;
166: PetscCall(PetscObjectQuery((PetscObject)dm, "DMDASNES_FDCOLORING", (PetscObject *)&fdcoloring));
167: if (!fdcoloring) {
168: ISColoring coloring;
170: PetscCall(DMCreateColoring(dm, dm->coloringtype, &coloring));
171: PetscCall(MatFDColoringCreate(B, coloring, &fdcoloring));
172: switch (dm->coloringtype) {
173: case IS_COLORING_GLOBAL:
174: PetscCall(MatFDColoringSetFunction(fdcoloring, (PetscErrorCode(*)(void))SNESComputeFunction_DMDA, dmdasnes));
175: break;
176: default:
177: SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_SUP, "No support for coloring type '%s'", ISColoringTypes[dm->coloringtype]);
178: }
179: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)fdcoloring, ((PetscObject)dm)->prefix));
180: PetscCall(MatFDColoringSetFromOptions(fdcoloring));
181: PetscCall(MatFDColoringSetUp(B, coloring, fdcoloring));
182: PetscCall(ISColoringDestroy(&coloring));
183: PetscCall(PetscObjectCompose((PetscObject)dm, "DMDASNES_FDCOLORING", (PetscObject)fdcoloring));
184: PetscCall(PetscObjectDereference((PetscObject)fdcoloring));
186: /* The following breaks an ugly reference counting loop that deserves a paragraph. MatFDColoringApply() will call
187: * VecDuplicate() with the state Vec and store inside the MatFDColoring. This Vec will duplicate the Vec, but the
188: * MatFDColoring is composed with the DM. We dereference the DM here so that the reference count will eventually
189: * drop to 0. Note the code in DMDestroy() that exits early for a negative reference count. That code path will be
190: * taken when the PetscObjectList for the Vec inside MatFDColoring is destroyed.
191: */
192: PetscCall(PetscObjectDereference((PetscObject)dm));
193: }
194: PetscCall(MatFDColoringApply(B, fdcoloring, X, snes));
195: }
196: /* This will be redundant if the user called both, but it's too common to forget. */
197: if (A != B) {
198: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
199: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
200: }
201: PetscFunctionReturn(PETSC_SUCCESS);
202: }
204: /*@C
205: DMDASNESSetFunctionLocal - set a local residual evaluation function for use with `DMDA`
207: Logically Collective
209: Input Parameters:
210: + dm - `DM` to associate callback with
211: . imode - `INSERT_VALUES` if local function computes owned part, `ADD_VALUES` if it contributes to ghosted part
212: . func - local residual evaluation
213: - ctx - optional context for local residual evaluation
215: Calling sequence of `func`:
216: $ PetscErrorCode func(DMDALocalInfo *info, void *x, void *f, void *ctx)
217: + info - `DMDALocalInfo` defining the subdomain to evaluate the residual on
218: . x - dimensional pointer to state at which to evaluate residual (e.g. PetscScalar *x or **x or ***x)
219: . f - dimensional pointer to residual, write the residual here (e.g. PetscScalar *f or **f or ***f)
220: - ctx - optional context passed above
222: Level: beginner
224: .seealso: `DMDA`, `DMDASNESSetJacobianLocal()`, `DMSNESSetFunction()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
225: @*/
226: PetscErrorCode DMDASNESSetFunctionLocal(DM dm, InsertMode imode, PetscErrorCode (*func)(DMDALocalInfo *, void *, void *, void *), void *ctx)
227: {
228: DMSNES sdm;
229: DMSNES_DA *dmdasnes;
231: PetscFunctionBegin;
233: PetscCall(DMGetDMSNESWrite(dm, &sdm));
234: PetscCall(DMDASNESGetContext(dm, sdm, &dmdasnes));
236: dmdasnes->residuallocalimode = imode;
237: dmdasnes->residuallocal = func;
238: dmdasnes->residuallocalctx = ctx;
240: PetscCall(DMSNESSetFunction(dm, SNESComputeFunction_DMDA, dmdasnes));
241: if (!sdm->ops->computejacobian) { /* Call us for the Jacobian too, can be overridden by the user. */
242: PetscCall(DMSNESSetJacobian(dm, SNESComputeJacobian_DMDA, dmdasnes));
243: }
244: PetscFunctionReturn(PETSC_SUCCESS);
245: }
247: /*@C
248: DMDASNESSetFunctionLocalVec - set a local residual evaluation function that operates on a local vector for `DMDA`
250: Logically Collective
252: Input Parameters:
253: + dm - `DM` to associate callback with
254: . imode - `INSERT_VALUES` if local function computes owned part, `ADD_VALUES` if it contributes to ghosted part
255: . func - local residual evaluation
256: - ctx - optional context for local residual evaluation
258: Calling sequence of `func`:
259: $ PetscErrorCode func(DMDALocalInfo *info, Vec x, Vec f, void *ctx),
260: + info - `DMDALocalInfo` defining the subdomain to evaluate the residual on
261: . x - state vector at which to evaluate residual
262: . f - residual vector
263: - ctx - optional context passed above
265: Level: beginner
267: .seealso: `DMDA`, `DMDASNESSetFunctionLocal()`, `DMDASNESSetJacobianLocalVec()`, `DMSNESSetFunction()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
268: @*/
269: PetscErrorCode DMDASNESSetFunctionLocalVec(DM dm, InsertMode imode, PetscErrorCode (*func)(DMDALocalInfo *, Vec, Vec, void *), void *ctx)
270: {
271: DMSNES sdm;
272: DMSNES_DA *dmdasnes;
274: PetscFunctionBegin;
276: PetscCall(DMGetDMSNESWrite(dm, &sdm));
277: PetscCall(DMDASNESGetContext(dm, sdm, &dmdasnes));
279: dmdasnes->residuallocalimode = imode;
280: dmdasnes->residuallocalvec = func;
281: dmdasnes->residuallocalctx = ctx;
283: PetscCall(DMSNESSetFunction(dm, SNESComputeFunction_DMDA, dmdasnes));
284: if (!sdm->ops->computejacobian) { /* Call us for the Jacobian too, can be overridden by the user. */
285: PetscCall(DMSNESSetJacobian(dm, SNESComputeJacobian_DMDA, dmdasnes));
286: }
287: PetscFunctionReturn(PETSC_SUCCESS);
288: }
290: /*@C
291: DMDASNESSetJacobianLocal - set a local Jacobian evaluation function for use with `DMDA`
293: Logically Collective
295: Input Parameters:
296: + dm - `DM` to associate callback with
297: . func - local Jacobian evaluation
298: - ctx - optional context for local Jacobian evaluation
300: Calling sequence of `func`:
301: $ PetscErrorCode func(DMDALocalInfo *info, void *x, Mat J, Mat M, void *ctx),
302: + info - `DMDALocalInfo` defining the subdomain to evaluate the Jacobian at
303: . x - dimensional pointer to state at which to evaluate Jacobian (e.g. PetscScalar *x or **x or ***x)
304: . J - Mat object for the Jacobian
305: . M - Mat object for the Jacobian preconditioner matrix, often `J`
306: - ctx - optional context passed above
308: Level: beginner
310: .seealso: `DMDA`, `DMDASNESSetFunctionLocal()`, `DMSNESSetJacobian()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
311: @*/
312: PetscErrorCode DMDASNESSetJacobianLocal(DM dm, PetscErrorCode (*func)(DMDALocalInfo *, void *, Mat, Mat, void *), void *ctx)
313: {
314: DMSNES sdm;
315: DMSNES_DA *dmdasnes;
317: PetscFunctionBegin;
319: PetscCall(DMGetDMSNESWrite(dm, &sdm));
320: PetscCall(DMDASNESGetContext(dm, sdm, &dmdasnes));
322: dmdasnes->jacobianlocal = func;
323: dmdasnes->jacobianlocalctx = ctx;
325: PetscCall(DMSNESSetJacobian(dm, SNESComputeJacobian_DMDA, dmdasnes));
326: PetscFunctionReturn(PETSC_SUCCESS);
327: }
329: /*@C
330: DMDASNESSetJacobianLocalVec - set a local Jacobian evaluation function that operates on a local vector with `DMDA`
332: Logically Collective
334: Input Parameters:
335: + dm - `DM` to associate callback with
336: . func - local Jacobian evaluation
337: - ctx - optional context for local Jacobian evaluation
339: Calling sequence of `func`:
340: $ PetscErrorCode (*func)(DMDALocalInfo *info, Vec x, Mat J, Mat M, void *ctx),
341: + info - `DMDALocalInfo` defining the subdomain to evaluate the Jacobian at
342: . x - state vector at which to evaluate Jacobian
343: . J - Mat object for the Jacobian
344: . M - Mat object for the Jacobian preconditioner matrix, often `J`
345: - ctx - optional context passed above
347: Level: beginner
349: .seealso: `DMDA`, `DMDASNESSetJacobianLocal()`, `DMDASNESSetFunctionLocalVec()`, `DMSNESSetJacobian()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
350: @*/
351: PetscErrorCode DMDASNESSetJacobianLocalVec(DM dm, PetscErrorCode (*func)(DMDALocalInfo *, Vec, Mat, Mat, void *), void *ctx)
352: {
353: DMSNES sdm;
354: DMSNES_DA *dmdasnes;
356: PetscFunctionBegin;
358: PetscCall(DMGetDMSNESWrite(dm, &sdm));
359: PetscCall(DMDASNESGetContext(dm, sdm, &dmdasnes));
361: dmdasnes->jacobianlocalvec = func;
362: dmdasnes->jacobianlocalctx = ctx;
364: PetscCall(DMSNESSetJacobian(dm, SNESComputeJacobian_DMDA, dmdasnes));
365: PetscFunctionReturn(PETSC_SUCCESS);
366: }
368: /*@C
369: DMDASNESSetObjectiveLocal - set a local residual evaluation function to used with a `DMDA`
371: Logically Collective
373: Input Parameters:
374: + dm - `DM` to associate callback with
375: . func - local objective evaluation
376: - ctx - optional context for local residual evaluation
378: Calling sequence of `func`:
379: $ PetscErrorCode func(DMDALocalInfo *info, void *x, PetscReal obj, void *ctx);
380: + info - defines the subdomain to evaluate the residual on
381: . x - dimensional pointer to state at which to evaluate residual
382: . ob - eventual objective value
383: - ctx - optional context passed above
385: Level: beginner
387: .seealso: `DMDA`, `DMSNESSetFunction()`, `DMDASNESSetJacobianLocal()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
388: @*/
389: PetscErrorCode DMDASNESSetObjectiveLocal(DM dm, DMDASNESObjective func, void *ctx)
390: {
391: DMSNES sdm;
392: DMSNES_DA *dmdasnes;
394: PetscFunctionBegin;
396: PetscCall(DMGetDMSNESWrite(dm, &sdm));
397: PetscCall(DMDASNESGetContext(dm, sdm, &dmdasnes));
399: dmdasnes->objectivelocal = func;
400: dmdasnes->objectivelocalctx = ctx;
402: PetscCall(DMSNESSetObjective(dm, SNESComputeObjective_DMDA, dmdasnes));
403: PetscFunctionReturn(PETSC_SUCCESS);
404: }
406: /*@C
407: DMDASNESSetObjectiveLocal - set a local residual evaluation function that operates on a local vector with `DMDA`
409: Logically Collective
411: Input Parameters:
412: + dm - `DM` to associate callback with
413: . func - local objective evaluation
414: - ctx - optional context for local residual evaluation
416: Calling sequence of `func`:
417: $ PetscErrorCode func(DMDALocalInfo *info, Vec x, PetscReal *ob, void *ctx);
418: + info - defines the subdomain to evaluate the residual on
419: . x - state vector at which to evaluate residual
420: . ob - eventual objective value
421: - ctx - optional context passed above
423: Level: beginner
425: .seealso: `DMDA`, `DMDASNESSetObjectiveLocal()`, `DMSNESSetFunction()`, `DMDASNESSetJacobianLocalVec()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
426: @*/
427: PetscErrorCode DMDASNESSetObjectiveLocalVec(DM dm, DMDASNESObjectiveVec func, void *ctx)
428: {
429: DMSNES sdm;
430: DMSNES_DA *dmdasnes;
432: PetscFunctionBegin;
434: PetscCall(DMGetDMSNESWrite(dm, &sdm));
435: PetscCall(DMDASNESGetContext(dm, sdm, &dmdasnes));
437: dmdasnes->objectivelocalvec = func;
438: dmdasnes->objectivelocalctx = ctx;
440: PetscCall(DMSNESSetObjective(dm, SNESComputeObjective_DMDA, dmdasnes));
441: PetscFunctionReturn(PETSC_SUCCESS);
442: }
444: static PetscErrorCode SNESComputePicard_DMDA(SNES snes, Vec X, Vec F, void *ctx)
445: {
446: DM dm;
447: DMSNES_DA *dmdasnes = (DMSNES_DA *)ctx;
448: DMDALocalInfo info;
449: Vec Xloc;
450: void *x, *f;
452: PetscFunctionBegin;
456: PetscCheck(dmdasnes->rhsplocal, PetscObjectComm((PetscObject)snes), PETSC_ERR_PLIB, "Corrupt context");
457: PetscCall(SNESGetDM(snes, &dm));
458: PetscCall(DMGetLocalVector(dm, &Xloc));
459: PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc));
460: PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc));
461: PetscCall(DMDAGetLocalInfo(dm, &info));
462: PetscCall(DMDAVecGetArray(dm, Xloc, &x));
463: switch (dmdasnes->residuallocalimode) {
464: case INSERT_VALUES: {
465: PetscCall(DMDAVecGetArray(dm, F, &f));
466: PetscCallBack("SNES Picard DMDA local callback function", (*dmdasnes->rhsplocal)(&info, x, f, dmdasnes->picardlocalctx));
467: PetscCall(DMDAVecRestoreArray(dm, F, &f));
468: } break;
469: case ADD_VALUES: {
470: Vec Floc;
471: PetscCall(DMGetLocalVector(dm, &Floc));
472: PetscCall(VecZeroEntries(Floc));
473: PetscCall(DMDAVecGetArray(dm, Floc, &f));
474: PetscCallBack("SNES Picard DMDA local callback function", (*dmdasnes->rhsplocal)(&info, x, f, dmdasnes->picardlocalctx));
475: PetscCall(DMDAVecRestoreArray(dm, Floc, &f));
476: PetscCall(VecZeroEntries(F));
477: PetscCall(DMLocalToGlobalBegin(dm, Floc, ADD_VALUES, F));
478: PetscCall(DMLocalToGlobalEnd(dm, Floc, ADD_VALUES, F));
479: PetscCall(DMRestoreLocalVector(dm, &Floc));
480: } break;
481: default:
482: SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_ARG_INCOMP, "Cannot use imode=%d", (int)dmdasnes->residuallocalimode);
483: }
484: PetscCall(DMDAVecRestoreArray(dm, Xloc, &x));
485: PetscCall(DMRestoreLocalVector(dm, &Xloc));
486: PetscFunctionReturn(PETSC_SUCCESS);
487: }
489: static PetscErrorCode SNESComputePicardJacobian_DMDA(SNES snes, Vec X, Mat A, Mat B, void *ctx)
490: {
491: DM dm;
492: DMSNES_DA *dmdasnes = (DMSNES_DA *)ctx;
493: DMDALocalInfo info;
494: Vec Xloc;
495: void *x;
497: PetscFunctionBegin;
498: PetscCheck(dmdasnes->jacobianplocal, PetscObjectComm((PetscObject)snes), PETSC_ERR_PLIB, "Corrupt context");
499: PetscCall(SNESGetDM(snes, &dm));
501: PetscCall(DMGetLocalVector(dm, &Xloc));
502: PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc));
503: PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc));
504: PetscCall(DMDAGetLocalInfo(dm, &info));
505: PetscCall(DMDAVecGetArray(dm, Xloc, &x));
506: PetscCallBack("SNES Picard DMDA local callback Jacobian", (*dmdasnes->jacobianplocal)(&info, x, A, B, dmdasnes->picardlocalctx));
507: PetscCall(DMDAVecRestoreArray(dm, Xloc, &x));
508: PetscCall(DMRestoreLocalVector(dm, &Xloc));
509: if (A != B) {
510: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
511: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
512: }
513: PetscFunctionReturn(PETSC_SUCCESS);
514: }
516: /*@C
517: DMDASNESSetPicardLocal - set a local right hand side and matrix evaluation function for Picard iteration with `DMDA`
519: Logically Collective
521: Input Parameters:
522: + dm - `DM` to associate callback with
523: . imode - `INSERT_VALUES` if local function computes owned part, `ADD_VALUES` if it contributes to ghosted part
524: . func - local residual evaluation
525: - ctx - optional context for local residual evaluation
527: Calling sequence of `func`:
528: $ PetscErrorCode func(DMDALocalInfo *info, void *x, void *f, void *ctx);
529: + info - defines the subdomain to evaluate the residual on
530: . x - dimensional pointer to state at which to evaluate residual
531: . f - dimensional pointer to residual, write the residual here
532: - ctx - optional context passed above
534: Level: beginner
536: Note:
537: The user must use `SNESSetFunction`(snes,`NULL`,`SNESPicardComputeFunction`,&user));
538: in their code before calling this routine.
540: .seealso: `DMDA`, `DMSNESSetFunction()`, `DMDASNESSetJacobian()`, `DMDACreate1d()`, `DMDACreate2d()`, `DMDACreate3d()`
541: @*/
542: PetscErrorCode DMDASNESSetPicardLocal(DM dm, InsertMode imode, PetscErrorCode (*func)(DMDALocalInfo *, void *, void *, void *), PetscErrorCode (*jac)(DMDALocalInfo *, void *, Mat, Mat, void *), void *ctx)
543: {
544: DMSNES sdm;
545: DMSNES_DA *dmdasnes;
547: PetscFunctionBegin;
549: PetscCall(DMGetDMSNESWrite(dm, &sdm));
550: PetscCall(DMDASNESGetContext(dm, sdm, &dmdasnes));
552: dmdasnes->residuallocalimode = imode;
553: dmdasnes->rhsplocal = func;
554: dmdasnes->jacobianplocal = jac;
555: dmdasnes->picardlocalctx = ctx;
557: PetscCall(DMSNESSetPicard(dm, SNESComputePicard_DMDA, SNESComputePicardJacobian_DMDA, dmdasnes));
558: PetscCall(DMSNESSetMFFunction(dm, SNESComputeFunction_DMDA, dmdasnes));
559: PetscFunctionReturn(PETSC_SUCCESS);
560: }