Actual source code: dmlocalts.c
petsc-3.5.4 2015-05-23
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: #if 0
98: Vec locF;
99: #endif
100: DMTS_Local *dmlocalts = (DMTS_Local *) ctx;
107: TSGetDM(ts, &dm);
108: DMGetLocalVector(dm, &locX);
109: VecZeroEntries(locX);
110: DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX);
111: DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX);
112: VecZeroEntries(F);
113: CHKMEMQ;
114: (*dmlocalts->rhsfunctionlocal)(dm, time, locX, F, dmlocalts->rhsfunctionlocalctx);
115: CHKMEMQ;
116: DMRestoreLocalVector(dm, &locX);
117: return(0);
118: }
122: static PetscErrorCode TSComputeIJacobian_DMLocal(TS ts, PetscReal time, Vec X, Vec X_t, PetscReal a, Mat A, Mat B, void *ctx)
123: {
124: DM dm;
125: Vec locX, locX_t;
126: DMTS_Local *dmlocalts = (DMTS_Local *) ctx;
130: TSGetDM(ts, &dm);
131: if (dmlocalts->ijacobianlocal) {
132: DMGetLocalVector(dm, &locX);
133: DMGetLocalVector(dm, &locX_t);
134: VecZeroEntries(locX);
135: VecZeroEntries(locX_t);
136: DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX);
137: DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX);
138: DMGlobalToLocalBegin(dm, X_t, INSERT_VALUES, locX_t);
139: DMGlobalToLocalEnd(dm, X_t, INSERT_VALUES, locX_t);
140: CHKMEMQ;
141: (*dmlocalts->ijacobianlocal)(dm, time, locX, locX_t, a, A, B, dmlocalts->ijacobianlocalctx);
142: CHKMEMQ;
143: DMRestoreLocalVector(dm, &locX);
144: DMRestoreLocalVector(dm, &locX_t);
145: } else {
146: MatFDColoring fdcoloring;
147: PetscObjectQuery((PetscObject) dm, "DMDASNES_FDCOLORING", (PetscObject *) &fdcoloring);
148: if (!fdcoloring) {
149: ISColoring coloring;
151: DMCreateColoring(dm, dm->coloringtype, &coloring);
152: MatFDColoringCreate(B, coloring, &fdcoloring);
153: ISColoringDestroy(&coloring);
154: switch (dm->coloringtype) {
155: case IS_COLORING_GLOBAL:
156: MatFDColoringSetFunction(fdcoloring, (PetscErrorCode (*)(void)) TSComputeIFunction_DMLocal, dmlocalts);
157: break;
158: default: SETERRQ1(PetscObjectComm((PetscObject) ts), PETSC_ERR_SUP, "No support for coloring type '%s'", ISColoringTypes[dm->coloringtype]);
159: }
160: PetscObjectSetOptionsPrefix((PetscObject) fdcoloring, ((PetscObject) dm)->prefix);
161: MatFDColoringSetFromOptions(fdcoloring);
162: MatFDColoringSetUp(B, coloring, fdcoloring);
163: PetscObjectCompose((PetscObject) dm, "DMDASNES_FDCOLORING", (PetscObject) fdcoloring);
164: PetscObjectDereference((PetscObject) fdcoloring);
166: /* The following breaks an ugly reference counting loop that deserves a paragraph. MatFDColoringApply() will call
167: * VecDuplicate() with the state Vec and store inside the MatFDColoring. This Vec will duplicate the Vec, but the
168: * MatFDColoring is composed with the DM. We dereference the DM here so that the reference count will eventually
169: * drop to 0. Note the code in DMDestroy() that exits early for a negative reference count. That code path will be
170: * taken when the PetscObjectList for the Vec inside MatFDColoring is destroyed.
171: */
172: PetscObjectDereference((PetscObject) dm);
173: }
174: MatFDColoringApply(B, fdcoloring, X, ts);
175: }
176: /* This will be redundant if the user called both, but it's too common to forget. */
177: if (A != B) {
178: MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
179: MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
180: }
181: return(0);
182: }
186: /*@C
187: DMTSSetIFunctionLocal - set a local implicit function evaluation function. This function is called with local vector
188: containing the local vector information PLUS ghost point information. It should compute a result for all local
189: elements and DMTS will automatically accumulate the overlapping values.
191: Logically Collective
193: Input Arguments:
194: + dm - DM to associate callback with
195: . func - local function evaluation
196: - ctx - context for function evaluation
198: Level: beginner
200: .seealso: DMTSSetIFunction(), DMTSSetIJacobianLocal()
201: @*/
202: PetscErrorCode DMTSSetIFunctionLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, Vec, void *), void *ctx)
203: {
204: DMTS tdm;
205: DMTS_Local *dmlocalts;
210: DMGetDMTSWrite(dm, &tdm);
211: DMLocalTSGetContext(dm, tdm, &dmlocalts);
213: dmlocalts->ifunctionlocal = func;
214: dmlocalts->ifunctionlocalctx = ctx;
216: DMTSSetIFunction(dm, TSComputeIFunction_DMLocal, dmlocalts);
217: if (!tdm->ops->ijacobian) { /* Call us for the Jacobian too, can be overridden by the user. */
218: DMTSSetIJacobian(dm, TSComputeIJacobian_DMLocal, dmlocalts);
219: }
220: return(0);
221: }
225: /*@C
226: DMTSSetIJacobianLocal - set a local Jacobian evaluation function
228: Logically Collective
230: Input Arguments:
231: + dm - DM to associate callback with
232: . func - local Jacobian evaluation
233: - ctx - optional context for local Jacobian evaluation
235: Level: beginner
237: .seealso: DMTSSetIFunctionLocal(), DMTSSetIJacobian(), DMTSSetIFunction()
238: @*/
239: PetscErrorCode DMTSSetIJacobianLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *), void *ctx)
240: {
241: DMTS tdm;
242: DMTS_Local *dmlocalts;
247: DMGetDMTSWrite(dm, &tdm);
248: DMLocalTSGetContext(dm, tdm, &dmlocalts);
250: dmlocalts->ijacobianlocal = func;
251: dmlocalts->ijacobianlocalctx = ctx;
253: DMTSSetIJacobian(dm, TSComputeIJacobian_DMLocal, dmlocalts);
254: return(0);
255: }
259: /*@C
260: DMTSSetRHSFunctionLocal - set a local rhs function evaluation function. This function is called with local vector
261: containing the local vector information PLUS ghost point information. It should compute a result for all local
262: elements and DMTS will automatically accumulate the overlapping values.
264: Logically Collective
266: Input Arguments:
267: + dm - DM to associate callback with
268: . func - local function evaluation
269: - ctx - context for function evaluation
271: Level: beginner
273: .seealso: DMTSSetRHSFunction(), DMTSSetIFunction(), DMTSSetIJacobianLocal()
274: @*/
275: PetscErrorCode DMTSSetRHSFunctionLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, void *), void *ctx)
276: {
277: DMTS tdm;
278: DMTS_Local *dmlocalts;
283: DMGetDMTSWrite(dm, &tdm);
284: DMLocalTSGetContext(dm, tdm, &dmlocalts);
286: dmlocalts->rhsfunctionlocal = func;
287: dmlocalts->rhsfunctionlocalctx = ctx;
289: DMTSSetRHSFunction(dm, TSComputeRHSFunction_DMLocal, dmlocalts);
290: return(0);
291: }