Actual source code: dmlocalsnes.c
1: #include <petsc/private/dmimpl.h>
2: #include <petsc/private/snesimpl.h>
4: typedef struct {
5: PetscErrorCode (*objectivelocal)(DM, Vec, PetscReal *, void *);
6: PetscErrorCode (*residuallocal)(DM, Vec, Vec, void *);
7: PetscErrorCode (*jacobianlocal)(DM, Vec, Mat, Mat, void *);
8: PetscErrorCode (*boundarylocal)(DM, Vec, void *);
9: void *objectivelocalctx;
10: void *residuallocalctx;
11: void *jacobianlocalctx;
12: void *boundarylocalctx;
13: } DMSNES_Local;
15: static PetscErrorCode DMSNESDestroy_DMLocal(DMSNES sdm)
16: {
17: PetscFunctionBegin;
18: PetscCall(PetscFree(sdm->data));
19: sdm->data = NULL;
20: PetscFunctionReturn(PETSC_SUCCESS);
21: }
23: static PetscErrorCode DMSNESDuplicate_DMLocal(DMSNES oldsdm, DMSNES sdm)
24: {
25: PetscFunctionBegin;
26: if (sdm->data != oldsdm->data) {
27: PetscCall(PetscFree(sdm->data));
28: PetscCall(PetscNew((DMSNES_Local **)&sdm->data));
29: if (oldsdm->data) PetscCall(PetscMemcpy(sdm->data, oldsdm->data, sizeof(DMSNES_Local)));
30: }
31: PetscFunctionReturn(PETSC_SUCCESS);
32: }
34: static PetscErrorCode DMLocalSNESGetContext(DM dm, DMSNES sdm, DMSNES_Local **dmlocalsnes)
35: {
36: PetscFunctionBegin;
37: *dmlocalsnes = NULL;
38: if (!sdm->data) {
39: PetscCall(PetscNew((DMSNES_Local **)&sdm->data));
41: sdm->ops->destroy = DMSNESDestroy_DMLocal;
42: sdm->ops->duplicate = DMSNESDuplicate_DMLocal;
43: }
44: *dmlocalsnes = (DMSNES_Local *)sdm->data;
45: PetscFunctionReturn(PETSC_SUCCESS);
46: }
48: static PetscErrorCode SNESComputeObjective_DMLocal(SNES snes, Vec X, PetscReal *obj, void *ctx)
49: {
50: DMSNES_Local *dmlocalsnes = (DMSNES_Local *)ctx;
51: DM dm;
52: Vec Xloc;
53: PetscBool transform;
55: PetscFunctionBegin;
58: PetscCall(SNESGetDM(snes, &dm));
59: PetscCall(DMGetLocalVector(dm, &Xloc));
60: PetscCall(VecZeroEntries(Xloc));
61: /* Non-conforming routines needs boundary values before G2L */
62: if (dmlocalsnes->boundarylocal) PetscCall((*dmlocalsnes->boundarylocal)(dm, Xloc, dmlocalsnes->boundarylocalctx));
63: PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc));
64: PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc));
65: /* Need to reset boundary values if we transformed */
66: PetscCall(DMHasBasisTransform(dm, &transform));
67: if (transform && dmlocalsnes->boundarylocal) PetscCall((*dmlocalsnes->boundarylocal)(dm, Xloc, dmlocalsnes->boundarylocalctx));
68: CHKMEMQ;
69: PetscCall((*dmlocalsnes->objectivelocal)(dm, Xloc, obj, dmlocalsnes->objectivelocalctx));
70: CHKMEMQ;
71: PetscCall(MPIU_Allreduce(MPI_IN_PLACE, obj, 1, MPIU_REAL, MPIU_SUM, PetscObjectComm((PetscObject)snes)));
72: PetscCall(DMRestoreLocalVector(dm, &Xloc));
73: PetscFunctionReturn(PETSC_SUCCESS);
74: }
76: static PetscErrorCode SNESComputeFunction_DMLocal(SNES snes, Vec X, Vec F, void *ctx)
77: {
78: DMSNES_Local *dmlocalsnes = (DMSNES_Local *)ctx;
79: DM dm;
80: Vec Xloc, Floc;
81: PetscBool transform;
83: PetscFunctionBegin;
87: PetscCall(SNESGetDM(snes, &dm));
88: PetscCall(DMGetLocalVector(dm, &Xloc));
89: PetscCall(DMGetLocalVector(dm, &Floc));
90: PetscCall(VecZeroEntries(Xloc));
91: PetscCall(VecZeroEntries(Floc));
92: /* Non-conforming routines needs boundary values before G2L */
93: if (dmlocalsnes->boundarylocal) PetscCall((*dmlocalsnes->boundarylocal)(dm, Xloc, dmlocalsnes->boundarylocalctx));
94: PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc));
95: PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc));
96: /* Need to reset boundary values if we transformed */
97: PetscCall(DMHasBasisTransform(dm, &transform));
98: if (transform && dmlocalsnes->boundarylocal) PetscCall((*dmlocalsnes->boundarylocal)(dm, Xloc, dmlocalsnes->boundarylocalctx));
99: CHKMEMQ;
100: PetscCall((*dmlocalsnes->residuallocal)(dm, Xloc, Floc, dmlocalsnes->residuallocalctx));
101: CHKMEMQ;
102: PetscCall(VecZeroEntries(F));
103: PetscCall(DMLocalToGlobalBegin(dm, Floc, ADD_VALUES, F));
104: PetscCall(DMLocalToGlobalEnd(dm, Floc, ADD_VALUES, F));
105: PetscCall(DMRestoreLocalVector(dm, &Floc));
106: PetscCall(DMRestoreLocalVector(dm, &Xloc));
107: {
108: char name[PETSC_MAX_PATH_LEN];
109: char oldname[PETSC_MAX_PATH_LEN];
110: const char *tmp;
111: PetscInt it;
113: PetscCall(SNESGetIterationNumber(snes, &it));
114: PetscCall(PetscSNPrintf(name, PETSC_MAX_PATH_LEN, "Solution, Iterate %d", (int)it));
115: PetscCall(PetscObjectGetName((PetscObject)X, &tmp));
116: PetscCall(PetscStrncpy(oldname, tmp, PETSC_MAX_PATH_LEN - 1));
117: PetscCall(PetscObjectSetName((PetscObject)X, name));
118: PetscCall(VecViewFromOptions(X, (PetscObject)snes, "-dmsnes_solution_vec_view"));
119: PetscCall(PetscObjectSetName((PetscObject)X, oldname));
120: PetscCall(PetscSNPrintf(name, PETSC_MAX_PATH_LEN, "Residual, Iterate %d", (int)it));
121: PetscCall(PetscObjectSetName((PetscObject)F, name));
122: PetscCall(VecViewFromOptions(F, (PetscObject)snes, "-dmsnes_residual_vec_view"));
123: }
124: PetscFunctionReturn(PETSC_SUCCESS);
125: }
127: static PetscErrorCode SNESComputeJacobian_DMLocal(SNES snes, Vec X, Mat A, Mat B, void *ctx)
128: {
129: DMSNES_Local *dmlocalsnes = (DMSNES_Local *)ctx;
130: DM dm;
131: Vec Xloc;
132: PetscBool transform;
134: PetscFunctionBegin;
135: PetscCall(SNESGetDM(snes, &dm));
136: if (dmlocalsnes->jacobianlocal) {
137: PetscCall(DMGetLocalVector(dm, &Xloc));
138: PetscCall(VecZeroEntries(Xloc));
139: /* Non-conforming routines needs boundary values before G2L */
140: if (dmlocalsnes->boundarylocal) PetscCall((*dmlocalsnes->boundarylocal)(dm, Xloc, dmlocalsnes->boundarylocalctx));
141: PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, Xloc));
142: PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, Xloc));
143: /* Need to reset boundary values if we transformed */
144: PetscCall(DMHasBasisTransform(dm, &transform));
145: if (transform && dmlocalsnes->boundarylocal) PetscCall((*dmlocalsnes->boundarylocal)(dm, Xloc, dmlocalsnes->boundarylocalctx));
146: CHKMEMQ;
147: PetscCall((*dmlocalsnes->jacobianlocal)(dm, Xloc, A, B, dmlocalsnes->jacobianlocalctx));
148: CHKMEMQ;
149: PetscCall(DMRestoreLocalVector(dm, &Xloc));
150: } else {
151: MatFDColoring fdcoloring;
152: PetscCall(PetscObjectQuery((PetscObject)dm, "DMDASNES_FDCOLORING", (PetscObject *)&fdcoloring));
153: if (!fdcoloring) {
154: ISColoring coloring;
156: PetscCall(DMCreateColoring(dm, dm->coloringtype, &coloring));
157: PetscCall(MatFDColoringCreate(B, coloring, &fdcoloring));
158: PetscCall(ISColoringDestroy(&coloring));
159: switch (dm->coloringtype) {
160: case IS_COLORING_GLOBAL:
161: PetscCall(MatFDColoringSetFunction(fdcoloring, (PetscErrorCode(*)(void))SNESComputeFunction_DMLocal, dmlocalsnes));
162: break;
163: default:
164: SETERRQ(PetscObjectComm((PetscObject)snes), PETSC_ERR_SUP, "No support for coloring type '%s'", ISColoringTypes[dm->coloringtype]);
165: }
166: PetscCall(PetscObjectSetOptionsPrefix((PetscObject)fdcoloring, ((PetscObject)dm)->prefix));
167: PetscCall(MatFDColoringSetFromOptions(fdcoloring));
168: PetscCall(MatFDColoringSetUp(B, coloring, fdcoloring));
169: PetscCall(PetscObjectCompose((PetscObject)dm, "DMDASNES_FDCOLORING", (PetscObject)fdcoloring));
170: PetscCall(PetscObjectDereference((PetscObject)fdcoloring));
172: /* The following breaks an ugly reference counting loop that deserves a paragraph. MatFDColoringApply() will call
173: * VecDuplicate() with the state Vec and store inside the MatFDColoring. This Vec will duplicate the Vec, but the
174: * MatFDColoring is composed with the DM. We dereference the DM here so that the reference count will eventually
175: * drop to 0. Note the code in DMDestroy() that exits early for a negative reference count. That code path will be
176: * taken when the PetscObjectList for the Vec inside MatFDColoring is destroyed.
177: */
178: PetscCall(PetscObjectDereference((PetscObject)dm));
179: }
180: PetscCall(MatFDColoringApply(B, fdcoloring, X, snes));
181: }
182: /* This will be redundant if the user called both, but it's too common to forget. */
183: if (A != B) {
184: PetscCall(MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY));
185: PetscCall(MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY));
186: }
187: PetscFunctionReturn(PETSC_SUCCESS);
188: }
190: /*@C
191: DMSNESSetObjectiveLocal - set a local objective evaluation function. This function is called with local vector
192: containing the local vector information PLUS ghost point information. It should compute a result for all local
193: elements and `DMSNES` will automatically accumulate the overlapping values.
195: Logically Collective
197: Input Parameters:
198: + dm - `DM` to associate callback with
199: . func - local objective evaluation
200: - ctx - optional context for local residual evaluation
202: Level: advanced
204: .seealso: `DMSNESSetFunctionLocal()`, `DMSNESSetJacobianLocal()`
205: @*/
206: PetscErrorCode DMSNESSetObjectiveLocal(DM dm, PetscErrorCode (*func)(DM, Vec, PetscReal *, void *), void *ctx)
207: {
208: DMSNES sdm;
209: DMSNES_Local *dmlocalsnes;
211: PetscFunctionBegin;
213: PetscCall(DMGetDMSNESWrite(dm, &sdm));
214: PetscCall(DMLocalSNESGetContext(dm, sdm, &dmlocalsnes));
216: dmlocalsnes->objectivelocal = func;
217: dmlocalsnes->objectivelocalctx = ctx;
219: PetscCall(DMSNESSetObjective(dm, SNESComputeObjective_DMLocal, dmlocalsnes));
220: PetscFunctionReturn(PETSC_SUCCESS);
221: }
223: /*@C
224: DMSNESSetFunctionLocal - set a local residual evaluation function. This function is called with local vector
225: containing the local vector information PLUS ghost point information. It should compute a result for all local
226: elements and `DMSNES` will automatically accumulate the overlapping values.
228: Logically Collective
230: Input Parameters:
231: + dm - `DM` to associate callback with
232: . func - local residual evaluation
233: - ctx - optional context for local residual evaluation
235: Calling sequence of `func`:
236: + dm - `DM` for the function
237: . x - vector to state at which to evaluate residual
238: . f - vector to hold the function evaluation
239: - ctx - optional context passed above
241: Level: advanced
243: .seealso: [](ch_snes), `DMSNESSetFunction()`, `DMSNESSetJacobianLocal()`
244: @*/
245: PetscErrorCode DMSNESSetFunctionLocal(DM dm, PetscErrorCode (*func)(DM dm, Vec x, Vec f, void *ctx), void *ctx)
246: {
247: DMSNES sdm;
248: DMSNES_Local *dmlocalsnes;
250: PetscFunctionBegin;
252: PetscCall(DMGetDMSNESWrite(dm, &sdm));
253: PetscCall(DMLocalSNESGetContext(dm, sdm, &dmlocalsnes));
255: dmlocalsnes->residuallocal = func;
256: dmlocalsnes->residuallocalctx = ctx;
258: PetscCall(DMSNESSetFunction(dm, SNESComputeFunction_DMLocal, dmlocalsnes));
259: if (!sdm->ops->computejacobian) { /* Call us for the Jacobian too, can be overridden by the user. */
260: PetscCall(DMSNESSetJacobian(dm, SNESComputeJacobian_DMLocal, dmlocalsnes));
261: }
262: PetscFunctionReturn(PETSC_SUCCESS);
263: }
265: /*@C
266: DMSNESSetBoundaryLocal - set a function to insert, for example, essential boundary conditions into a ghosted solution vector
268: Logically Collective
270: Input Parameters:
271: + dm - `DM` to associate callback with
272: . func - local boundary value evaluation
273: - ctx - optional context for local boundary value evaluation
275: Calling sequence of `func`:
276: + dm - the `DM` context
277: . X - ghosted solution vector, appropriate locations (such as essential boundary condition nodes) should be filled
278: - ctx - option context passed in `DMSNESSetBoundaryLocal()`
280: Level: advanced
282: .seealso: [](ch_snes), `DMSNESSetObjectiveLocal()`, `DMSNESSetFunctionLocal()`, `DMSNESSetJacobianLocal()`
283: @*/
284: PetscErrorCode DMSNESSetBoundaryLocal(DM dm, PetscErrorCode (*func)(DM dm, Vec X, void *ctx), void *ctx)
285: {
286: DMSNES sdm;
287: DMSNES_Local *dmlocalsnes;
289: PetscFunctionBegin;
291: PetscCall(DMGetDMSNESWrite(dm, &sdm));
292: PetscCall(DMLocalSNESGetContext(dm, sdm, &dmlocalsnes));
294: dmlocalsnes->boundarylocal = func;
295: dmlocalsnes->boundarylocalctx = ctx;
296: PetscFunctionReturn(PETSC_SUCCESS);
297: }
299: /*@C
300: DMSNESSetJacobianLocal - set a local Jacobian evaluation function
302: Logically Collective
304: Input Parameters:
305: + dm - `DM` to associate callback with
306: . func - local Jacobian evaluation
307: - ctx - optional context for local Jacobian evaluation
309: Calling sequence of `func`:
310: + dm - the `DM` context
311: . X - current solution vector (ghosted or not?)
312: . J - the Jacobian
313: . Jp - approximate Jacobian used to compute the preconditioner, often `J`
314: - ctx - a user provided context
316: Level: advanced
318: .seealso: [](ch_snes), `DMSNESSetObjectiveLocal()`, `DMSNESSetFunctionLocal()`, `DMSNESSetBoundaryLocal()`
319: @*/
320: PetscErrorCode DMSNESSetJacobianLocal(DM dm, PetscErrorCode (*func)(DM dm, Vec X, Mat J, Mat Jp, void *ctx), void *ctx)
321: {
322: DMSNES sdm;
323: DMSNES_Local *dmlocalsnes;
325: PetscFunctionBegin;
327: PetscCall(DMGetDMSNESWrite(dm, &sdm));
328: PetscCall(DMLocalSNESGetContext(dm, sdm, &dmlocalsnes));
330: dmlocalsnes->jacobianlocal = func;
331: dmlocalsnes->jacobianlocalctx = ctx;
333: PetscCall(DMSNESSetJacobian(dm, SNESComputeJacobian_DMLocal, dmlocalsnes));
334: PetscFunctionReturn(PETSC_SUCCESS);
335: }
337: /*@C
338: DMSNESGetObjectiveLocal - get the local objective evaluation function information set with `DMSNESSetObjectiveLocal()`.
340: Not Collective
342: Input Parameter:
343: . dm - `DM` with the associated callback
345: Output Parameters:
346: + func - local objective evaluation
347: - ctx - context for local residual evaluation
349: Level: beginner
351: .seealso: `DMSNESSetObjective()`, `DMSNESSetObjectiveLocal()`, `DMSNESSetFunctionLocal()`
352: @*/
353: PetscErrorCode DMSNESGetObjectiveLocal(DM dm, PetscErrorCode (**func)(DM, Vec, PetscReal *, void *), void **ctx)
354: {
355: DMSNES sdm;
356: DMSNES_Local *dmlocalsnes;
358: PetscFunctionBegin;
360: PetscCall(DMGetDMSNES(dm, &sdm));
361: PetscCall(DMLocalSNESGetContext(dm, sdm, &dmlocalsnes));
362: if (func) *func = dmlocalsnes->objectivelocal;
363: if (ctx) *ctx = dmlocalsnes->objectivelocalctx;
364: PetscFunctionReturn(PETSC_SUCCESS);
365: }
367: /*@C
368: DMSNESGetFunctionLocal - get the local residual evaluation function information set with `DMSNESSetFunctionLocal()`.
370: Not Collective
372: Input Parameter:
373: . dm - `DM` with the associated callback
375: Output Parameters:
376: + func - local residual evaluation
377: - ctx - context for local residual evaluation
379: Level: beginner
381: .seealso: [](ch_snes), `DMSNESSetFunction()`, `DMSNESSetFunctionLocal()`, `DMSNESSetJacobianLocal()`
382: @*/
383: PetscErrorCode DMSNESGetFunctionLocal(DM dm, PetscErrorCode (**func)(DM, Vec, Vec, void *), void **ctx)
384: {
385: DMSNES sdm;
386: DMSNES_Local *dmlocalsnes;
388: PetscFunctionBegin;
390: PetscCall(DMGetDMSNES(dm, &sdm));
391: PetscCall(DMLocalSNESGetContext(dm, sdm, &dmlocalsnes));
392: if (func) *func = dmlocalsnes->residuallocal;
393: if (ctx) *ctx = dmlocalsnes->residuallocalctx;
394: PetscFunctionReturn(PETSC_SUCCESS);
395: }
397: /*@C
398: DMSNESGetBoundaryLocal - get the local boundary value function set with `DMSNESSetBoundaryLocal()`.
400: Not Collective
402: Input Parameter:
403: . dm - `DM` with the associated callback
405: Output Parameters:
406: + func - local boundary value evaluation
407: - ctx - context for local boundary value evaluation
409: Level: intermediate
411: .seealso: [](ch_snes), `DMSNESSetFunctionLocal()`, `DMSNESSetBoundaryLocal()`, `DMSNESSetJacobianLocal()`
412: @*/
413: PetscErrorCode DMSNESGetBoundaryLocal(DM dm, PetscErrorCode (**func)(DM, Vec, void *), void **ctx)
414: {
415: DMSNES sdm;
416: DMSNES_Local *dmlocalsnes;
418: PetscFunctionBegin;
420: PetscCall(DMGetDMSNES(dm, &sdm));
421: PetscCall(DMLocalSNESGetContext(dm, sdm, &dmlocalsnes));
422: if (func) *func = dmlocalsnes->boundarylocal;
423: if (ctx) *ctx = dmlocalsnes->boundarylocalctx;
424: PetscFunctionReturn(PETSC_SUCCESS);
425: }
427: /*@C
428: DMSNESGetJacobianLocal - the local Jacobian evaluation function set with `DMSNESSetJacobianLocal()`.
430: Logically Collective
432: Input Parameter:
433: . dm - `DM` with the associated callback
435: Output Parameters:
436: + func - local Jacobian evaluation
437: - ctx - context for local Jacobian evaluation
439: Level: beginner
441: .seealso: [](ch_snes), `DMSNESSetJacobianLocal()`, `DMSNESSetJacobian()`
442: @*/
443: PetscErrorCode DMSNESGetJacobianLocal(DM dm, PetscErrorCode (**func)(DM, Vec, Mat, Mat, void *), void **ctx)
444: {
445: DMSNES sdm;
446: DMSNES_Local *dmlocalsnes;
448: PetscFunctionBegin;
450: PetscCall(DMGetDMSNES(dm, &sdm));
451: PetscCall(DMLocalSNESGetContext(dm, sdm, &dmlocalsnes));
452: if (func) *func = dmlocalsnes->jacobianlocal;
453: if (ctx) *ctx = dmlocalsnes->jacobianlocalctx;
454: PetscFunctionReturn(PETSC_SUCCESS);
455: }