Actual source code: dmlocalts.c
petsc-3.6.4 2016-04-12
1: #include <petsc/private/dmimpl.h>
2: #include <petsc/private/tsimpl.h> /*I "petscts.h" I*/
4: typedef struct {
5: PetscErrorCode (*ifunctionlocal)(DM,PetscReal,Vec,Vec,Vec,void*);
6: PetscErrorCode (*ijacobianlocal)(DM,PetscReal,Vec,Vec,PetscReal,Mat,Mat,void*);
7: PetscErrorCode (*rhsfunctionlocal)(DM,PetscReal,Vec,Vec,void*);
8: void *ifunctionlocalctx;
9: void *ijacobianlocalctx;
10: void *rhsfunctionlocalctx;
11: } DMTS_Local;
15: static PetscErrorCode DMTSDestroy_DMLocal(DMTS tdm)
16: {
20: PetscFree(tdm->data);
21: return(0);
22: }
26: static PetscErrorCode DMTSDuplicate_DMLocal(DMTS oldtdm, DMTS tdm)
27: {
31: PetscNewLog(tdm, (DMTS_Local **) &tdm->data);
32: if (oldtdm->data) {PetscMemcpy(tdm->data, oldtdm->data, sizeof(DMTS_Local));}
33: return(0);
34: }
38: static PetscErrorCode DMLocalTSGetContext(DM dm, DMTS tdm, DMTS_Local **dmlocalts)
39: {
43: *dmlocalts = NULL;
44: if (!tdm->data) {
45: PetscNewLog(dm, (DMTS_Local **) &tdm->data);
47: tdm->ops->destroy = DMTSDestroy_DMLocal;
48: tdm->ops->duplicate = DMTSDuplicate_DMLocal;
49: }
50: *dmlocalts = (DMTS_Local *) tdm->data;
51: return(0);
52: }
56: static PetscErrorCode TSComputeIFunction_DMLocal(TS ts, PetscReal time, Vec X, Vec X_t, Vec F, void *ctx)
57: {
58: DM dm;
59: Vec locX, locX_t, locF;
60: DMTS_Local *dmlocalts = (DMTS_Local *) ctx;
68: TSGetDM(ts, &dm);
69: DMGetLocalVector(dm, &locX);
70: DMGetLocalVector(dm, &locX_t);
71: DMGetLocalVector(dm, &locF);
72: VecZeroEntries(locX);
73: VecZeroEntries(locX_t);
74: DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX);
75: DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX);
76: DMGlobalToLocalBegin(dm, X_t, INSERT_VALUES, locX_t);
77: DMGlobalToLocalEnd(dm, X_t, INSERT_VALUES, locX_t);
78: VecZeroEntries(locF);
79: CHKMEMQ;
80: (*dmlocalts->ifunctionlocal)(dm, time, locX, locX_t, locF, dmlocalts->ifunctionlocalctx);
81: CHKMEMQ;
82: VecZeroEntries(F);
83: DMLocalToGlobalBegin(dm, locF, ADD_VALUES, F);
84: DMLocalToGlobalEnd(dm, locF, ADD_VALUES, F);
85: DMRestoreLocalVector(dm, &locX);
86: DMRestoreLocalVector(dm, &locX_t);
87: DMRestoreLocalVector(dm, &locF);
88: return(0);
89: }
93: static PetscErrorCode TSComputeRHSFunction_DMLocal(TS ts, PetscReal time, Vec X, Vec F, void *ctx)
94: {
95: DM dm;
96: Vec locX;
97: DMTS_Local *dmlocalts = (DMTS_Local *) ctx;
104: TSGetDM(ts, &dm);
105: DMGetLocalVector(dm, &locX);
106: VecZeroEntries(locX);
107: DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX);
108: DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX);
109: VecZeroEntries(F);
110: CHKMEMQ;
111: (*dmlocalts->rhsfunctionlocal)(dm, time, locX, F, dmlocalts->rhsfunctionlocalctx);
112: CHKMEMQ;
113: DMRestoreLocalVector(dm, &locX);
114: return(0);
115: }
119: static PetscErrorCode TSComputeIJacobian_DMLocal(TS ts, PetscReal time, Vec X, Vec X_t, PetscReal a, Mat A, Mat B, void *ctx)
120: {
121: DM dm;
122: Vec locX, locX_t;
123: DMTS_Local *dmlocalts = (DMTS_Local *) ctx;
127: TSGetDM(ts, &dm);
128: if (dmlocalts->ijacobianlocal) {
129: DMGetLocalVector(dm, &locX);
130: DMGetLocalVector(dm, &locX_t);
131: VecZeroEntries(locX);
132: VecZeroEntries(locX_t);
133: DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX);
134: DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX);
135: DMGlobalToLocalBegin(dm, X_t, INSERT_VALUES, locX_t);
136: DMGlobalToLocalEnd(dm, X_t, INSERT_VALUES, locX_t);
137: CHKMEMQ;
138: (*dmlocalts->ijacobianlocal)(dm, time, locX, locX_t, a, A, B, dmlocalts->ijacobianlocalctx);
139: CHKMEMQ;
140: DMRestoreLocalVector(dm, &locX);
141: DMRestoreLocalVector(dm, &locX_t);
142: } else {
143: MatFDColoring fdcoloring;
144: PetscObjectQuery((PetscObject) dm, "DMDASNES_FDCOLORING", (PetscObject *) &fdcoloring);
145: if (!fdcoloring) {
146: ISColoring coloring;
148: DMCreateColoring(dm, dm->coloringtype, &coloring);
149: MatFDColoringCreate(B, coloring, &fdcoloring);
150: ISColoringDestroy(&coloring);
151: switch (dm->coloringtype) {
152: case IS_COLORING_GLOBAL:
153: MatFDColoringSetFunction(fdcoloring, (PetscErrorCode (*)(void)) TSComputeIFunction_DMLocal, dmlocalts);
154: break;
155: default: SETERRQ1(PetscObjectComm((PetscObject) ts), PETSC_ERR_SUP, "No support for coloring type '%s'", ISColoringTypes[dm->coloringtype]);
156: }
157: PetscObjectSetOptionsPrefix((PetscObject) fdcoloring, ((PetscObject) dm)->prefix);
158: MatFDColoringSetFromOptions(fdcoloring);
159: MatFDColoringSetUp(B, coloring, fdcoloring);
160: PetscObjectCompose((PetscObject) dm, "DMDASNES_FDCOLORING", (PetscObject) fdcoloring);
161: PetscObjectDereference((PetscObject) fdcoloring);
163: /* The following breaks an ugly reference counting loop that deserves a paragraph. MatFDColoringApply() will call
164: * VecDuplicate() with the state Vec and store inside the MatFDColoring. This Vec will duplicate the Vec, but the
165: * MatFDColoring is composed with the DM. We dereference the DM here so that the reference count will eventually
166: * drop to 0. Note the code in DMDestroy() that exits early for a negative reference count. That code path will be
167: * taken when the PetscObjectList for the Vec inside MatFDColoring is destroyed.
168: */
169: PetscObjectDereference((PetscObject) dm);
170: }
171: MatFDColoringApply(B, fdcoloring, X, ts);
172: }
173: /* This will be redundant if the user called both, but it's too common to forget. */
174: if (A != B) {
175: MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
176: MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
177: }
178: return(0);
179: }
183: /*@C
184: DMTSSetIFunctionLocal - set a local implicit function evaluation function. This function is called with local vector
185: containing the local vector information PLUS ghost point information. It should compute a result for all local
186: elements and DMTS will automatically accumulate the overlapping values.
188: Logically Collective
190: Input Arguments:
191: + dm - DM to associate callback with
192: . func - local function evaluation
193: - ctx - context for function evaluation
195: Level: beginner
197: .seealso: DMTSSetIFunction(), DMTSSetIJacobianLocal()
198: @*/
199: PetscErrorCode DMTSSetIFunctionLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, Vec, void *), void *ctx)
200: {
201: DMTS tdm;
202: DMTS_Local *dmlocalts;
207: DMGetDMTSWrite(dm, &tdm);
208: DMLocalTSGetContext(dm, tdm, &dmlocalts);
210: dmlocalts->ifunctionlocal = func;
211: dmlocalts->ifunctionlocalctx = ctx;
213: DMTSSetIFunction(dm, TSComputeIFunction_DMLocal, dmlocalts);
214: if (!tdm->ops->ijacobian) { /* Call us for the Jacobian too, can be overridden by the user. */
215: DMTSSetIJacobian(dm, TSComputeIJacobian_DMLocal, dmlocalts);
216: }
217: return(0);
218: }
222: /*@C
223: DMTSSetIJacobianLocal - set a local Jacobian evaluation function
225: Logically Collective
227: Input Arguments:
228: + dm - DM to associate callback with
229: . func - local Jacobian evaluation
230: - ctx - optional context for local Jacobian evaluation
232: Level: beginner
234: .seealso: DMTSSetIFunctionLocal(), DMTSSetIJacobian(), DMTSSetIFunction()
235: @*/
236: PetscErrorCode DMTSSetIJacobianLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *), void *ctx)
237: {
238: DMTS tdm;
239: DMTS_Local *dmlocalts;
244: DMGetDMTSWrite(dm, &tdm);
245: DMLocalTSGetContext(dm, tdm, &dmlocalts);
247: dmlocalts->ijacobianlocal = func;
248: dmlocalts->ijacobianlocalctx = ctx;
250: DMTSSetIJacobian(dm, TSComputeIJacobian_DMLocal, dmlocalts);
251: return(0);
252: }
256: /*@C
257: DMTSSetRHSFunctionLocal - set a local rhs function evaluation function. This function is called with local vector
258: containing the local vector information PLUS ghost point information. It should compute a result for all local
259: elements and DMTS will automatically accumulate the overlapping values.
261: Logically Collective
263: Input Arguments:
264: + dm - DM to associate callback with
265: . func - local function evaluation
266: - ctx - context for function evaluation
268: Level: beginner
270: .seealso: DMTSSetRHSFunction(), DMTSSetIFunction(), DMTSSetIJacobianLocal()
271: @*/
272: PetscErrorCode DMTSSetRHSFunctionLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, void *), void *ctx)
273: {
274: DMTS tdm;
275: DMTS_Local *dmlocalts;
280: DMGetDMTSWrite(dm, &tdm);
281: DMLocalTSGetContext(dm, tdm, &dmlocalts);
283: dmlocalts->rhsfunctionlocal = func;
284: dmlocalts->rhsfunctionlocalctx = ctx;
286: DMTSSetRHSFunction(dm, TSComputeRHSFunction_DMLocal, dmlocalts);
287: return(0);
288: }