Actual source code: dmlocalts.c

petsc-master 2019-12-10
Report Typos and Errors
  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: } DMTS_Local;

 15: static PetscErrorCode DMTSDestroy_DMLocal(DMTS tdm)
 16: {

 20:   PetscFree(tdm->data);
 21:   return(0);
 22: }

 24: static PetscErrorCode DMTSDuplicate_DMLocal(DMTS oldtdm, DMTS tdm)
 25: {

 29:   PetscNewLog(tdm, (DMTS_Local **) &tdm->data);
 30:   if (oldtdm->data) {PetscMemcpy(tdm->data, oldtdm->data, sizeof(DMTS_Local));}
 31:   return(0);
 32: }

 34: static PetscErrorCode DMLocalTSGetContext(DM dm, DMTS tdm, DMTS_Local **dmlocalts)
 35: {

 39:   *dmlocalts = NULL;
 40:   if (!tdm->data) {
 41:     PetscNewLog(dm, (DMTS_Local **) &tdm->data);

 43:     tdm->ops->destroy   = DMTSDestroy_DMLocal;
 44:     tdm->ops->duplicate = DMTSDuplicate_DMLocal;
 45:   }
 46:   *dmlocalts = (DMTS_Local *) tdm->data;
 47:   return(0);
 48: }

 50: static PetscErrorCode TSComputeIFunction_DMLocal(TS ts, PetscReal time, Vec X, Vec X_t, Vec F, void *ctx)
 51: {
 52:   DM             dm;
 53:   Vec            locX, locX_t, locF;
 54:   DMTS_Local    *dmlocalts = (DMTS_Local *) ctx;

 62:   TSGetDM(ts, &dm);
 63:   DMGetLocalVector(dm, &locX);
 64:   DMGetLocalVector(dm, &locX_t);
 65:   DMGetLocalVector(dm, &locF);
 66:   VecZeroEntries(locX);
 67:   VecZeroEntries(locX_t);
 68:   if (dmlocalts->boundarylocal) {(*dmlocalts->boundarylocal)(dm, time, locX, locX_t,dmlocalts->boundarylocalctx);}
 69:   DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX);
 70:   DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX);
 71:   DMGlobalToLocalBegin(dm, X_t, INSERT_VALUES, locX_t);
 72:   DMGlobalToLocalEnd(dm, X_t, INSERT_VALUES, locX_t);
 73:   VecZeroEntries(locF);
 74:   CHKMEMQ;
 75:   (*dmlocalts->ifunctionlocal)(dm, time, locX, locX_t, locF, dmlocalts->ifunctionlocalctx);
 76:   CHKMEMQ;
 77:   VecZeroEntries(F);
 78:   DMLocalToGlobalBegin(dm, locF, ADD_VALUES, F);
 79:   DMLocalToGlobalEnd(dm, locF, ADD_VALUES, F);
 80:   DMRestoreLocalVector(dm, &locX);
 81:   DMRestoreLocalVector(dm, &locX_t);
 82:   DMRestoreLocalVector(dm, &locF);
 83:   return(0);
 84: }

 86: static PetscErrorCode TSComputeRHSFunction_DMLocal(TS ts, PetscReal time, Vec X, Vec F, void *ctx)
 87: {
 88:   DM             dm;
 89:   Vec            locX;
 90:   DMTS_Local    *dmlocalts = (DMTS_Local *) ctx;

 97:   TSGetDM(ts, &dm);
 98:   DMGetLocalVector(dm, &locX);
 99:   VecZeroEntries(locX);
100:   if (dmlocalts->boundarylocal) {(*dmlocalts->boundarylocal)(dm,time,locX,NULL,dmlocalts->boundarylocalctx);}
101:   DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX);
102:   DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX);
103:   VecZeroEntries(F);
104:   CHKMEMQ;
105:   (*dmlocalts->rhsfunctionlocal)(dm, time, locX, F, dmlocalts->rhsfunctionlocalctx);
106:   CHKMEMQ;
107:   DMRestoreLocalVector(dm, &locX);
108:   return(0);
109: }

111: static PetscErrorCode TSComputeIJacobian_DMLocal(TS ts, PetscReal time, Vec X, Vec X_t, PetscReal a, Mat A, Mat B, void *ctx)
112: {
113:   DM             dm;
114:   Vec            locX, locX_t;
115:   DMTS_Local    *dmlocalts = (DMTS_Local *) ctx;

119:   TSGetDM(ts, &dm);
120:   if (dmlocalts->ijacobianlocal) {
121:     DMGetLocalVector(dm, &locX);
122:     DMGetLocalVector(dm, &locX_t);
123:     VecZeroEntries(locX);
124:     if (dmlocalts->boundarylocal) {(*dmlocalts->boundarylocal)(dm,time,locX,locX_t,dmlocalts->boundarylocalctx);}
125:     VecZeroEntries(locX_t);
126:     DMGlobalToLocalBegin(dm, X, INSERT_VALUES, locX);
127:     DMGlobalToLocalEnd(dm, X, INSERT_VALUES, locX);
128:     DMGlobalToLocalBegin(dm, X_t, INSERT_VALUES, locX_t);
129:     DMGlobalToLocalEnd(dm, X_t, INSERT_VALUES, locX_t);
130:     CHKMEMQ;
131:     (*dmlocalts->ijacobianlocal)(dm, time, locX, locX_t, a, A, B, dmlocalts->ijacobianlocalctx);
132:     CHKMEMQ;
133:     DMRestoreLocalVector(dm, &locX);
134:     DMRestoreLocalVector(dm, &locX_t);
135:   } else {
136:     MatFDColoring fdcoloring;
137:     PetscObjectQuery((PetscObject) dm, "DMDASNES_FDCOLORING", (PetscObject *) &fdcoloring);
138:     if (!fdcoloring) {
139:       ISColoring coloring;

141:       DMCreateColoring(dm, dm->coloringtype, &coloring);
142:       MatFDColoringCreate(B, coloring, &fdcoloring);
143:       ISColoringDestroy(&coloring);
144:       switch (dm->coloringtype) {
145:       case IS_COLORING_GLOBAL:
146:         MatFDColoringSetFunction(fdcoloring, (PetscErrorCode (*)(void)) TSComputeIFunction_DMLocal, dmlocalts);
147:         break;
148:       default: SETERRQ1(PetscObjectComm((PetscObject) ts), PETSC_ERR_SUP, "No support for coloring type '%s'", ISColoringTypes[dm->coloringtype]);
149:       }
150:       PetscObjectSetOptionsPrefix((PetscObject) fdcoloring, ((PetscObject) dm)->prefix);
151:       MatFDColoringSetFromOptions(fdcoloring);
152:       MatFDColoringSetUp(B, coloring, fdcoloring);
153:       PetscObjectCompose((PetscObject) dm, "DMDASNES_FDCOLORING", (PetscObject) fdcoloring);
154:       PetscObjectDereference((PetscObject) fdcoloring);

156:       /* The following breaks an ugly reference counting loop that deserves a paragraph. MatFDColoringApply() will call
157:        * VecDuplicate() with the state Vec and store inside the MatFDColoring. This Vec will duplicate the Vec, but the
158:        * MatFDColoring is composed with the DM. We dereference the DM here so that the reference count will eventually
159:        * drop to 0. Note the code in DMDestroy() that exits early for a negative reference count. That code path will be
160:        * taken when the PetscObjectList for the Vec inside MatFDColoring is destroyed.
161:        */
162:       PetscObjectDereference((PetscObject) dm);
163:     }
164:     MatFDColoringApply(B, fdcoloring, X, ts);
165:   }
166:   /* This will be redundant if the user called both, but it's too common to forget. */
167:   if (A != B) {
168:     MatAssemblyBegin(A, MAT_FINAL_ASSEMBLY);
169:     MatAssemblyEnd(A, MAT_FINAL_ASSEMBLY);
170:   }
171:   return(0);
172: }

174: /*@C
175:   DMTSSetBoundaryLocal - set the function for essential boundary data for a local implicit function evaluation.
176:     It should set the essential boundary data for the local portion of the solution X, as well its time derivative X_t (if it is not NULL).
177:     Vectors are initialized to zero before this function, so it is only needed for non homogeneous data.

179:   Note that this function is somewhat optional: boundary data could potentially be inserted by a function passed to
180:   DMTSSetIFunctionLocal().  The use case for this function is for discretizations with constraints (see
181:   DMGetDefaultConstraints()): this function inserts boundary values before constraint interpolation.

183:   Logically Collective

185:   Input Arguments:
186: + dm   - DM to associate callback with
187: . func - local function evaluation
188: - ctx  - context for function evaluation

190:   Level: intermediate

192: .seealso: DMTSSetIFunction(), DMTSSetIJacobianLocal()
193: @*/
194: PetscErrorCode DMTSSetBoundaryLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, void *), void *ctx)
195: {
196:   DMTS           tdm;
197:   DMTS_Local    *dmlocalts;

202:   DMGetDMTSWrite(dm, &tdm);
203:   DMLocalTSGetContext(dm, tdm, &dmlocalts);

205:   dmlocalts->boundarylocal    = func;
206:   dmlocalts->boundarylocalctx = ctx;

208:   return(0);
209: }

211: /*@C
212:   DMTSSetIFunctionLocal - set a local implicit function evaluation function. This function is called with local vector
213:       containing the local vector information PLUS ghost point information. It should compute a result for all local
214:       elements and DMTS will automatically accumulate the overlapping values.

216:   Logically Collective

218:   Input Arguments:
219: + dm   - DM to associate callback with
220: . func - local function evaluation
221: - ctx  - context for function evaluation

223:   Level: beginner

225: .seealso: DMTSSetIFunction(), DMTSSetIJacobianLocal()
226: @*/
227: PetscErrorCode DMTSSetIFunctionLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, Vec, void *), void *ctx)
228: {
229:   DMTS           tdm;
230:   DMTS_Local    *dmlocalts;

235:   DMGetDMTSWrite(dm, &tdm);
236:   DMLocalTSGetContext(dm, tdm, &dmlocalts);

238:   dmlocalts->ifunctionlocal    = func;
239:   dmlocalts->ifunctionlocalctx = ctx;

241:   DMTSSetIFunction(dm, TSComputeIFunction_DMLocal, dmlocalts);
242:   if (!tdm->ops->ijacobian) {  /* Call us for the Jacobian too, can be overridden by the user. */
243:     DMTSSetIJacobian(dm, TSComputeIJacobian_DMLocal, dmlocalts);
244:   }
245:   return(0);
246: }

248: /*@C
249:   DMTSSetIJacobianLocal - set a local Jacobian evaluation function

251:   Logically Collective

253:   Input Arguments:
254: + dm - DM to associate callback with
255: . func - local Jacobian evaluation
256: - ctx - optional context for local Jacobian evaluation

258:   Level: beginner

260: .seealso: DMTSSetIFunctionLocal(), DMTSSetIJacobian(), DMTSSetIFunction()
261: @*/
262: PetscErrorCode DMTSSetIJacobianLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, PetscReal, Mat, Mat, void *), void *ctx)
263: {
264:   DMTS           tdm;
265:   DMTS_Local    *dmlocalts;

270:   DMGetDMTSWrite(dm, &tdm);
271:   DMLocalTSGetContext(dm, tdm, &dmlocalts);

273:   dmlocalts->ijacobianlocal    = func;
274:   dmlocalts->ijacobianlocalctx = ctx;

276:   DMTSSetIJacobian(dm, TSComputeIJacobian_DMLocal, dmlocalts);
277:   return(0);
278: }

280: /*@C
281:   DMTSSetRHSFunctionLocal - set a local rhs function evaluation function. This function is called with local vector
282:       containing the local vector information PLUS ghost point information. It should compute a result for all local
283:       elements and DMTS will automatically accumulate the overlapping values.

285:   Logically Collective

287:   Input Arguments:
288: + dm   - DM to associate callback with
289: . func - local function evaluation
290: - ctx  - context for function evaluation

292:   Level: beginner

294: .seealso: DMTSSetRHSFunction(), DMTSSetIFunction(), DMTSSetIJacobianLocal()
295: @*/
296: PetscErrorCode DMTSSetRHSFunctionLocal(DM dm, PetscErrorCode (*func)(DM, PetscReal, Vec, Vec, void *), void *ctx)
297: {
298:   DMTS           tdm;
299:   DMTS_Local    *dmlocalts;

304:   DMGetDMTSWrite(dm, &tdm);
305:   DMLocalTSGetContext(dm, tdm, &dmlocalts);

307:   dmlocalts->rhsfunctionlocal    = func;
308:   dmlocalts->rhsfunctionlocalctx = ctx;

310:   DMTSSetRHSFunction(dm, TSComputeRHSFunction_DMLocal, dmlocalts);
311:   return(0);
312: }