Actual source code: dmlocalts.c
1: #include <petsc/private/dmimpl.h>
2: #include <petsc/private/tsimpl.h>
4: typedef struct {
5: PetscErrorCode (*boundarylocal)(DM, PetscReal, Vec, Vec, void *);
6: PetscErrorCode (*ifunctionlocal)(DM, PetscReal, Vec, Vec, Vec, void *);
7: PetscErrorCode (*ijacobianlocal)(DM, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *);
8: PetscErrorCode (*rhsfunctionlocal)(DM, PetscReal, Vec, Vec, void *);
9: void *boundarylocalctx;
10: void *ifunctionlocalctx;
11: void *ijacobianlocalctx;
12: void *rhsfunctionlocalctx;
13: Vec lumpedmassinv;
14: Mat mass;
15: KSP kspmass;
16: } DMTS_Local;
18: static PetscErrorCode DMTSDestroy_DMLocal(DMTS tdm)
19: {
20: PetscFunctionBegin;
21: PetscCall(PetscFree(tdm->data));
22: PetscFunctionReturn(PETSC_SUCCESS);
23: }
25: static PetscErrorCode DMTSDuplicate_DMLocal(DMTS oldtdm, DMTS tdm)
26: {
27: PetscFunctionBegin;
28: PetscCall(PetscNew((DMTS_Local **)&tdm->data));
29: if (oldtdm->data) PetscCall(PetscMemcpy(tdm->data, oldtdm->data, sizeof(DMTS_Local)));
30: PetscFunctionReturn(PETSC_SUCCESS);
31: }
33: static PetscErrorCode DMLocalTSGetContext(DM dm, DMTS tdm, DMTS_Local **dmlocalts)
34: {
35: PetscFunctionBegin;
36: *dmlocalts = NULL;
37: if (!tdm->data) {
38: PetscCall(PetscNew((DMTS_Local **)&tdm->data));
40: tdm->ops->destroy = DMTSDestroy_DMLocal;
41: tdm->ops->duplicate = DMTSDuplicate_DMLocal;
42: }
43: *dmlocalts = (DMTS_Local *)tdm->data;
44: PetscFunctionReturn(PETSC_SUCCESS);
45: }
47: static PetscErrorCode TSComputeIFunction_DMLocal(TS ts, PetscReal time, Vec X, Vec X_t, Vec F, PetscCtx ctx)
48: {
49: DM dm;
50: Vec locX, locX_t, locF;
51: DMTS_Local *dmlocalts = (DMTS_Local *)ctx;
53: PetscFunctionBegin;
58: PetscCall(TSGetDM(ts, &dm));
59: PetscCall(DMGetLocalVector(dm, &locX));
60: PetscCall(DMGetLocalVector(dm, &locX_t));
61: PetscCall(DMGetLocalVector(dm, &locF));
62: PetscCall(VecZeroEntries(locX));
63: PetscCall(VecZeroEntries(locX_t));
64: if (dmlocalts->boundarylocal) PetscCall((*dmlocalts->boundarylocal)(dm, time, locX, locX_t, dmlocalts->boundarylocalctx));
65: PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX));
66: PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX));
67: PetscCall(DMGlobalToLocalBegin(dm, X_t, INSERT_VALUES, locX_t));
68: PetscCall(DMGlobalToLocalEnd(dm, X_t, INSERT_VALUES, locX_t));
69: PetscCall(VecZeroEntries(locF));
70: CHKMEMQ;
71: PetscCall((*dmlocalts->ifunctionlocal)(dm, time, locX, locX_t, locF, dmlocalts->ifunctionlocalctx));
72: CHKMEMQ;
73: PetscCall(VecZeroEntries(F));
74: PetscCall(DMLocalToGlobalBegin(dm, locF, ADD_VALUES, F));
75: PetscCall(DMLocalToGlobalEnd(dm, locF, ADD_VALUES, F));
76: PetscCall(DMRestoreLocalVector(dm, &locX));
77: PetscCall(DMRestoreLocalVector(dm, &locX_t));
78: PetscCall(DMRestoreLocalVector(dm, &locF));
80: /* remove nullspace from residual */
81: {
82: MatNullSpace nullsp;
83: PetscCall(PetscObjectQuery((PetscObject)dm, "__dmtsnullspace", (PetscObject *)&nullsp));
84: if (nullsp) PetscCall(MatNullSpaceRemove(nullsp, F));
85: }
86: PetscFunctionReturn(PETSC_SUCCESS);
87: }
89: static PetscErrorCode TSComputeRHSFunction_DMLocal(TS ts, PetscReal time, Vec X, Vec F, PetscCtx ctx)
90: {
91: DM dm;
92: Vec locX, locF;
93: DMTS_Local *dmlocalts = (DMTS_Local *)ctx;
95: PetscFunctionBegin;
99: PetscCall(TSGetDM(ts, &dm));
100: PetscCall(DMGetLocalVector(dm, &locX));
101: PetscCall(DMGetLocalVector(dm, &locF));
102: PetscCall(VecZeroEntries(locX));
103: if (dmlocalts->boundarylocal) PetscCall((*dmlocalts->boundarylocal)(dm, time, locX, NULL, dmlocalts->boundarylocalctx));
104: PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX));
105: PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX));
106: PetscCall(VecZeroEntries(locF));
107: CHKMEMQ;
108: PetscCall((*dmlocalts->rhsfunctionlocal)(dm, time, locX, locF, dmlocalts->rhsfunctionlocalctx));
109: CHKMEMQ;
110: PetscCall(VecZeroEntries(F));
111: PetscCall(DMLocalToGlobalBegin(dm, locF, ADD_VALUES, F));
112: PetscCall(DMLocalToGlobalEnd(dm, locF, ADD_VALUES, F));
113: if (dmlocalts->lumpedmassinv) PetscCall(VecPointwiseMult(F, dmlocalts->lumpedmassinv, F));
114: else if (dmlocalts->kspmass) {
115: Vec tmp;
117: PetscCall(DMGetGlobalVector(dm, &tmp));
118: PetscCall(KSPSolve(dmlocalts->kspmass, F, tmp));
119: PetscCall(VecCopy(tmp, F));
120: PetscCall(DMRestoreGlobalVector(dm, &tmp));
121: }
122: PetscCall(DMRestoreLocalVector(dm, &locX));
123: PetscCall(DMRestoreLocalVector(dm, &locF));
124: PetscFunctionReturn(PETSC_SUCCESS);
125: }
127: static PetscErrorCode TSComputeIJacobian_DMLocal(TS ts, PetscReal time, Vec X, Vec X_t, PetscReal a, Mat A, Mat B, PetscCtx ctx)
128: {
129: DM dm;
130: Vec locX, locX_t;
131: DMTS_Local *dmlocalts = (DMTS_Local *)ctx;
133: PetscFunctionBegin;
134: PetscCall(TSGetDM(ts, &dm));
135: if (dmlocalts->ijacobianlocal) {
136: PetscCall(DMGetLocalVector(dm, &locX));
137: PetscCall(DMGetLocalVector(dm, &locX_t));
138: PetscCall(VecZeroEntries(locX));
139: PetscCall(VecZeroEntries(locX_t));
140: if (dmlocalts->boundarylocal) PetscCall((*dmlocalts->boundarylocal)(dm, time, locX, locX_t, dmlocalts->boundarylocalctx));
141: PetscCall(DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX));
142: PetscCall(DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX));
143: PetscCall(DMGlobalToLocalBegin(dm, X_t, INSERT_VALUES, locX_t));
144: PetscCall(DMGlobalToLocalEnd(dm, X_t, INSERT_VALUES, locX_t));
145: CHKMEMQ;
146: PetscCall((*dmlocalts->ijacobianlocal)(dm, time, locX, locX_t, a, A, B, dmlocalts->ijacobianlocalctx));
147: CHKMEMQ;
148: PetscCall(DMRestoreLocalVector(dm, &locX));
149: PetscCall(DMRestoreLocalVector(dm, &locX_t));
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, (MatFDColoringFn *)(PetscVoidFn *)TSComputeIFunction_DMLocal, dmlocalts));
162: break;
163: default:
164: SETERRQ(PetscObjectComm((PetscObject)ts), 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, ts));
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: DMTSSetBoundaryLocal - set the function for essential boundary data for a local implicit function evaluation.
193: Logically Collective
195: Input Parameters:
196: + dm - `DM` to associate callback with
197: . func - local function evaluation
198: - ctx - context for function evaluation
200: Calling sequence of `func`:
201: + dm - the `DM`
202: . t - the current time
203: . u - the current solution
204: . f - output, the computed right hand side function
205: - ctx - the application context for the function
207: Level: intermediate
209: Notes:
210: `func` should set the essential boundary data for the local portion of the solution, as
211: well its time derivative (if it is not `NULL`).
213: Vectors are initialized to zero before this function, so it is only needed for non
214: homogeneous data.
216: This function is somewhat optional: boundary data could potentially be inserted by a function
217: passed to `DMTSSetIFunctionLocal()`. The use case for this function is for discretizations
218: with constraints (see `DMGetDefaultConstraints()`): this function inserts boundary values
219: before constraint interpolation.
221: .seealso: [](ch_ts), `DM`, `TS`, `DMTSSetIFunction()`, `DMTSSetIJacobianLocal()`
222: @*/
223: PetscErrorCode DMTSSetBoundaryLocal(DM dm, PetscErrorCode (*func)(DM dm, PetscReal t, Vec u, Vec f, PetscCtx ctx), PetscCtx ctx)
224: {
225: DMTS tdm;
226: DMTS_Local *dmlocalts;
228: PetscFunctionBegin;
230: PetscCall(DMGetDMTSWrite(dm, &tdm));
231: PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts));
233: dmlocalts->boundarylocal = func;
234: dmlocalts->boundarylocalctx = ctx;
235: PetscFunctionReturn(PETSC_SUCCESS);
236: }
238: /*@C
239: DMTSGetIFunctionLocal - get the local implicit function evaluation function. This function is called with local vector
240: containing the local vector information PLUS ghost point information. It should compute a result for all local
241: elements and `DM` will automatically accumulate the overlapping values.
243: Logically Collective
245: Input Parameter:
246: . dm - `DM` to associate callback with
248: Output Parameters:
249: + func - local function evaluation
250: - ctx - context for function evaluation
252: Calling sequence of `func`:
253: + dm - the `DM`
254: . t - the current time
255: . u - the current solution
256: . udot - the derivative of `u`
257: . F - output, the computed implicit function
258: - ctx - the application context for the function
260: Level: beginner
262: .seealso: [](ch_ts), `DM`, `DMTSSetIFunctionLocal()`, `DMTSSetIFunction()`, `DMTSSetIJacobianLocal()`
263: @*/
264: PetscErrorCode DMTSGetIFunctionLocal(DM dm, PetscErrorCode (**func)(DM dm, PetscReal t, Vec u, Vec udot, Vec F, PetscCtx ctx), PetscCtxRt ctx)
265: {
266: DMTS tdm;
267: DMTS_Local *dmlocalts;
269: PetscFunctionBegin;
271: PetscCall(DMGetDMTS(dm, &tdm));
272: PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts));
273: if (func) {
274: PetscAssertPointer(func, 2);
275: *func = dmlocalts->ifunctionlocal;
276: }
277: if (ctx) {
278: PetscAssertPointer(ctx, 3);
279: *(void **)ctx = dmlocalts->ifunctionlocalctx;
280: }
281: PetscFunctionReturn(PETSC_SUCCESS);
282: }
284: /*@C
285: DMTSSetIFunctionLocal - set a local implicit function evaluation function. This function is called with local vector
286: containing the local vector information PLUS ghost point information. It should compute a result for all local
287: elements and `DM` will automatically accumulate the overlapping values.
289: Logically Collective
291: Input Parameters:
292: + dm - `DM` to associate callback with
293: . func - local function evaluation
294: - ctx - context for function evaluation
296: Calling sequence of `func`:
297: + dm - the `DM`
298: . t - the current time
299: . u - the current solution
300: . udot - the derivative of `u`
301: . F - output, the computed implicit function
302: - ctx - the application context for the function
304: Level: beginner
306: .seealso: [](ch_ts), `DM`, `DMTSGetIFunctionLocal()`, `DMTSSetIFunction()`, `DMTSSetIJacobianLocal()`
307: @*/
308: PetscErrorCode DMTSSetIFunctionLocal(DM dm, PetscErrorCode (*func)(DM dm, PetscReal t, Vec u, Vec udot, Vec F, PetscCtx ctx), PetscCtx ctx)
309: {
310: DMTS tdm;
311: DMTS_Local *dmlocalts;
313: PetscFunctionBegin;
315: PetscCall(DMGetDMTSWrite(dm, &tdm));
316: PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts));
318: dmlocalts->ifunctionlocal = func;
319: dmlocalts->ifunctionlocalctx = ctx;
321: PetscCall(DMTSSetIFunction(dm, TSComputeIFunction_DMLocal, dmlocalts));
322: if (!tdm->ops->ijacobian) { /* Call us for the Jacobian too, can be overridden by the user. */
323: PetscCall(DMTSSetIJacobian(dm, TSComputeIJacobian_DMLocal, dmlocalts));
324: }
325: PetscFunctionReturn(PETSC_SUCCESS);
326: }
328: /*@C
329: DMTSGetIJacobianLocal - get a local Jacobian evaluation function
331: Logically Collective
333: Input Parameter:
334: . dm - `DM` to associate callback with
336: Output Parameters:
337: + func - local Jacobian evaluation
338: - ctx - optional context for local Jacobian evaluation
340: Calling sequence of `func`:
341: + dm - the `DM`
342: . t - the current time
343: . u - the current solution
344: . udot - the derivative of `u`
345: . shift - the shift factoring arising from the implicit time-step
346: . J - output, the Jacobian
347: . Jpre - output, matrix from which to compute the preconditioner for `J`, often the same as `J`
348: - ctx - the application context for the function
350: Level: beginner
352: .seealso: [](ch_ts), `DM`, `DMTSSetIJacobianLocal()`, `DMTSSetIFunctionLocal()`, `DMTSSetIJacobian()`, `DMTSSetIFunction()`
353: @*/
354: PetscErrorCode DMTSGetIJacobianLocal(DM dm, PetscErrorCode (**func)(DM dm, PetscReal t, Vec u, Vec udot, PetscReal shift, Mat J, Mat Jpre, PetscCtx ctx), PetscCtxRt ctx)
355: {
356: DMTS tdm;
357: DMTS_Local *dmlocalts;
359: PetscFunctionBegin;
361: PetscCall(DMGetDMTS(dm, &tdm));
362: PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts));
363: if (func) {
364: PetscAssertPointer(func, 2);
365: *func = dmlocalts->ijacobianlocal;
366: }
367: if (ctx) {
368: PetscAssertPointer(ctx, 3);
369: *(void **)ctx = dmlocalts->ijacobianlocalctx;
370: }
371: PetscFunctionReturn(PETSC_SUCCESS);
372: }
374: /*@C
375: DMTSSetIJacobianLocal - set a local Jacobian evaluation function
377: Logically Collective
379: Input Parameters:
380: + dm - `DM` to associate callback with
381: . func - local Jacobian evaluation
382: - ctx - optional context for local Jacobian evaluation
384: Calling sequence of `func`:
385: + dm - the `DM`
386: . t - the current time
387: . u - the current solution
388: . udot - the derivative of `u`
389: . shift - the shift factoring arising from the implicit time-step
390: . J - output, the Jacobian
391: . Jpre - output, matrix from which to compute the preconditioner for `J`, often the same as `J`
392: - ctx - the application context for the function
394: Level: beginner
396: .seealso: [](ch_ts), `DM`, `DMTSGetIJacobianLocal()`, `DMTSSetIFunctionLocal()`, `DMTSSetIJacobian()`, `DMTSSetIFunction()`
397: @*/
398: PetscErrorCode DMTSSetIJacobianLocal(DM dm, PetscErrorCode (*func)(DM dm, PetscReal t, Vec u, Vec udot, PetscReal shift, Mat J, Mat Jpre, PetscCtx ctx), PetscCtx ctx)
399: {
400: DMTS tdm;
401: DMTS_Local *dmlocalts;
403: PetscFunctionBegin;
405: PetscCall(DMGetDMTSWrite(dm, &tdm));
406: PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts));
408: dmlocalts->ijacobianlocal = func;
409: dmlocalts->ijacobianlocalctx = ctx;
411: PetscCall(DMTSSetIJacobian(dm, TSComputeIJacobian_DMLocal, dmlocalts));
412: PetscFunctionReturn(PETSC_SUCCESS);
413: }
415: /*@C
416: DMTSGetRHSFunctionLocal - get a local rhs function evaluation function. This function is called with local vector
417: containing the local vector information PLUS ghost point information. It should compute a result for all local
418: elements and `DM` will automatically accumulate the overlapping values.
420: Logically Collective
422: Input Parameter:
423: . dm - `DM` to associate callback with
425: Output Parameters:
426: + func - local function evaluation
427: - ctx - context for function evaluation
429: Calling sequence of `func`:
430: + dm - the `DM`
431: . t - the current time
432: . u - the current solution
433: . udot - output, the evaluated right hand side
434: - ctx - the application context for the function
436: Level: beginner
438: .seealso: [](ch_ts), `DM`, `DMTSSetRHSFunctionLocal()`, `DMTSSetRHSFunction()`, `DMTSSetIFunction()`, `DMTSSetIJacobianLocal()`
439: @*/
440: PetscErrorCode DMTSGetRHSFunctionLocal(DM dm, PetscErrorCode (**func)(DM dm, PetscReal t, Vec u, Vec udot, PetscCtx ctx), PetscCtxRt ctx)
441: {
442: DMTS tdm;
443: DMTS_Local *dmlocalts;
445: PetscFunctionBegin;
447: PetscCall(DMGetDMTS(dm, &tdm));
448: PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts));
449: if (func) {
450: PetscAssertPointer(func, 2);
451: *func = dmlocalts->rhsfunctionlocal;
452: }
453: if (ctx) {
454: PetscAssertPointer(ctx, 3);
455: *(void **)ctx = dmlocalts->rhsfunctionlocalctx;
456: }
457: PetscFunctionReturn(PETSC_SUCCESS);
458: }
460: /*@C
461: DMTSSetRHSFunctionLocal - set a local rhs function evaluation function. This function is called with local vector
462: containing the local vector information PLUS ghost point information. It should compute a result for all local
463: elements and `DM` will automatically accumulate the overlapping values.
465: Logically Collective
467: Input Parameters:
468: + dm - `DM` to associate callback with
469: . func - local function evaluation
470: - ctx - context for function evaluation
472: Calling sequence of `func`:
473: + dm - the `DM`
474: . t - the current time
475: . u - the current solution
476: . f - output, the evaluated right hand side
477: - ctx - the application context for the function
479: Level: beginner
481: .seealso: [](ch_ts), `DM`, `DMTSGetRHSFunctionLocal()`, `DMTSSetRHSFunction()`, `DMTSSetIFunction()`, `DMTSSetIJacobianLocal()`
482: @*/
483: PetscErrorCode DMTSSetRHSFunctionLocal(DM dm, PetscErrorCode (*func)(DM dm, PetscReal t, Vec u, Vec f, PetscCtx ctx), PetscCtx ctx)
484: {
485: DMTS tdm;
486: DMTS_Local *dmlocalts;
488: PetscFunctionBegin;
490: PetscCall(DMGetDMTSWrite(dm, &tdm));
491: PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts));
493: dmlocalts->rhsfunctionlocal = func;
494: dmlocalts->rhsfunctionlocalctx = ctx;
496: PetscCall(DMTSSetRHSFunction(dm, TSComputeRHSFunction_DMLocal, dmlocalts));
497: PetscFunctionReturn(PETSC_SUCCESS);
498: }
500: /*@
501: DMTSCreateRHSMassMatrix - This creates the mass matrix associated with the given `DM`, and a solver to invert it, and stores them in the `DM` context.
503: Collective
505: Input Parameter:
506: . dm - `DM` providing the mass matrix
508: Level: developer
510: Note:
511: The idea here is that an explicit system can be given a mass matrix, based on the `DM`, which is inverted on the RHS at each step.
513: .seealso: [](ch_ts), `DM`, `DMTSCreateRHSMassMatrixLumped()`, `DMTSDestroyRHSMassMatrix()`, `DMCreateMassMatrix()`, `DMTS`
514: @*/
515: PetscErrorCode DMTSCreateRHSMassMatrix(DM dm)
516: {
517: DMTS tdm;
518: DMTS_Local *dmlocalts;
519: const char *prefix;
521: PetscFunctionBegin;
523: PetscCall(DMGetDMTSWrite(dm, &tdm));
524: PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts));
525: PetscCall(DMCreateMassMatrix(dm, dm, &dmlocalts->mass));
526: PetscCall(KSPCreate(PetscObjectComm((PetscObject)dm), &dmlocalts->kspmass));
527: PetscCall(PetscObjectGetOptionsPrefix((PetscObject)dm, &prefix));
528: PetscCall(KSPSetOptionsPrefix(dmlocalts->kspmass, prefix));
529: PetscCall(KSPAppendOptionsPrefix(dmlocalts->kspmass, "mass_"));
530: PetscCall(KSPSetFromOptions(dmlocalts->kspmass));
531: PetscCall(KSPSetOperators(dmlocalts->kspmass, dmlocalts->mass, dmlocalts->mass));
532: PetscFunctionReturn(PETSC_SUCCESS);
533: }
535: /*@
536: DMTSCreateRHSMassMatrixLumped - This creates the lumped mass matrix associated with the given `DM`, and a solver to invert it, and stores them in the `DM` context.
538: Collective
540: Input Parameter:
541: . dm - `DM` providing the mass matrix
543: Level: developer
545: Note:
546: The idea here is that an explicit system can be given a mass matrix, based on the `DM`, which is inverted on the RHS at each step.
547: Since the matrix is lumped, inversion is trivial.
549: .seealso: [](ch_ts), `DM`, `DMTSCreateRHSMassMatrix()`, `DMTSDestroyRHSMassMatrix()`, `DMCreateMassMatrix()`, `DMTS`
550: @*/
551: PetscErrorCode DMTSCreateRHSMassMatrixLumped(DM dm)
552: {
553: DMTS tdm;
554: DMTS_Local *dmlocalts;
556: PetscFunctionBegin;
558: PetscCall(DMGetDMTSWrite(dm, &tdm));
559: PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts));
560: PetscCall(DMCreateMassMatrixLumped(dm, NULL, &dmlocalts->lumpedmassinv));
561: PetscCall(VecReciprocal(dmlocalts->lumpedmassinv));
562: PetscCall(VecViewFromOptions(dmlocalts->lumpedmassinv, NULL, "-lumped_mass_inv_view"));
563: PetscFunctionReturn(PETSC_SUCCESS);
564: }
566: /*@
567: DMTSDestroyRHSMassMatrix - Destroys the mass matrix and solver stored in the `DM` context, if they exist.
569: Logically Collective
571: Input Parameter:
572: . dm - `DM` providing the mass matrix
574: Level: developer
576: .seealso: [](ch_ts), `DM`, `DMTSCreateRHSMassMatrixLumped()`, `DMCreateMassMatrix()`, `DMTS`
577: @*/
578: PetscErrorCode DMTSDestroyRHSMassMatrix(DM dm)
579: {
580: DMTS tdm;
581: DMTS_Local *dmlocalts;
583: PetscFunctionBegin;
585: PetscCall(DMGetDMTSWrite(dm, &tdm));
586: PetscCall(DMLocalTSGetContext(dm, tdm, &dmlocalts));
587: PetscCall(VecDestroy(&dmlocalts->lumpedmassinv));
588: PetscCall(MatDestroy(&dmlocalts->mass));
589: PetscCall(KSPDestroy(&dmlocalts->kspmass));
590: PetscFunctionReturn(PETSC_SUCCESS);
591: }