Actual source code: dtds.c

  1: #include <petsc/private/petscdsimpl.h>

  3: PetscClassId PETSCDS_CLASSID = 0;

  5: PetscFunctionList PetscDSList              = NULL;
  6: PetscBool         PetscDSRegisterAllCalled = PETSC_FALSE;

  8: /* A PetscDS (Discrete System) encodes a set of equations posed in a discrete space, which represents a set of
  9:    nonlinear continuum equations. The equations can have multiple fields, each field having a different
 10:    discretization. In addition, different pieces of the domain can have different field combinations and equations.

 12:    The DS provides the user a description of the approximation space on any given cell. It also gives pointwise
 13:    functions representing the equations.

 15:    Each field is associated with a label, marking the cells on which it is supported. Note that a field can be
 16:    supported on the closure of a cell not in the label due to overlap of the boundary of neighboring cells. The DM
 17:    then creates a DS for each set of cells with identical approximation spaces. When assembling, the user asks for
 18:    the space associated with a given cell. DMPlex uses the labels associated with each DS in the default integration loop.
 19: */

 21: /*@C
 22:   PetscDSRegister - Adds a new `PetscDS` implementation

 24:   Not Collective

 26:   Input Parameters:
 27: + name        - The name of a new user-defined creation routine
 28: - create_func - The creation routine itself

 30:   Sample usage:
 31: .vb
 32:     PetscDSRegister("my_ds", MyPetscDSCreate);
 33: .ve

 35:   Then, your PetscDS type can be chosen with the procedural interface via
 36: .vb
 37:     PetscDSCreate(MPI_Comm, PetscDS *);
 38:     PetscDSSetType(PetscDS, "my_ds");
 39: .ve
 40:    or at runtime via the option
 41: .vb
 42:     -petscds_type my_ds
 43: .ve

 45:   Level: advanced

 47:   Note:
 48:   `PetscDSRegister()` may be called multiple times to add several user-defined `PetscDSs`

 50:   Fortran Note:
 51:   Not available from Fortran

 53: .seealso: `PetscDSType`, `PetscDS`, `PetscDSRegisterAll()`, `PetscDSRegisterDestroy()`
 54: @*/
 55: PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS))
 56: {
 57:   PetscFunctionListAdd(&PetscDSList, sname, function);
 58:   return 0;
 59: }

 61: /*@C
 62:   PetscDSSetType - Builds a particular `PetscDS`

 64:   Collective on prob

 66:   Input Parameters:
 67: + prob - The `PetscDS` object
 68: - name - The `PetscDSType`

 70:   Options Database Key:
 71: . -petscds_type <type> - Sets the PetscDS type; use -help for a list of available types

 73:   Level: intermediate

 75:   Fortran Note:
 76:   Not available from Fortran

 78: .seealso: `PetscDSType`, `PetscDS`, `PetscDSGetType()`, `PetscDSCreate()`
 79: @*/
 80: PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name)
 81: {
 82:   PetscErrorCode (*r)(PetscDS);
 83:   PetscBool match;

 86:   PetscObjectTypeCompare((PetscObject)prob, name, &match);
 87:   if (match) return 0;

 89:   PetscDSRegisterAll();
 90:   PetscFunctionListFind(PetscDSList, name, &r);

 93:   PetscTryTypeMethod(prob, destroy);
 94:   prob->ops->destroy = NULL;

 96:   (*r)(prob);
 97:   PetscObjectChangeTypeName((PetscObject)prob, name);
 98:   return 0;
 99: }

101: /*@C
102:   PetscDSGetType - Gets the `PetscDSType` name (as a string) from the `PetscDS`

104:   Not Collective

106:   Input Parameter:
107: . prob  - The `PetscDS`

109:   Output Parameter:
110: . name - The `PetscDSType` name

112:   Level: intermediate

114:   Fortran Note:
115:   Not available from Fortran

117: .seealso: `PetscDSType`, `PetscDS`, `PetscDSSetType()`, `PetscDSCreate()`
118: @*/
119: PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name)
120: {
123:   PetscDSRegisterAll();
124:   *name = ((PetscObject)prob)->type_name;
125:   return 0;
126: }

128: static PetscErrorCode PetscDSView_Ascii(PetscDS ds, PetscViewer viewer)
129: {
130:   PetscViewerFormat  format;
131:   const PetscScalar *constants;
132:   PetscInt           Nf, numConstants, f;

134:   PetscDSGetNumFields(ds, &Nf);
135:   PetscViewerGetFormat(viewer, &format);
136:   PetscViewerASCIIPrintf(viewer, "Discrete System with %" PetscInt_FMT " fields\n", Nf);
137:   PetscViewerASCIIPushTab(viewer);
138:   PetscViewerASCIIPrintf(viewer, "  cell total dim %" PetscInt_FMT " total comp %" PetscInt_FMT "\n", ds->totDim, ds->totComp);
139:   if (ds->isCohesive) PetscViewerASCIIPrintf(viewer, "  cohesive cell\n");
140:   for (f = 0; f < Nf; ++f) {
141:     DSBoundary      b;
142:     PetscObject     obj;
143:     PetscClassId    id;
144:     PetscQuadrature q;
145:     const char     *name;
146:     PetscInt        Nc, Nq, Nqc;

148:     PetscDSGetDiscretization(ds, f, &obj);
149:     PetscObjectGetClassId(obj, &id);
150:     PetscObjectGetName(obj, &name);
151:     PetscViewerASCIIPrintf(viewer, "Field %s", name ? name : "<unknown>");
152:     PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
153:     if (id == PETSCFE_CLASSID) {
154:       PetscFEGetNumComponents((PetscFE)obj, &Nc);
155:       PetscFEGetQuadrature((PetscFE)obj, &q);
156:       PetscViewerASCIIPrintf(viewer, " FEM");
157:     } else if (id == PETSCFV_CLASSID) {
158:       PetscFVGetNumComponents((PetscFV)obj, &Nc);
159:       PetscFVGetQuadrature((PetscFV)obj, &q);
160:       PetscViewerASCIIPrintf(viewer, " FVM");
161:     } else SETERRQ(PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f);
162:     if (Nc > 1) PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " components", Nc);
163:     else PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " component ", Nc);
164:     if (ds->implicit[f]) PetscViewerASCIIPrintf(viewer, " (implicit)");
165:     else PetscViewerASCIIPrintf(viewer, " (explicit)");
166:     if (q) {
167:       PetscQuadratureGetData(q, NULL, &Nqc, &Nq, NULL, NULL);
168:       PetscViewerASCIIPrintf(viewer, " (Nq %" PetscInt_FMT " Nqc %" PetscInt_FMT ")", Nq, Nqc);
169:     }
170:     PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT "-jet", ds->jetDegree[f]);
171:     PetscViewerASCIIPrintf(viewer, "\n");
172:     PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
173:     PetscViewerASCIIPushTab(viewer);
174:     if (id == PETSCFE_CLASSID) PetscFEView((PetscFE)obj, viewer);
175:     else if (id == PETSCFV_CLASSID) PetscFVView((PetscFV)obj, viewer);
176:     PetscViewerASCIIPopTab(viewer);

178:     for (b = ds->boundary; b; b = b->next) {
179:       char    *name;
180:       PetscInt c, i;

182:       if (b->field != f) continue;
183:       PetscViewerASCIIPushTab(viewer);
184:       PetscViewerASCIIPrintf(viewer, "Boundary %s (%s) %s\n", b->name, b->lname, DMBoundaryConditionTypes[b->type]);
185:       if (!b->Nc) {
186:         PetscViewerASCIIPrintf(viewer, "  all components\n");
187:       } else {
188:         PetscViewerASCIIPrintf(viewer, "  components: ");
189:         PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
190:         for (c = 0; c < b->Nc; ++c) {
191:           if (c > 0) PetscViewerASCIIPrintf(viewer, ", ");
192:           PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->comps[c]);
193:         }
194:         PetscViewerASCIIPrintf(viewer, "\n");
195:         PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
196:       }
197:       PetscViewerASCIIPrintf(viewer, "  values: ");
198:       PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
199:       for (i = 0; i < b->Nv; ++i) {
200:         if (i > 0) PetscViewerASCIIPrintf(viewer, ", ");
201:         PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->values[i]);
202:       }
203:       PetscViewerASCIIPrintf(viewer, "\n");
204:       PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
205:       if (b->func) {
206:         PetscDLAddr(b->func, &name);
207:         if (name) PetscViewerASCIIPrintf(viewer, "  func: %s\n", name);
208:         else PetscViewerASCIIPrintf(viewer, "  func: %p\n", b->func);
209:         PetscFree(name);
210:       }
211:       if (b->func_t) {
212:         PetscDLAddr(b->func_t, &name);
213:         if (name) PetscViewerASCIIPrintf(viewer, "  func_t: %s\n", name);
214:         else PetscViewerASCIIPrintf(viewer, "  func_t: %p\n", b->func_t);
215:         PetscFree(name);
216:       }
217:       PetscWeakFormView(b->wf, viewer);
218:       PetscViewerASCIIPopTab(viewer);
219:     }
220:   }
221:   PetscDSGetConstants(ds, &numConstants, &constants);
222:   if (numConstants) {
223:     PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " constants\n", numConstants);
224:     PetscViewerASCIIPushTab(viewer);
225:     for (f = 0; f < numConstants; ++f) PetscViewerASCIIPrintf(viewer, "%g\n", (double)PetscRealPart(constants[f]));
226:     PetscViewerASCIIPopTab(viewer);
227:   }
228:   PetscWeakFormView(ds->wf, viewer);
229:   PetscViewerASCIIPopTab(viewer);
230:   return 0;
231: }

233: /*@C
234:    PetscDSViewFromOptions - View a `PetscDS` based on values in the options database

236:    Collective on PetscDS

238:    Input Parameters:
239: +  A - the `PetscDS` object
240: .  obj - Optional object
241: -  name - command line option

243:    Level: intermediate

245: .seealso: `PetscDSType`, `PetscDS`, `PetscDSView()`, `PetscObjectViewFromOptions()`, `PetscDSCreate()`
246: @*/
247: PetscErrorCode PetscDSViewFromOptions(PetscDS A, PetscObject obj, const char name[])
248: {
250:   PetscObjectViewFromOptions((PetscObject)A, obj, name);
251:   return 0;
252: }

254: /*@C
255:   PetscDSView - Views a `PetscDS`

257:   Collective on prob

259:   Input Parameters:
260: + prob - the `PetscDS` object to view
261: - v  - the viewer

263:   Level: developer

265: .seealso: `PetscDSType`, `PetscDS`, `PetscViewer`, `PetscDSDestroy()`
266: @*/
267: PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v)
268: {
269:   PetscBool iascii;

272:   if (!v) PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)prob), &v);
274:   PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &iascii);
275:   if (iascii) PetscDSView_Ascii(prob, v);
276:   PetscTryTypeMethod(prob, view, v);
277:   return 0;
278: }

280: /*@
281:   PetscDSSetFromOptions - sets parameters in a `PetscDS` from the options database

283:   Collective on prob

285:   Input Parameter:
286: . prob - the `PetscDS` object to set options for

288:   Options Database Keys:
289: + -petscds_type <type>     - Set the `PetscDS` type
290: . -petscds_view <view opt> - View the `PetscDS`
291: . -petscds_jac_pre         - Turn formation of a separate Jacobian preconditioner on or off
292: . -bc_<name> <ids>         - Specify a list of label ids for a boundary condition
293: - -bc_<name>_comp <comps>  - Specify a list of field components to constrain for a boundary condition

295:   Level: intermediate

297: .seealso: `PetscDS`, `PetscDSView()`
298: @*/
299: PetscErrorCode PetscDSSetFromOptions(PetscDS prob)
300: {
301:   DSBoundary  b;
302:   const char *defaultType;
303:   char        name[256];
304:   PetscBool   flg;

307:   if (!((PetscObject)prob)->type_name) {
308:     defaultType = PETSCDSBASIC;
309:   } else {
310:     defaultType = ((PetscObject)prob)->type_name;
311:   }
312:   PetscDSRegisterAll();

314:   PetscObjectOptionsBegin((PetscObject)prob);
315:   for (b = prob->boundary; b; b = b->next) {
316:     char      optname[1024];
317:     PetscInt  ids[1024], len = 1024;
318:     PetscBool flg;

320:     PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name);
321:     PetscMemzero(ids, sizeof(ids));
322:     PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg);
323:     if (flg) {
324:       b->Nv = len;
325:       PetscFree(b->values);
326:       PetscMalloc1(len, &b->values);
327:       PetscArraycpy(b->values, ids, len);
328:       PetscWeakFormRewriteKeys(b->wf, b->label, len, b->values);
329:     }
330:     len = 1024;
331:     PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name);
332:     PetscMemzero(ids, sizeof(ids));
333:     PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg);
334:     if (flg) {
335:       b->Nc = len;
336:       PetscFree(b->comps);
337:       PetscMalloc1(len, &b->comps);
338:       PetscArraycpy(b->comps, ids, len);
339:     }
340:   }
341:   PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg);
342:   if (flg) {
343:     PetscDSSetType(prob, name);
344:   } else if (!((PetscObject)prob)->type_name) {
345:     PetscDSSetType(prob, defaultType);
346:   }
347:   PetscOptionsBool("-petscds_jac_pre", "Discrete System", "PetscDSUseJacobianPreconditioner", prob->useJacPre, &prob->useJacPre, &flg);
348:   PetscTryTypeMethod(prob, setfromoptions);
349:   /* process any options handlers added with PetscObjectAddOptionsHandler() */
350:   PetscObjectProcessOptionsHandlers((PetscObject)prob, PetscOptionsObject);
351:   PetscOptionsEnd();
352:   if (prob->Nf) PetscDSViewFromOptions(prob, NULL, "-petscds_view");
353:   return 0;
354: }

356: /*@C
357:   PetscDSSetUp - Construct data structures for the `PetscDS`

359:   Collective on prob

361:   Input Parameter:
362: . prob - the `PetscDS` object to setup

364:   Level: developer

366: .seealso: `PetscDS`, `PetscDSView()`, `PetscDSDestroy()`
367: @*/
368: PetscErrorCode PetscDSSetUp(PetscDS prob)
369: {
370:   const PetscInt Nf   = prob->Nf;
371:   PetscBool      hasH = PETSC_FALSE;
372:   PetscInt       dim, dimEmbed, NbMax = 0, NcMax = 0, NqMax = 0, NsMax = 1, f;

375:   if (prob->setup) return 0;
376:   /* Calculate sizes */
377:   PetscDSGetSpatialDimension(prob, &dim);
378:   PetscDSGetCoordinateDimension(prob, &dimEmbed);
379:   prob->totDim = prob->totComp = 0;
380:   PetscMalloc2(Nf, &prob->Nc, Nf, &prob->Nb);
381:   PetscCalloc2(Nf + 1, &prob->off, Nf + 1, &prob->offDer);
382:   PetscCalloc6(Nf + 1, &prob->offCohesive[0], Nf + 1, &prob->offCohesive[1], Nf + 1, &prob->offCohesive[2], Nf + 1, &prob->offDerCohesive[0], Nf + 1, &prob->offDerCohesive[1], Nf + 1, &prob->offDerCohesive[2]);
383:   PetscMalloc2(Nf, &prob->T, Nf, &prob->Tf);
384:   for (f = 0; f < Nf; ++f) {
385:     PetscObject     obj;
386:     PetscClassId    id;
387:     PetscQuadrature q  = NULL;
388:     PetscInt        Nq = 0, Nb, Nc;

390:     PetscDSGetDiscretization(prob, f, &obj);
391:     if (prob->jetDegree[f] > 1) hasH = PETSC_TRUE;
392:     if (!obj) {
393:       /* Empty mesh */
394:       Nb = Nc    = 0;
395:       prob->T[f] = prob->Tf[f] = NULL;
396:     } else {
397:       PetscObjectGetClassId(obj, &id);
398:       if (id == PETSCFE_CLASSID) {
399:         PetscFE fe = (PetscFE)obj;

401:         PetscFEGetQuadrature(fe, &q);
402:         PetscFEGetDimension(fe, &Nb);
403:         PetscFEGetNumComponents(fe, &Nc);
404:         PetscFEGetCellTabulation(fe, prob->jetDegree[f], &prob->T[f]);
405:         PetscFEGetFaceTabulation(fe, prob->jetDegree[f], &prob->Tf[f]);
406:       } else if (id == PETSCFV_CLASSID) {
407:         PetscFV fv = (PetscFV)obj;

409:         PetscFVGetQuadrature(fv, &q);
410:         PetscFVGetNumComponents(fv, &Nc);
411:         Nb = Nc;
412:         PetscFVGetCellTabulation(fv, &prob->T[f]);
413:         /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */
414:       } else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f);
415:     }
416:     prob->Nc[f]                    = Nc;
417:     prob->Nb[f]                    = Nb;
418:     prob->off[f + 1]               = Nc + prob->off[f];
419:     prob->offDer[f + 1]            = Nc * dim + prob->offDer[f];
420:     prob->offCohesive[0][f + 1]    = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[0][f];
421:     prob->offDerCohesive[0][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[0][f];
422:     prob->offCohesive[1][f]        = (prob->cohesive[f] ? 0 : Nc) + prob->offCohesive[0][f];
423:     prob->offDerCohesive[1][f]     = (prob->cohesive[f] ? 0 : Nc) * dimEmbed + prob->offDerCohesive[0][f];
424:     prob->offCohesive[2][f + 1]    = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[2][f];
425:     prob->offDerCohesive[2][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[2][f];
426:     if (q) PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL);
427:     NqMax = PetscMax(NqMax, Nq);
428:     NbMax = PetscMax(NbMax, Nb);
429:     NcMax = PetscMax(NcMax, Nc);
430:     prob->totDim += Nb;
431:     prob->totComp += Nc;
432:     /* There are two faces for all fields on a cohesive cell, except for cohesive fields */
433:     if (prob->isCohesive && !prob->cohesive[f]) prob->totDim += Nb;
434:   }
435:   prob->offCohesive[1][Nf]    = prob->offCohesive[0][Nf];
436:   prob->offDerCohesive[1][Nf] = prob->offDerCohesive[0][Nf];
437:   /* Allocate works space */
438:   NsMax = 2; /* A non-cohesive discretizations can be used on a cohesive cell, so we need this extra workspace for all DS */
439:   PetscMalloc3(NsMax * prob->totComp, &prob->u, NsMax * prob->totComp, &prob->u_t, NsMax * prob->totComp * dimEmbed + (hasH ? NsMax * prob->totComp * dimEmbed * dimEmbed : 0), &prob->u_x);
440:   PetscMalloc5(dimEmbed, &prob->x, NbMax * NcMax, &prob->basisReal, NbMax * NcMax * dimEmbed, &prob->basisDerReal, NbMax * NcMax, &prob->testReal, NbMax * NcMax * dimEmbed, &prob->testDerReal);
441:   PetscCall(PetscMalloc6(NsMax * NqMax * NcMax, &prob->f0, NsMax * NqMax * NcMax * dimEmbed, &prob->f1, NsMax * NsMax * NqMax * NcMax * NcMax, &prob->g0, NsMax * NsMax * NqMax * NcMax * NcMax * dimEmbed, &prob->g1, NsMax * NsMax * NqMax * NcMax * NcMax * dimEmbed,
442:                          &prob->g2, NsMax * NsMax * NqMax * NcMax * NcMax * dimEmbed * dimEmbed, &prob->g3));
443:   PetscTryTypeMethod(prob, setup);
444:   prob->setup = PETSC_TRUE;
445:   return 0;
446: }

448: static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob)
449: {
450:   PetscFree2(prob->Nc, prob->Nb);
451:   PetscFree2(prob->off, prob->offDer);
452:   PetscFree6(prob->offCohesive[0], prob->offCohesive[1], prob->offCohesive[2], prob->offDerCohesive[0], prob->offDerCohesive[1], prob->offDerCohesive[2]);
453:   PetscFree2(prob->T, prob->Tf);
454:   PetscFree3(prob->u, prob->u_t, prob->u_x);
455:   PetscFree5(prob->x, prob->basisReal, prob->basisDerReal, prob->testReal, prob->testDerReal);
456:   PetscFree6(prob->f0, prob->f1, prob->g0, prob->g1, prob->g2, prob->g3);
457:   return 0;
458: }

460: static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew)
461: {
462:   PetscObject          *tmpd;
463:   PetscBool            *tmpi;
464:   PetscInt             *tmpk;
465:   PetscBool            *tmpc;
466:   PetscPointFunc       *tmpup;
467:   PetscSimplePointFunc *tmpexactSol, *tmpexactSol_t;
468:   void                **tmpexactCtx, **tmpexactCtx_t;
469:   void                **tmpctx;
470:   PetscInt              Nf = prob->Nf, f;

472:   if (Nf >= NfNew) return 0;
473:   prob->setup = PETSC_FALSE;
474:   PetscDSDestroyStructs_Static(prob);
475:   PetscMalloc4(NfNew, &tmpd, NfNew, &tmpi, NfNew, &tmpc, NfNew, &tmpk);
476:   for (f = 0; f < Nf; ++f) {
477:     tmpd[f] = prob->disc[f];
478:     tmpi[f] = prob->implicit[f];
479:     tmpc[f] = prob->cohesive[f];
480:     tmpk[f] = prob->jetDegree[f];
481:   }
482:   for (f = Nf; f < NfNew; ++f) {
483:     tmpd[f] = NULL;
484:     tmpi[f] = PETSC_TRUE, tmpc[f] = PETSC_FALSE;
485:     tmpk[f] = 1;
486:   }
487:   PetscFree4(prob->disc, prob->implicit, prob->cohesive, prob->jetDegree);
488:   PetscWeakFormSetNumFields(prob->wf, NfNew);
489:   prob->Nf        = NfNew;
490:   prob->disc      = tmpd;
491:   prob->implicit  = tmpi;
492:   prob->cohesive  = tmpc;
493:   prob->jetDegree = tmpk;
494:   PetscCalloc2(NfNew, &tmpup, NfNew, &tmpctx);
495:   for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f];
496:   for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f];
497:   for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL;
498:   for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL;
499:   PetscFree2(prob->update, prob->ctx);
500:   prob->update = tmpup;
501:   prob->ctx    = tmpctx;
502:   PetscCalloc4(NfNew, &tmpexactSol, NfNew, &tmpexactCtx, NfNew, &tmpexactSol_t, NfNew, &tmpexactCtx_t);
503:   for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f];
504:   for (f = 0; f < Nf; ++f) tmpexactCtx[f] = prob->exactCtx[f];
505:   for (f = 0; f < Nf; ++f) tmpexactSol_t[f] = prob->exactSol_t[f];
506:   for (f = 0; f < Nf; ++f) tmpexactCtx_t[f] = prob->exactCtx_t[f];
507:   for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL;
508:   for (f = Nf; f < NfNew; ++f) tmpexactCtx[f] = NULL;
509:   for (f = Nf; f < NfNew; ++f) tmpexactSol_t[f] = NULL;
510:   for (f = Nf; f < NfNew; ++f) tmpexactCtx_t[f] = NULL;
511:   PetscFree4(prob->exactSol, prob->exactCtx, prob->exactSol_t, prob->exactCtx_t);
512:   prob->exactSol   = tmpexactSol;
513:   prob->exactCtx   = tmpexactCtx;
514:   prob->exactSol_t = tmpexactSol_t;
515:   prob->exactCtx_t = tmpexactCtx_t;
516:   return 0;
517: }

519: /*@
520:   PetscDSDestroy - Destroys a PetscDS object

522:   Collective on prob

524:   Input Parameter:
525: . prob - the PetscDS object to destroy

527:   Level: developer

529: .seealso: `PetscDSView()`
530: @*/
531: PetscErrorCode PetscDSDestroy(PetscDS *ds)
532: {
533:   PetscInt f;

535:   if (!*ds) return 0;

538:   if (--((PetscObject)(*ds))->refct > 0) {
539:     *ds = NULL;
540:     return 0;
541:   }
542:   ((PetscObject)(*ds))->refct = 0;
543:   if ((*ds)->subprobs) {
544:     PetscInt dim, d;

546:     PetscDSGetSpatialDimension(*ds, &dim);
547:     for (d = 0; d < dim; ++d) PetscDSDestroy(&(*ds)->subprobs[d]);
548:   }
549:   PetscFree((*ds)->subprobs);
550:   PetscDSDestroyStructs_Static(*ds);
551:   for (f = 0; f < (*ds)->Nf; ++f) PetscObjectDereference((*ds)->disc[f]);
552:   PetscFree4((*ds)->disc, (*ds)->implicit, (*ds)->cohesive, (*ds)->jetDegree);
553:   PetscWeakFormDestroy(&(*ds)->wf);
554:   PetscFree2((*ds)->update, (*ds)->ctx);
555:   PetscFree4((*ds)->exactSol, (*ds)->exactCtx, (*ds)->exactSol_t, (*ds)->exactCtx_t);
556:   PetscTryTypeMethod((*ds), destroy);
557:   PetscDSDestroyBoundary(*ds);
558:   PetscFree((*ds)->constants);
559:   PetscHeaderDestroy(ds);
560:   return 0;
561: }

563: /*@
564:   PetscDSCreate - Creates an empty `PetscDS` object. The type can then be set with `PetscDSSetType()`.

566:   Collective

568:   Input Parameter:
569: . comm - The communicator for the `PetscDS` object

571:   Output Parameter:
572: . ds   - The `PetscDS` object

574:   Level: beginner

576: .seealso: `PetscDS`, `PetscDSSetType()`, `PETSCDSBASIC`, `PetscDSType`
577: @*/
578: PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *ds)
579: {
580:   PetscDS p;

583:   *ds = NULL;
584:   PetscDSInitializePackage();

586:   PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView);

588:   p->Nf           = 0;
589:   p->setup        = PETSC_FALSE;
590:   p->numConstants = 0;
591:   p->constants    = NULL;
592:   p->dimEmbed     = -1;
593:   p->useJacPre    = PETSC_TRUE;
594:   PetscWeakFormCreate(comm, &p->wf);

596:   *ds = p;
597:   return 0;
598: }

600: /*@
601:   PetscDSGetNumFields - Returns the number of fields in the `PetscDS`

603:   Not collective

605:   Input Parameter:
606: . prob - The PetscDS object

608:   Output Parameter:
609: . Nf - The number of fields

611:   Level: beginner

613: .seealso: `PetscDS`, `PetscDSGetSpatialDimension()`, `PetscDSCreate()`
614: @*/
615: PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf)
616: {
619:   *Nf = prob->Nf;
620:   return 0;
621: }

623: /*@
624:   PetscDSGetSpatialDimension - Returns the spatial dimension of the `PetscDS`, meaning the topological dimension of the discretizations

626:   Not collective

628:   Input Parameter:
629: . prob - The `PetscDS` object

631:   Output Parameter:
632: . dim - The spatial dimension

634:   Level: beginner

636: .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
637: @*/
638: PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim)
639: {
642:   *dim = 0;
643:   if (prob->Nf) {
644:     PetscObject  obj;
645:     PetscClassId id;

647:     PetscDSGetDiscretization(prob, 0, &obj);
648:     if (obj) {
649:       PetscObjectGetClassId(obj, &id);
650:       if (id == PETSCFE_CLASSID) PetscFEGetSpatialDimension((PetscFE)obj, dim);
651:       else if (id == PETSCFV_CLASSID) PetscFVGetSpatialDimension((PetscFV)obj, dim);
652:       else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
653:     }
654:   }
655:   return 0;
656: }

658: /*@
659:   PetscDSGetCoordinateDimension - Returns the coordinate dimension of the `PetscDS`, meaning the dimension of the space into which the discretiaztions are embedded

661:   Not collective

663:   Input Parameter:
664: . prob - The `PetscDS` object

666:   Output Parameter:
667: . dimEmbed - The coordinate dimension

669:   Level: beginner

671: .seealso: `PetscDS`, `PetscDSSetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
672: @*/
673: PetscErrorCode PetscDSGetCoordinateDimension(PetscDS prob, PetscInt *dimEmbed)
674: {
678:   *dimEmbed = prob->dimEmbed;
679:   return 0;
680: }

682: /*@
683:   PetscDSSetCoordinateDimension - Set the coordinate dimension of the `PetscDS`, meaning the dimension of the space into which the discretiaztions are embedded

685:   Logically collective on prob

687:   Input Parameters:
688: + prob - The `PetscDS` object
689: - dimEmbed - The coordinate dimension

691:   Level: beginner

693: .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
694: @*/
695: PetscErrorCode PetscDSSetCoordinateDimension(PetscDS prob, PetscInt dimEmbed)
696: {
699:   prob->dimEmbed = dimEmbed;
700:   return 0;
701: }

703: /*@
704:   PetscDSIsCohesive - Returns the flag indicating that this `PetscDS` is for a cohesive cell

706:   Not collective

708:   Input Parameter:
709: . ds - The `PetscDS` object

711:   Output Parameter:
712: . isCohesive - The flag

714:   Level: developer

716: .seealso: `PetscDS`, `PetscDSGetNumCohesive()`, `PetscDSGetCohesive()`, `PetscDSSetCohesive()`, `PetscDSCreate()`
717: @*/
718: PetscErrorCode PetscDSIsCohesive(PetscDS ds, PetscBool *isCohesive)
719: {
722:   *isCohesive = ds->isCohesive;
723:   return 0;
724: }

726: /*@
727:   PetscDSGetNumCohesive - Returns the numer of cohesive fields, meaning those defined on the interior of a cohesive cell

729:   Not collective

731:   Input Parameter:
732: . ds - The `PetscDS` object

734:   Output Parameter:
735: . numCohesive - The number of cohesive fields

737:   Level: developer

739: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSCreate()`
740: @*/
741: PetscErrorCode PetscDSGetNumCohesive(PetscDS ds, PetscInt *numCohesive)
742: {
743:   PetscInt f;

747:   *numCohesive = 0;
748:   for (f = 0; f < ds->Nf; ++f) *numCohesive += ds->cohesive[f] ? 1 : 0;
749:   return 0;
750: }

752: /*@
753:   PetscDSGetCohesive - Returns the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell

755:   Not collective

757:   Input Parameters:
758: + ds - The `PetscDS` object
759: - f  - The field index

761:   Output Parameter:
762: . isCohesive - The flag

764:   Level: developer

766: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
767: @*/
768: PetscErrorCode PetscDSGetCohesive(PetscDS ds, PetscInt f, PetscBool *isCohesive)
769: {
773:   *isCohesive = ds->cohesive[f];
774:   return 0;
775: }

777: /*@
778:   PetscDSSetCohesive - Set the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell

780:   Not collective

782:   Input Parameters:
783: + ds - The `PetscDS` object
784: . f  - The field index
785: - isCohesive - The flag for a cohesive field

787:   Level: developer

789: .seealso: `PetscDS`, `PetscDSGetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
790: @*/
791: PetscErrorCode PetscDSSetCohesive(PetscDS ds, PetscInt f, PetscBool isCohesive)
792: {
793:   PetscInt i;

797:   ds->cohesive[f] = isCohesive;
798:   ds->isCohesive  = PETSC_FALSE;
799:   for (i = 0; i < ds->Nf; ++i) ds->isCohesive = ds->isCohesive || ds->cohesive[f] ? PETSC_TRUE : PETSC_FALSE;
800:   return 0;
801: }

803: /*@
804:   PetscDSGetTotalDimension - Returns the total size of the approximation space for this system

806:   Not collective

808:   Input Parameter:
809: . prob - The `PetscDS` object

811:   Output Parameter:
812: . dim - The total problem dimension

814:   Level: beginner

816: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
817: @*/
818: PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim)
819: {
821:   PetscDSSetUp(prob);
823:   *dim = prob->totDim;
824:   return 0;
825: }

827: /*@
828:   PetscDSGetTotalComponents - Returns the total number of components in this system

830:   Not collective

832:   Input Parameter:
833: . prob - The `PetscDS` object

835:   Output Parameter:
836: . dim - The total number of components

838:   Level: beginner

840: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
841: @*/
842: PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc)
843: {
845:   PetscDSSetUp(prob);
847:   *Nc = prob->totComp;
848:   return 0;
849: }

851: /*@
852:   PetscDSGetDiscretization - Returns the discretization object for the given field

854:   Not collective

856:   Input Parameters:
857: + prob - The `PetscDS` object
858: - f - The field number

860:   Output Parameter:
861: . disc - The discretization object

863:   Level: beginner

865: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
866: @*/
867: PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc)
868: {
873:   *disc = prob->disc[f];
874:   return 0;
875: }

877: /*@
878:   PetscDSSetDiscretization - Sets the discretization object for the given field

880:   Not collective

882:   Input Parameters:
883: + prob - The `PetscDS` object
884: . f - The field number
885: - disc - The discretization object

887:   Level: beginner

889: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSGetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
890: @*/
891: PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc)
892: {
896:   PetscDSEnlarge_Static(prob, f + 1);
897:   PetscObjectDereference(prob->disc[f]);
898:   prob->disc[f] = disc;
899:   PetscObjectReference(disc);
900:   if (disc) {
901:     PetscClassId id;

903:     PetscObjectGetClassId(disc, &id);
904:     if (id == PETSCFE_CLASSID) {
905:       PetscDSSetImplicit(prob, f, PETSC_TRUE);
906:     } else if (id == PETSCFV_CLASSID) {
907:       PetscDSSetImplicit(prob, f, PETSC_FALSE);
908:     }
909:     PetscDSSetJetDegree(prob, f, 1);
910:   }
911:   return 0;
912: }

914: /*@
915:   PetscDSGetWeakForm - Returns the weak form object

917:   Not collective

919:   Input Parameter:
920: . ds - The `PetscDS` object

922:   Output Parameter:
923: . wf - The weak form object

925:   Level: beginner

927: .seealso: `PetscWeakForm`, `PetscDSSetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
928: @*/
929: PetscErrorCode PetscDSGetWeakForm(PetscDS ds, PetscWeakForm *wf)
930: {
933:   *wf = ds->wf;
934:   return 0;
935: }

937: /*@
938:   PetscDSSetWeakForm - Sets the weak form object

940:   Not collective

942:   Input Parameters:
943: + ds - The `PetscDS` object
944: - wf - The weak form object

946:   Level: beginner

948: .seealso: `PetscWeakForm`, `PetscDSGetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
949: @*/
950: PetscErrorCode PetscDSSetWeakForm(PetscDS ds, PetscWeakForm wf)
951: {
954:   PetscObjectDereference((PetscObject)ds->wf);
955:   ds->wf = wf;
956:   PetscObjectReference((PetscObject)wf);
957:   PetscWeakFormSetNumFields(wf, ds->Nf);
958:   return 0;
959: }

961: /*@
962:   PetscDSAddDiscretization - Adds a discretization object

964:   Not collective

966:   Input Parameters:
967: + prob - The `PetscDS` object
968: - disc - The boundary discretization object

970:   Level: beginner

972: .seealso: `PetscWeakForm`, `PetscDSGetDiscretization()`, `PetscDSSetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
973: @*/
974: PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc)
975: {
976:   PetscDSSetDiscretization(prob, prob->Nf, disc);
977:   return 0;
978: }

980: /*@
981:   PetscDSGetQuadrature - Returns the quadrature, which must agree for all fields in the `PetscDS`

983:   Not collective

985:   Input Parameter:
986: . prob - The `PetscDS` object

988:   Output Parameter:
989: . q - The quadrature object

991:   Level: intermediate

993: .seealso: `PetscDS`, `PetscQuadrature`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
994: @*/
995: PetscErrorCode PetscDSGetQuadrature(PetscDS prob, PetscQuadrature *q)
996: {
997:   PetscObject  obj;
998:   PetscClassId id;

1000:   *q = NULL;
1001:   if (!prob->Nf) return 0;
1002:   PetscDSGetDiscretization(prob, 0, &obj);
1003:   PetscObjectGetClassId(obj, &id);
1004:   if (id == PETSCFE_CLASSID) PetscFEGetQuadrature((PetscFE)obj, q);
1005:   else if (id == PETSCFV_CLASSID) PetscFVGetQuadrature((PetscFV)obj, q);
1006:   else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
1007:   return 0;
1008: }

1010: /*@
1011:   PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for `TSIMEX`

1013:   Not collective

1015:   Input Parameters:
1016: + prob - The `PetscDS` object
1017: - f - The field number

1019:   Output Parameter:
1020: . implicit - The flag indicating what kind of solve to use for this field

1022:   Level: developer

1024: .seealso: `TSIMEX`, `PetscDS`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1025: @*/
1026: PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit)
1027: {
1031:   *implicit = prob->implicit[f];
1032:   return 0;
1033: }

1035: /*@
1036:   PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for `TSIMEX`

1038:   Not collective

1040:   Input Parameters:
1041: + prob - The `PetscDS` object
1042: . f - The field number
1043: - implicit - The flag indicating what kind of solve to use for this field

1045:   Level: developer

1047: .seealso: `TSIMEX`, `PetscDSGetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1048: @*/
1049: PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit)
1050: {
1053:   prob->implicit[f] = implicit;
1054:   return 0;
1055: }

1057: /*@
1058:   PetscDSGetJetDegree - Returns the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate.

1060:   Not collective

1062:   Input Parameters:
1063: + ds - The `PetscDS` object
1064: - f  - The field number

1066:   Output Parameter:
1067: . k  - The highest derivative we need to tabulate

1069:   Level: developer

1071: .seealso: `PetscDS`, `PetscDSSetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1072: @*/
1073: PetscErrorCode PetscDSGetJetDegree(PetscDS ds, PetscInt f, PetscInt *k)
1074: {
1078:   *k = ds->jetDegree[f];
1079:   return 0;
1080: }

1082: /*@
1083:   PetscDSSetJetDegree - Set the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate.

1085:   Not collective

1087:   Input Parameters:
1088: + ds - The `PetscDS` object
1089: . f  - The field number
1090: - k  - The highest derivative we need to tabulate

1092:   Level: developer

1094: .seealso: ``PetscDS`, PetscDSGetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1095: @*/
1096: PetscErrorCode PetscDSSetJetDegree(PetscDS ds, PetscInt f, PetscInt k)
1097: {
1100:   ds->jetDegree[f] = k;
1101:   return 0;
1102: }

1104: PetscErrorCode PetscDSGetObjective(PetscDS ds, PetscInt f, void (**obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
1105: {
1106:   PetscPointFunc *tmp;
1107:   PetscInt        n;

1112:   PetscWeakFormGetObjective(ds->wf, NULL, 0, f, 0, &n, &tmp);
1113:   *obj = tmp ? tmp[0] : NULL;
1114:   return 0;
1115: }

1117: PetscErrorCode PetscDSSetObjective(PetscDS ds, PetscInt f, void (*obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
1118: {
1122:   PetscWeakFormSetIndexObjective(ds->wf, NULL, 0, f, 0, 0, obj);
1123:   return 0;
1124: }

1126: /*@C
1127:   PetscDSGetResidual - Get the pointwise residual function for a given test field

1129:   Not collective

1131:   Input Parameters:
1132: + ds - The `PetscDS`
1133: - f  - The test field number

1135:   Output Parameters:
1136: + f0 - integrand for the test function term
1137: - f1 - integrand for the test function gradient term

1139:   Calling sequence for the callbacks f0 and f1:
1140: .vb
1141:   f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1142:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1143:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1144:      PetscReal t, const PetscReal x[], PetscScalar f0[])
1145: .ve
1146: + dim - the spatial dimension
1147: . Nf - the number of fields
1148: . uOff - the offset into u[] and u_t[] for each field
1149: . uOff_x - the offset into u_x[] for each field
1150: . u - each field evaluated at the current point
1151: . u_t - the time derivative of each field evaluated at the current point
1152: . u_x - the gradient of each field evaluated at the current point
1153: . aOff - the offset into a[] and a_t[] for each auxiliary field
1154: . aOff_x - the offset into a_x[] for each auxiliary field
1155: . a - each auxiliary field evaluated at the current point
1156: . a_t - the time derivative of each auxiliary field evaluated at the current point
1157: . a_x - the gradient of auxiliary each field evaluated at the current point
1158: . t - current time
1159: . x - coordinates of the current point
1160: . numConstants - number of constant parameters
1161: . constants - constant parameters
1162: - f0 - output values at the current point

1164:   Level: intermediate

1166:   Note:
1167:   We are using a first order FEM model for the weak form:  \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t)

1169: .seealso: `PetscDS`, `PetscDSSetResidual()`
1170: @*/
1171: PetscErrorCode PetscDSGetResidual(PetscDS ds, PetscInt f, void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1172: {
1173:   PetscPointFunc *tmp0, *tmp1;
1174:   PetscInt        n0, n1;

1178:   PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1);
1179:   *f0 = tmp0 ? tmp0[0] : NULL;
1180:   *f1 = tmp1 ? tmp1[0] : NULL;
1181:   return 0;
1182: }

1184: /*@C
1185:   PetscDSSetResidual - Set the pointwise residual function for a given test field

1187:   Not collective

1189:   Input Parameters:
1190: + ds - The `PetscDS`
1191: . f  - The test field number
1192: . f0 - integrand for the test function term
1193: - f1 - integrand for the test function gradient term

1195:   Calling sequence for the callbacks f0 and f1:
1196: .vb
1197:   f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1198:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1199:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1200:      PetscReal t, const PetscReal x[], PetscScalar f0[])
1201: .ve
1202: + dim - the spatial dimension
1203: . Nf - the number of fields
1204: . uOff - the offset into u[] and u_t[] for each field
1205: . uOff_x - the offset into u_x[] for each field
1206: . u - each field evaluated at the current point
1207: . u_t - the time derivative of each field evaluated at the current point
1208: . u_x - the gradient of each field evaluated at the current point
1209: . aOff - the offset into a[] and a_t[] for each auxiliary field
1210: . aOff_x - the offset into a_x[] for each auxiliary field
1211: . a - each auxiliary field evaluated at the current point
1212: . a_t - the time derivative of each auxiliary field evaluated at the current point
1213: . a_x - the gradient of auxiliary each field evaluated at the current point
1214: . t - current time
1215: . x - coordinates of the current point
1216: . numConstants - number of constant parameters
1217: . constants - constant parameters
1218: - f0 - output values at the current point

1220:   Level: intermediate

1222:   Note:
1223:   We are using a first order FEM model for the weak form:  \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t)

1225: .seealso: `PetscDS`, `PetscDSGetResidual()`
1226: @*/
1227: PetscErrorCode PetscDSSetResidual(PetscDS ds, PetscInt f, void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1228: {
1233:   PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1);
1234:   return 0;
1235: }

1237: /*@C
1238:   PetscDSGetRHSResidual - Get the pointwise RHS residual function for explicit timestepping for a given test field

1240:   Not collective

1242:   Input Parameters:
1243: + ds - The `PetscDS`
1244: - f  - The test field number

1246:   Output Parameters:
1247: + f0 - integrand for the test function term
1248: - f1 - integrand for the test function gradient term

1250:   Calling sequence for the callbacks f0 and f1:
1251: .vb
1252:   f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1253:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1254:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1255:     PetscReal t, const PetscReal x[], PetscScalar f0[])
1256: .ve
1257: + dim - the spatial dimension
1258: . Nf - the number of fields
1259: . uOff - the offset into u[] and u_t[] for each field
1260: . uOff_x - the offset into u_x[] for each field
1261: . u - each field evaluated at the current point
1262: . u_t - the time derivative of each field evaluated at the current point
1263: . u_x - the gradient of each field evaluated at the current point
1264: . aOff - the offset into a[] and a_t[] for each auxiliary field
1265: . aOff_x - the offset into a_x[] for each auxiliary field
1266: . a - each auxiliary field evaluated at the current point
1267: . a_t - the time derivative of each auxiliary field evaluated at the current point
1268: . a_x - the gradient of auxiliary each field evaluated at the current point
1269: . t - current time
1270: . x - coordinates of the current point
1271: . numConstants - number of constant parameters
1272: . constants - constant parameters
1273: - f0 - output values at the current point

1275:   Level: intermediate

1277:   Note:
1278:   We are using a first order FEM model for the weak form: \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t)

1280: .seealso: `PetscDS`, `PetscDSSetRHSResidual()`
1281: @*/
1282: PetscErrorCode PetscDSGetRHSResidual(PetscDS ds, PetscInt f, void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1283: {
1284:   PetscPointFunc *tmp0, *tmp1;
1285:   PetscInt        n0, n1;

1289:   PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 100, &n0, &tmp0, &n1, &tmp1);
1290:   *f0 = tmp0 ? tmp0[0] : NULL;
1291:   *f1 = tmp1 ? tmp1[0] : NULL;
1292:   return 0;
1293: }

1295: /*@C
1296:   PetscDSSetRHSResidual - Set the pointwise residual function for explicit timestepping for a given test field

1298:   Not collective

1300:   Input Parameters:
1301: + ds - The `PetscDS`
1302: . f  - The test field number
1303: . f0 - integrand for the test function term
1304: - f1 - integrand for the test function gradient term

1306:   Clling sequence for the callbacks f0 and f1:
1307: .vb
1308:   f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1309:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1310:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1311:      PetscReal t, const PetscReal x[], PetscScalar f0[])
1312: .ve
1313: + dim - the spatial dimension
1314: . Nf - the number of fields
1315: . uOff - the offset into u[] and u_t[] for each field
1316: . uOff_x - the offset into u_x[] for each field
1317: . u - each field evaluated at the current point
1318: . u_t - the time derivative of each field evaluated at the current point
1319: . u_x - the gradient of each field evaluated at the current point
1320: . aOff - the offset into a[] and a_t[] for each auxiliary field
1321: . aOff_x - the offset into a_x[] for each auxiliary field
1322: . a - each auxiliary field evaluated at the current point
1323: . a_t - the time derivative of each auxiliary field evaluated at the current point
1324: . a_x - the gradient of auxiliary each field evaluated at the current point
1325: . t - current time
1326: . x - coordinates of the current point
1327: . numConstants - number of constant parameters
1328: . constants - constant parameters
1329: - f0 - output values at the current point

1331:   Level: intermediate

1333:   Note:
1334:   We are using a first order FEM model for the weak form:  \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t)

1336: .seealso: `PetscDS`, `PetscDSGetResidual()`
1337: @*/
1338: PetscErrorCode PetscDSSetRHSResidual(PetscDS ds, PetscInt f, void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1339: {
1344:   PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 100, 0, f0, 0, f1);
1345:   return 0;
1346: }

1348: /*@C
1349:   PetscDSHasJacobian - Checks that the Jacobian functions have been set

1351:   Not collective

1353:   Input Parameter:
1354: . prob - The `PetscDS`

1356:   Output Parameter:
1357: . hasJac - flag that pointwise function for the Jacobian has been set

1359:   Level: intermediate

1361: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1362: @*/
1363: PetscErrorCode PetscDSHasJacobian(PetscDS ds, PetscBool *hasJac)
1364: {
1366:   PetscWeakFormHasJacobian(ds->wf, hasJac);
1367:   return 0;
1368: }

1370: /*@C
1371:   PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field

1373:   Not collective

1375:   Input Parameters:
1376: + ds - The `PetscDS`
1377: . f  - The test field number
1378: - g  - The field number

1380:   Output Parameters:
1381: + g0 - integrand for the test and basis function term
1382: . g1 - integrand for the test function and basis function gradient term
1383: . g2 - integrand for the test function gradient and basis function term
1384: - g3 - integrand for the test function gradient and basis function gradient term

1386:   Calling sequence for the callbacks g0, g1, g2 and g3:
1387: .vb
1388:   g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1389:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1390:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1391:      PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1392: .ve
1393: + dim - the spatial dimension
1394: . Nf - the number of fields
1395: . uOff - the offset into u[] and u_t[] for each field
1396: . uOff_x - the offset into u_x[] for each field
1397: . u - each field evaluated at the current point
1398: . u_t - the time derivative of each field evaluated at the current point
1399: . u_x - the gradient of each field evaluated at the current point
1400: . aOff - the offset into a[] and a_t[] for each auxiliary field
1401: . aOff_x - the offset into a_x[] for each auxiliary field
1402: . a - each auxiliary field evaluated at the current point
1403: . a_t - the time derivative of each auxiliary field evaluated at the current point
1404: . a_x - the gradient of auxiliary each field evaluated at the current point
1405: . t - current time
1406: . u_tShift - the multiplier a for dF/dU_t
1407: . x - coordinates of the current point
1408: . numConstants - number of constant parameters
1409: . constants - constant parameters
1410: - g0 - output values at the current point

1412:   Level: intermediate

1414:   Note:
1415:   We are using a first order FEM model for the weak form:
1416:   \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi

1418: .seealso: `PetscDS`, `PetscDSSetJacobian()`
1419: @*/
1420: PetscErrorCode PetscDSGetJacobian(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1421: {
1422:   PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1423:   PetscInt       n0, n1, n2, n3;

1428:   PetscWeakFormGetJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3);
1429:   *g0 = tmp0 ? tmp0[0] : NULL;
1430:   *g1 = tmp1 ? tmp1[0] : NULL;
1431:   *g2 = tmp2 ? tmp2[0] : NULL;
1432:   *g3 = tmp3 ? tmp3[0] : NULL;
1433:   return 0;
1434: }

1436: /*@C
1437:   PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields

1439:   Not collective

1441:   Input Parameters:
1442: + ds - The `PetscDS`
1443: . f  - The test field number
1444: . g  - The field number
1445: . g0 - integrand for the test and basis function term
1446: . g1 - integrand for the test function and basis function gradient term
1447: . g2 - integrand for the test function gradient and basis function term
1448: - g3 - integrand for the test function gradient and basis function gradient term

1450:   Calling sequence for the callbacks g0, g1, g2 and g3:
1451: .vb
1452:   g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1453:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1454:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1455:      PetscReal t, const PetscReal x[], PetscScalar g0[])
1456: .ve
1457: + dim - the spatial dimension
1458: . Nf - the number of fields
1459: . uOff - the offset into u[] and u_t[] for each field
1460: . uOff_x - the offset into u_x[] for each field
1461: . u - each field evaluated at the current point
1462: . u_t - the time derivative of each field evaluated at the current point
1463: . u_x - the gradient of each field evaluated at the current point
1464: . aOff - the offset into a[] and a_t[] for each auxiliary field
1465: . aOff_x - the offset into a_x[] for each auxiliary field
1466: . a - each auxiliary field evaluated at the current point
1467: . a_t - the time derivative of each auxiliary field evaluated at the current point
1468: . a_x - the gradient of auxiliary each field evaluated at the current point
1469: . t - current time
1470: . u_tShift - the multiplier a for dF/dU_t
1471: . x - coordinates of the current point
1472: . numConstants - number of constant parameters
1473: . constants - constant parameters
1474: - g0 - output values at the current point

1476:   Level: intermediate

1478:   Note:
1479:    We are using a first order FEM model for the weak form:
1480:   \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi

1482: .seealso: `PetscDS`, `PetscDSGetJacobian()`
1483: @*/
1484: PetscErrorCode PetscDSSetJacobian(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1485: {
1493:   PetscWeakFormSetIndexJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3);
1494:   return 0;
1495: }

1497: /*@C
1498:   PetscDSUseJacobianPreconditioner - Set whether to construct a Jacobian preconditioner

1500:   Not collective

1502:   Input Parameters:
1503: + prob - The `PetscDS`
1504: - useJacPre - flag that enables construction of a Jacobian preconditioner

1506:   Level: intermediate

1508: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1509: @*/
1510: PetscErrorCode PetscDSUseJacobianPreconditioner(PetscDS prob, PetscBool useJacPre)
1511: {
1513:   prob->useJacPre = useJacPre;
1514:   return 0;
1515: }

1517: /*@C
1518:   PetscDSHasJacobianPreconditioner - Checks if a Jacobian preconditioner matrix has been set

1520:   Not collective

1522:   Input Parameter:
1523: . prob - The `PetscDS`

1525:   Output Parameter:
1526: . hasJacPre - flag that pointwise function for Jacobian preconditioner matrix has been set

1528:   Level: intermediate

1530: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1531: @*/
1532: PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS ds, PetscBool *hasJacPre)
1533: {
1535:   *hasJacPre = PETSC_FALSE;
1536:   if (!ds->useJacPre) return 0;
1537:   PetscWeakFormHasJacobianPreconditioner(ds->wf, hasJacPre);
1538:   return 0;
1539: }

1541: /*@C
1542:   PetscDSGetJacobianPreconditioner - Get the pointwise Jacobian preconditioner function for given test and basis field. If this is missing,
1543:    the system matrix is used to build the preconditioner.

1545:   Not collective

1547:   Input Parameters:
1548: + ds - The `PetscDS`
1549: . f  - The test field number
1550: - g  - The field number

1552:   Output Parameters:
1553: + g0 - integrand for the test and basis function term
1554: . g1 - integrand for the test function and basis function gradient term
1555: . g2 - integrand for the test function gradient and basis function term
1556: - g3 - integrand for the test function gradient and basis function gradient term

1558:   Calling sequence for the callbacks g0, g1, g2 and g3:
1559: .vb
1560:   g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1561:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1562:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1563:      PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1564: .ve
1565: + dim - the spatial dimension
1566: . Nf - the number of fields
1567: . uOff - the offset into u[] and u_t[] for each field
1568: . uOff_x - the offset into u_x[] for each field
1569: . u - each field evaluated at the current point
1570: . u_t - the time derivative of each field evaluated at the current point
1571: . u_x - the gradient of each field evaluated at the current point
1572: . aOff - the offset into a[] and a_t[] for each auxiliary field
1573: . aOff_x - the offset into a_x[] for each auxiliary field
1574: . a - each auxiliary field evaluated at the current point
1575: . a_t - the time derivative of each auxiliary field evaluated at the current point
1576: . a_x - the gradient of auxiliary each field evaluated at the current point
1577: . t - current time
1578: . u_tShift - the multiplier a for dF/dU_t
1579: . x - coordinates of the current point
1580: . numConstants - number of constant parameters
1581: . constants - constant parameters
1582: - g0 - output values at the current point

1584:   Level: intermediate

1586:   Note:
1587:   We are using a first order FEM model for the weak form:
1588:   \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi

1590: .seealso: `PetscDS`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1591: @*/
1592: PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1593: {
1594:   PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1595:   PetscInt       n0, n1, n2, n3;

1600:   PetscWeakFormGetJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3);
1601:   *g0 = tmp0 ? tmp0[0] : NULL;
1602:   *g1 = tmp1 ? tmp1[0] : NULL;
1603:   *g2 = tmp2 ? tmp2[0] : NULL;
1604:   *g3 = tmp3 ? tmp3[0] : NULL;
1605:   return 0;
1606: }

1608: /*@C
1609:   PetscDSSetJacobianPreconditioner - Set the pointwise Jacobian preconditioner function for given test and basis fields.
1610:   If this is missing, the system matrix is used to build the preconditioner.

1612:   Not collective

1614:   Input Parameters:
1615: + ds - The `PetscDS`
1616: . f  - The test field number
1617: . g  - The field number
1618: . g0 - integrand for the test and basis function term
1619: . g1 - integrand for the test function and basis function gradient term
1620: . g2 - integrand for the test function gradient and basis function term
1621: - g3 - integrand for the test function gradient and basis function gradient term

1623:   Calling sequence for the callbacks g0, g1, g2 and g3:
1624: .vb
1625:   g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1626:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1627:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1628:      PetscReal t, const PetscReal x[], PetscScalar g0[])
1629: .ve
1630: + dim - the spatial dimension
1631: . Nf - the number of fields
1632: . uOff - the offset into u[] and u_t[] for each field
1633: . uOff_x - the offset into u_x[] for each field
1634: . u - each field evaluated at the current point
1635: . u_t - the time derivative of each field evaluated at the current point
1636: . u_x - the gradient of each field evaluated at the current point
1637: . aOff - the offset into a[] and a_t[] for each auxiliary field
1638: . aOff_x - the offset into a_x[] for each auxiliary field
1639: . a - each auxiliary field evaluated at the current point
1640: . a_t - the time derivative of each auxiliary field evaluated at the current point
1641: . a_x - the gradient of auxiliary each field evaluated at the current point
1642: . t - current time
1643: . u_tShift - the multiplier a for dF/dU_t
1644: . x - coordinates of the current point
1645: . numConstants - number of constant parameters
1646: . constants - constant parameters
1647: - g0 - output values at the current point

1649:   Level: intermediate

1651:   Note:
1652:   We are using a first order FEM model for the weak form:
1653:   \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi

1655: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobian()`
1656: @*/
1657: PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1658: {
1666:   PetscWeakFormSetIndexJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3);
1667:   return 0;
1668: }

1670: /*@C
1671:   PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, dF/du_t, has been set

1673:   Not collective

1675:   Input Parameter:
1676: . ds - The `PetscDS`

1678:   Output Parameter:
1679: . hasDynJac - flag that pointwise function for dynamic Jacobian has been set

1681:   Level: intermediate

1683: .seealso: `PetscDS`, `PetscDSGetDynamicJacobian()`, `PetscDSSetDynamicJacobian()`, `PetscDSGetJacobian()`
1684: @*/
1685: PetscErrorCode PetscDSHasDynamicJacobian(PetscDS ds, PetscBool *hasDynJac)
1686: {
1688:   PetscWeakFormHasDynamicJacobian(ds->wf, hasDynJac);
1689:   return 0;
1690: }

1692: /*@C
1693:   PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, dF/du_t, function for given test and basis field

1695:   Not collective

1697:   Input Parameters:
1698: + ds - The `PetscDS`
1699: . f  - The test field number
1700: - g  - The field number

1702:   Output Parameters:
1703: + g0 - integrand for the test and basis function term
1704: . g1 - integrand for the test function and basis function gradient term
1705: . g2 - integrand for the test function gradient and basis function term
1706: - g3 - integrand for the test function gradient and basis function gradient term

1708:    Calling sequence for the callbacks g0, g1, g2 and g3:
1709: .vb
1710:   g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1711:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1712:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1713:      PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1714: .ve
1715: + dim - the spatial dimension
1716: . Nf - the number of fields
1717: . uOff - the offset into u[] and u_t[] for each field
1718: . uOff_x - the offset into u_x[] for each field
1719: . u - each field evaluated at the current point
1720: . u_t - the time derivative of each field evaluated at the current point
1721: . u_x - the gradient of each field evaluated at the current point
1722: . aOff - the offset into a[] and a_t[] for each auxiliary field
1723: . aOff_x - the offset into a_x[] for each auxiliary field
1724: . a - each auxiliary field evaluated at the current point
1725: . a_t - the time derivative of each auxiliary field evaluated at the current point
1726: . a_x - the gradient of auxiliary each field evaluated at the current point
1727: . t - current time
1728: . u_tShift - the multiplier a for dF/dU_t
1729: . x - coordinates of the current point
1730: . numConstants - number of constant parameters
1731: . constants - constant parameters
1732: - g0 - output values at the current point

1734:   Level: intermediate

1736:   Note:
1737:   We are using a first order FEM model for the weak form:
1738:   \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi

1740: .seealso: `PetscDS`, `PetscDSSetJacobian()`
1741: @*/
1742: PetscErrorCode PetscDSGetDynamicJacobian(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1743: {
1744:   PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1745:   PetscInt       n0, n1, n2, n3;

1750:   PetscWeakFormGetDynamicJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3);
1751:   *g0 = tmp0 ? tmp0[0] : NULL;
1752:   *g1 = tmp1 ? tmp1[0] : NULL;
1753:   *g2 = tmp2 ? tmp2[0] : NULL;
1754:   *g3 = tmp3 ? tmp3[0] : NULL;
1755:   return 0;
1756: }

1758: /*@C
1759:   PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, dF/du_t, function for given test and basis fields

1761:   Not collective

1763:   Input Parameters:
1764: + ds - The `PetscDS`
1765: . f  - The test field number
1766: . g  - The field number
1767: . g0 - integrand for the test and basis function term
1768: . g1 - integrand for the test function and basis function gradient term
1769: . g2 - integrand for the test function gradient and basis function term
1770: - g3 - integrand for the test function gradient and basis function gradient term

1772:   Calling sequence for the callbacks g0, g1, g2 and g3:
1773: .vb
1774:    g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1775:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1776:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1777:      PetscReal t, const PetscReal x[], PetscScalar g0[])
1778: .ve
1779: + dim - the spatial dimension
1780: . Nf - the number of fields
1781: . uOff - the offset into u[] and u_t[] for each field
1782: . uOff_x - the offset into u_x[] for each field
1783: . u - each field evaluated at the current point
1784: . u_t - the time derivative of each field evaluated at the current point
1785: . u_x - the gradient of each field evaluated at the current point
1786: . aOff - the offset into a[] and a_t[] for each auxiliary field
1787: . aOff_x - the offset into a_x[] for each auxiliary field
1788: . a - each auxiliary field evaluated at the current point
1789: . a_t - the time derivative of each auxiliary field evaluated at the current point
1790: . a_x - the gradient of auxiliary each field evaluated at the current point
1791: . t - current time
1792: . u_tShift - the multiplier a for dF/dU_t
1793: . x - coordinates of the current point
1794: . numConstants - number of constant parameters
1795: . constants - constant parameters
1796: - g0 - output values at the current point

1798:   Level: intermediate

1800:   Note:
1801:   We are using a first order FEM model for the weak form:
1802:   \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi

1804: .seealso: `PetscDS`, `PetscDSGetJacobian()`
1805: @*/
1806: PetscErrorCode PetscDSSetDynamicJacobian(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1807: {
1815:   PetscWeakFormSetIndexDynamicJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3);
1816:   return 0;
1817: }

1819: /*@C
1820:   PetscDSGetRiemannSolver - Returns the Riemann solver for the given field

1822:   Not collective

1824:   Input Parameters:
1825: + ds - The `PetscDS` object
1826: - f  - The field number

1828:   Output Parameter:
1829: . r    - Riemann solver

1831:   Calling sequence for r:
1832: .vb
1833:   r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)
1834: .ve
1835: + dim  - The spatial dimension
1836: . Nf   - The number of fields
1837: . x    - The coordinates at a point on the interface
1838: . n    - The normal vector to the interface
1839: . uL   - The state vector to the left of the interface
1840: . uR   - The state vector to the right of the interface
1841: . flux - output array of flux through the interface
1842: . numConstants - number of constant parameters
1843: . constants - constant parameters
1844: - ctx  - optional user context

1846:   Level: intermediate

1848: .seealso: `PetscDS`, `PetscDSSetRiemannSolver()`
1849: @*/
1850: PetscErrorCode PetscDSGetRiemannSolver(PetscDS ds, PetscInt f, void (**r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscInt numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx))
1851: {
1852:   PetscRiemannFunc *tmp;
1853:   PetscInt          n;

1858:   PetscWeakFormGetRiemannSolver(ds->wf, NULL, 0, f, 0, &n, &tmp);
1859:   *r = tmp ? tmp[0] : NULL;
1860:   return 0;
1861: }

1863: /*@C
1864:   PetscDSSetRiemannSolver - Sets the Riemann solver for the given field

1866:   Not collective

1868:   Input Parameters:
1869: + ds - The `PetscDS` object
1870: . f  - The field number
1871: - r  - Riemann solver

1873:   Calling sequence for r:
1874: .vb
1875:    r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)
1876: .ve
1877: + dim  - The spatial dimension
1878: . Nf   - The number of fields
1879: . x    - The coordinates at a point on the interface
1880: . n    - The normal vector to the interface
1881: . uL   - The state vector to the left of the interface
1882: . uR   - The state vector to the right of the interface
1883: . flux - output array of flux through the interface
1884: . numConstants - number of constant parameters
1885: . constants - constant parameters
1886: - ctx  - optional user context

1888:   Level: intermediate

1890: .seealso: `PetscDS`, `PetscDSGetRiemannSolver()`
1891: @*/
1892: PetscErrorCode PetscDSSetRiemannSolver(PetscDS ds, PetscInt f, void (*r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscInt numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx))
1893: {
1897:   PetscWeakFormSetIndexRiemannSolver(ds->wf, NULL, 0, f, 0, 0, r);
1898:   return 0;
1899: }

1901: /*@C
1902:   PetscDSGetUpdate - Get the pointwise update function for a given field

1904:   Not collective

1906:   Input Parameters:
1907: + ds - The `PetscDS`
1908: - f  - The field number

1910:   Output Parameter:
1911: . update - update function

1913:   Calling sequence for the callback update:
1914: .vb
1915:   update(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1916:          const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1917:          const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1918:          PetscReal t, const PetscReal x[], PetscScalar uNew[])
1919: .ve
1920: + dim - the spatial dimension
1921: . Nf - the number of fields
1922: . uOff - the offset into u[] and u_t[] for each field
1923: . uOff_x - the offset into u_x[] for each field
1924: . u - each field evaluated at the current point
1925: . u_t - the time derivative of each field evaluated at the current point
1926: . u_x - the gradient of each field evaluated at the current point
1927: . aOff - the offset into a[] and a_t[] for each auxiliary field
1928: . aOff_x - the offset into a_x[] for each auxiliary field
1929: . a - each auxiliary field evaluated at the current point
1930: . a_t - the time derivative of each auxiliary field evaluated at the current point
1931: . a_x - the gradient of auxiliary each field evaluated at the current point
1932: . t - current time
1933: . x - coordinates of the current point
1934: - uNew - new value for field at the current point

1936:   Level: intermediate

1938: .seealso: `PetscDS`, `PetscDSSetUpdate()`, `PetscDSSetResidual()`
1939: @*/
1940: PetscErrorCode PetscDSGetUpdate(PetscDS ds, PetscInt f, void (**update)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
1941: {
1944:   if (update) {
1946:     *update = ds->update[f];
1947:   }
1948:   return 0;
1949: }

1951: /*@C
1952:   PetscDSSetUpdate - Set the pointwise update function for a given field

1954:   Not collective

1956:   Input Parameters:
1957: + ds     - The `PetscDS`
1958: . f      - The field number
1959: - update - update function

1961:   Calling sequence for the callback update:
1962: .vb
1963:   update(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1964:          const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1965:          const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1966:          PetscReal t, const PetscReal x[], PetscScalar uNew[])
1967: .ve
1968: + dim - the spatial dimension
1969: . Nf - the number of fields
1970: . uOff - the offset into u[] and u_t[] for each field
1971: . uOff_x - the offset into u_x[] for each field
1972: . u - each field evaluated at the current point
1973: . u_t - the time derivative of each field evaluated at the current point
1974: . u_x - the gradient of each field evaluated at the current point
1975: . aOff - the offset into a[] and a_t[] for each auxiliary field
1976: . aOff_x - the offset into a_x[] for each auxiliary field
1977: . a - each auxiliary field evaluated at the current point
1978: . a_t - the time derivative of each auxiliary field evaluated at the current point
1979: . a_x - the gradient of auxiliary each field evaluated at the current point
1980: . t - current time
1981: . x - coordinates of the current point
1982: - uNew - new field values at the current point

1984:   Level: intermediate

1986: .seealso: `PetscDS`, `PetscDSGetResidual()`
1987: @*/
1988: PetscErrorCode PetscDSSetUpdate(PetscDS ds, PetscInt f, void (*update)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
1989: {
1993:   PetscDSEnlarge_Static(ds, f + 1);
1994:   ds->update[f] = update;
1995:   return 0;
1996: }

1998: PetscErrorCode PetscDSGetContext(PetscDS ds, PetscInt f, void *ctx)
1999: {
2003:   *(void **)ctx = ds->ctx[f];
2004:   return 0;
2005: }

2007: PetscErrorCode PetscDSSetContext(PetscDS ds, PetscInt f, void *ctx)
2008: {
2011:   PetscDSEnlarge_Static(ds, f + 1);
2012:   ds->ctx[f] = ctx;
2013:   return 0;
2014: }

2016: /*@C
2017:   PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field

2019:   Not collective

2021:   Input Parameters:
2022: + ds - The PetscDS
2023: - f  - The test field number

2025:   Output Parameters:
2026: + f0 - boundary integrand for the test function term
2027: - f1 - boundary integrand for the test function gradient term

2029:   Calling sequence for the callbacks f0 and f1:
2030: .vb
2031:   f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2032:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2033:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2034:      PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[])
2035: .ve
2036: + dim - the spatial dimension
2037: . Nf - the number of fields
2038: . uOff - the offset into u[] and u_t[] for each field
2039: . uOff_x - the offset into u_x[] for each field
2040: . u - each field evaluated at the current point
2041: . u_t - the time derivative of each field evaluated at the current point
2042: . u_x - the gradient of each field evaluated at the current point
2043: . aOff - the offset into a[] and a_t[] for each auxiliary field
2044: . aOff_x - the offset into a_x[] for each auxiliary field
2045: . a - each auxiliary field evaluated at the current point
2046: . a_t - the time derivative of each auxiliary field evaluated at the current point
2047: . a_x - the gradient of auxiliary each field evaluated at the current point
2048: . t - current time
2049: . x - coordinates of the current point
2050: . n - unit normal at the current point
2051: . numConstants - number of constant parameters
2052: . constants - constant parameters
2053: - f0 - output values at the current point

2055:   Level: intermediate

2057:   Note:
2058:   We are using a first order FEM model for the weak form:
2059:   \int_\Gamma \phi {\vec f}_0(u, u_t, \nabla u, x, t) \cdot \hat n + \nabla\phi \cdot {\overleftrightarrow f}_1(u, u_t, \nabla u, x, t) \cdot \hat n

2061: .seealso: `PetscDS`, `PetscDSSetBdResidual()`
2062: @*/
2063: PetscErrorCode PetscDSGetBdResidual(PetscDS ds, PetscInt f, void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
2064: {
2065:   PetscBdPointFunc *tmp0, *tmp1;
2066:   PetscInt          n0, n1;

2070:   PetscWeakFormGetBdResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1);
2071:   *f0 = tmp0 ? tmp0[0] : NULL;
2072:   *f1 = tmp1 ? tmp1[0] : NULL;
2073:   return 0;
2074: }

2076: /*@C
2077:   PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field

2079:   Not collective

2081:   Input Parameters:
2082: + ds - The `PetscDS`
2083: . f  - The test field number
2084: . f0 - boundary integrand for the test function term
2085: - f1 - boundary integrand for the test function gradient term

2087:   Calling sequence for the callbacks f0 and f1:
2088: .vb
2089:   f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2090:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2091:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2092:      PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[])
2093: .ve
2094: + dim - the spatial dimension
2095: . Nf - the number of fields
2096: . uOff - the offset into u[] and u_t[] for each field
2097: . uOff_x - the offset into u_x[] for each field
2098: . u - each field evaluated at the current point
2099: . u_t - the time derivative of each field evaluated at the current point
2100: . u_x - the gradient of each field evaluated at the current point
2101: . aOff - the offset into a[] and a_t[] for each auxiliary field
2102: . aOff_x - the offset into a_x[] for each auxiliary field
2103: . a - each auxiliary field evaluated at the current point
2104: . a_t - the time derivative of each auxiliary field evaluated at the current point
2105: . a_x - the gradient of auxiliary each field evaluated at the current point
2106: . t - current time
2107: . x - coordinates of the current point
2108: . n - unit normal at the current point
2109: . numConstants - number of constant parameters
2110: . constants - constant parameters
2111: - f0 - output values at the current point

2113:   Level: intermediate

2115:   Note:
2116:   We are using a first order FEM model for the weak form:
2117:   \int_\Gamma \phi {\vec f}_0(u, u_t, \nabla u, x, t) \cdot \hat n + \nabla\phi \cdot {\overleftrightarrow f}_1(u, u_t, \nabla u, x, t) \cdot \hat n

2119: .seealso: `PetscDS`, `PetscDSGetBdResidual()`
2120: @*/
2121: PetscErrorCode PetscDSSetBdResidual(PetscDS ds, PetscInt f, void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
2122: {
2125:   PetscWeakFormSetIndexBdResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1);
2126:   return 0;
2127: }

2129: /*@
2130:   PetscDSHasBdJacobian - Indicates that boundary Jacobian functions have been set

2132:   Not collective

2134:   Input Parameter:
2135: . ds - The `PetscDS`

2137:   Output Parameter:
2138: . hasBdJac - flag that pointwise function for the boundary Jacobian has been set

2140:   Level: intermediate

2142: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2143: @*/
2144: PetscErrorCode PetscDSHasBdJacobian(PetscDS ds, PetscBool *hasBdJac)
2145: {
2148:   PetscWeakFormHasBdJacobian(ds->wf, hasBdJac);
2149:   return 0;
2150: }

2152: /*@C
2153:   PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field

2155:   Not collective

2157:   Input Parameters:
2158: + ds - The `PetscDS`
2159: . f  - The test field number
2160: - g  - The field number

2162:   Output Parameters:
2163: + g0 - integrand for the test and basis function term
2164: . g1 - integrand for the test function and basis function gradient term
2165: . g2 - integrand for the test function gradient and basis function term
2166: - g3 - integrand for the test function gradient and basis function gradient term

2168:   Calling sequence for the callbacks g0, g1, g2 and g3:
2169: .vb
2170:   g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2171:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2172:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2173:      PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[])
2174: .ve
2175: + dim - the spatial dimension
2176: . Nf - the number of fields
2177: . uOff - the offset into u[] and u_t[] for each field
2178: . uOff_x - the offset into u_x[] for each field
2179: . u - each field evaluated at the current point
2180: . u_t - the time derivative of each field evaluated at the current point
2181: . u_x - the gradient of each field evaluated at the current point
2182: . aOff - the offset into a[] and a_t[] for each auxiliary field
2183: . aOff_x - the offset into a_x[] for each auxiliary field
2184: . a - each auxiliary field evaluated at the current point
2185: . a_t - the time derivative of each auxiliary field evaluated at the current point
2186: . a_x - the gradient of auxiliary each field evaluated at the current point
2187: . t - current time
2188: . u_tShift - the multiplier a for dF/dU_t
2189: . x - coordinates of the current point
2190: . n - normal at the current point
2191: . numConstants - number of constant parameters
2192: . constants - constant parameters
2193: - g0 - output values at the current point

2195:   Level: intermediate

2197:   Note:
2198:   We are using a first order FEM model for the weak form:
2199:   \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi

2201: .seealso: `PetscDS`, `PetscDSSetBdJacobian()`
2202: @*/
2203: PetscErrorCode PetscDSGetBdJacobian(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2204: {
2205:   PetscBdPointJac *tmp0, *tmp1, *tmp2, *tmp3;
2206:   PetscInt         n0, n1, n2, n3;

2211:   PetscWeakFormGetBdJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3);
2212:   *g0 = tmp0 ? tmp0[0] : NULL;
2213:   *g1 = tmp1 ? tmp1[0] : NULL;
2214:   *g2 = tmp2 ? tmp2[0] : NULL;
2215:   *g3 = tmp3 ? tmp3[0] : NULL;
2216:   return 0;
2217: }

2219: /*@C
2220:   PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field

2222:   Not collective

2224:   Input Parameters:
2225: + ds - The PetscDS
2226: . f  - The test field number
2227: . g  - The field number
2228: . g0 - integrand for the test and basis function term
2229: . g1 - integrand for the test function and basis function gradient term
2230: . g2 - integrand for the test function gradient and basis function term
2231: - g3 - integrand for the test function gradient and basis function gradient term

2233:   Calling sequence for the callbacks g0, g1, g2 and g3:
2234: .vb
2235:   g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2236:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2237:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2238:      PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[])
2239: .ve
2240: + dim - the spatial dimension
2241: . Nf - the number of fields
2242: . uOff - the offset into u[] and u_t[] for each field
2243: . uOff_x - the offset into u_x[] for each field
2244: . u - each field evaluated at the current point
2245: . u_t - the time derivative of each field evaluated at the current point
2246: . u_x - the gradient of each field evaluated at the current point
2247: . aOff - the offset into a[] and a_t[] for each auxiliary field
2248: . aOff_x - the offset into a_x[] for each auxiliary field
2249: . a - each auxiliary field evaluated at the current point
2250: . a_t - the time derivative of each auxiliary field evaluated at the current point
2251: . a_x - the gradient of auxiliary each field evaluated at the current point
2252: . t - current time
2253: . u_tShift - the multiplier a for dF/dU_t
2254: . x - coordinates of the current point
2255: . n - normal at the current point
2256: . numConstants - number of constant parameters
2257: . constants - constant parameters
2258: - g0 - output values at the current point

2260:   Level: intermediate

2262:   Note:
2263:   We are using a first order FEM model for the weak form:
2264:   \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi

2266: .seealso: `PetscDS`, `PetscDSGetBdJacobian()`
2267: @*/
2268: PetscErrorCode PetscDSSetBdJacobian(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2269: {
2277:   PetscWeakFormSetIndexBdJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3);
2278:   return 0;
2279: }

2281: /*@
2282:   PetscDSHasBdJacobianPreconditioner - Signals that boundary Jacobian preconditioner functions have been set

2284:   Not collective

2286:   Input Parameter:
2287: . ds - The `PetscDS`

2289:   Output Parameter:
2290: . hasBdJac - flag that pointwise function for the boundary Jacobian preconditioner has been set

2292:   Level: intermediate

2294: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2295: @*/
2296: PetscErrorCode PetscDSHasBdJacobianPreconditioner(PetscDS ds, PetscBool *hasBdJacPre)
2297: {
2300:   PetscWeakFormHasBdJacobianPreconditioner(ds->wf, hasBdJacPre);
2301:   return 0;
2302: }

2304: /*@C
2305:   PetscDSGetBdJacobianPreconditioner - Get the pointwise boundary Jacobian preconditioner function for given test and basis field

2307:   Not collective

2309:   Input Parameters:
2310: + ds - The `PetscDS`
2311: . f  - The test field number
2312: - g  - The field number

2314:   Output Parameters:
2315: + g0 - integrand for the test and basis function term
2316: . g1 - integrand for the test function and basis function gradient term
2317: . g2 - integrand for the test function gradient and basis function term
2318: - g3 - integrand for the test function gradient and basis function gradient term

2320:    Calling sequence for the callbacks g0, g1, g2 and g3:
2321: .vb
2322:   g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2323:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2324:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2325:      PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[])
2326: .ve
2327: + dim - the spatial dimension
2328: . Nf - the number of fields
2329: . NfAux - the number of auxiliary fields
2330: . uOff - the offset into u[] and u_t[] for each field
2331: . uOff_x - the offset into u_x[] for each field
2332: . u - each field evaluated at the current point
2333: . u_t - the time derivative of each field evaluated at the current point
2334: . u_x - the gradient of each field evaluated at the current point
2335: . aOff - the offset into a[] and a_t[] for each auxiliary field
2336: . aOff_x - the offset into a_x[] for each auxiliary field
2337: . a - each auxiliary field evaluated at the current point
2338: . a_t - the time derivative of each auxiliary field evaluated at the current point
2339: . a_x - the gradient of auxiliary each field evaluated at the current point
2340: . t - current time
2341: . u_tShift - the multiplier a for dF/dU_t
2342: . x - coordinates of the current point
2343: . n - normal at the current point
2344: . numConstants - number of constant parameters
2345: . constants - constant parameters
2346: - g0 - output values at the current point

2348:   Level: intermediate

2350:   Note:
2351:   We are using a first order FEM model for the weak form:
2352:   \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi

2354:   Fortran Note:
2355:   This is not yet available in Fortran.

2357: .seealso: `PetscDS`, `PetscDSSetBdJacobianPreconditioner()`
2358: @*/
2359: PetscErrorCode PetscDSGetBdJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2360: {
2361:   PetscBdPointJac *tmp0, *tmp1, *tmp2, *tmp3;
2362:   PetscInt         n0, n1, n2, n3;

2367:   PetscWeakFormGetBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3);
2368:   *g0 = tmp0 ? tmp0[0] : NULL;
2369:   *g1 = tmp1 ? tmp1[0] : NULL;
2370:   *g2 = tmp2 ? tmp2[0] : NULL;
2371:   *g3 = tmp3 ? tmp3[0] : NULL;
2372:   return 0;
2373: }

2375: /*@C
2376:   PetscDSSetBdJacobianPreconditioner - Set the pointwise boundary Jacobian preconditioner function for given test and basis field

2378:   Not collective

2380:   Input Parameters:
2381: + ds - The `PetscDS`
2382: . f  - The test field number
2383: . g  - The field number
2384: . g0 - integrand for the test and basis function term
2385: . g1 - integrand for the test function and basis function gradient term
2386: . g2 - integrand for the test function gradient and basis function term
2387: - g3 - integrand for the test function gradient and basis function gradient term

2389:    Calling sequence for the callbacks g0, g1, g2 and g3:
2390: .vb
2391:   g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2392:      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2393:      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2394:      PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[])
2395: .ve
2396: + dim - the spatial dimension
2397: . Nf - the number of fields
2398: . NfAux - the number of auxiliary fields
2399: . uOff - the offset into u[] and u_t[] for each field
2400: . uOff_x - the offset into u_x[] for each field
2401: . u - each field evaluated at the current point
2402: . u_t - the time derivative of each field evaluated at the current point
2403: . u_x - the gradient of each field evaluated at the current point
2404: . aOff - the offset into a[] and a_t[] for each auxiliary field
2405: . aOff_x - the offset into a_x[] for each auxiliary field
2406: . a - each auxiliary field evaluated at the current point
2407: . a_t - the time derivative of each auxiliary field evaluated at the current point
2408: . a_x - the gradient of auxiliary each field evaluated at the current point
2409: . t - current time
2410: . u_tShift - the multiplier a for dF/dU_t
2411: . x - coordinates of the current point
2412: . n - normal at the current point
2413: . numConstants - number of constant parameters
2414: . constants - constant parameters
2415: - g0 - output values at the current point

2417:   Level: intermediate

2419:   Note:
2420:   We are using a first order FEM model for the weak form:
2421:   \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi

2423:   Fortran Note:
2424:   This is not yet available in Fortran.

2426: .seealso: `PetscDS`, `PetscDSGetBdJacobianPreconditioner()`
2427: @*/
2428: PetscErrorCode PetscDSSetBdJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]), void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]), void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2429: {
2437:   PetscWeakFormSetIndexBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3);
2438:   return 0;
2439: }

2441: /*@C
2442:   PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field

2444:   Not collective

2446:   Input Parameters:
2447: + prob - The PetscDS
2448: - f    - The test field number

2450:   Output Parameters:
2451: + exactSol - exact solution for the test field
2452: - exactCtx - exact solution context

2454:   Calling sequence for the solution functions:
2455: .vb
2456:   sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2457: .ve
2458: + dim - the spatial dimension
2459: . t - current time
2460: . x - coordinates of the current point
2461: . Nc - the number of field components
2462: . u - the solution field evaluated at the current point
2463: - ctx - a user context

2465:   Level: intermediate

2467: .seealso: `PetscDS`, `PetscDSSetExactSolution()`, `PetscDSGetExactSolutionTimeDerivative()`
2468: @*/
2469: PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx)
2470: {
2473:   if (sol) {
2475:     *sol = prob->exactSol[f];
2476:   }
2477:   if (ctx) {
2479:     *ctx = prob->exactCtx[f];
2480:   }
2481:   return 0;
2482: }

2484: /*@C
2485:   PetscDSSetExactSolution - Set the pointwise exact solution function for a given test field

2487:   Not collective

2489:   Input Parameters:
2490: + prob - The `PetscDS`
2491: . f    - The test field number
2492: . sol  - solution function for the test fields
2493: - ctx  - solution context or NULL

2495:   Calling sequence for solution functions:
2496: .vb
2497:   sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2498: .ve
2499: + dim - the spatial dimension
2500: . t - current time
2501: . x - coordinates of the current point
2502: . Nc - the number of field components
2503: . u - the solution field evaluated at the current point
2504: - ctx - a user context

2506:   Level: intermediate

2508: .seealso: `PetscDS`, `PetscDSGetExactSolution()`
2509: @*/
2510: PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx)
2511: {
2514:   PetscDSEnlarge_Static(prob, f + 1);
2515:   if (sol) {
2517:     prob->exactSol[f] = sol;
2518:   }
2519:   if (ctx) {
2521:     prob->exactCtx[f] = ctx;
2522:   }
2523:   return 0;
2524: }

2526: /*@C
2527:   PetscDSGetExactSolutionTimeDerivative - Get the pointwise time derivative of the exact solution function for a given test field

2529:   Not collective

2531:   Input Parameters:
2532: + prob - The `PetscDS`
2533: - f    - The test field number

2535:   Output Parameters:
2536: + exactSol - time derivative of the exact solution for the test field
2537: - exactCtx - time derivative of the exact solution context

2539:   Calling sequence for the solution functions:
2540: .vb
2541:   sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2542: .ve
2543: + dim - the spatial dimension
2544: . t - current time
2545: . x - coordinates of the current point
2546: . Nc - the number of field components
2547: . u - the solution field evaluated at the current point
2548: - ctx - a user context

2550:   Level: intermediate

2552: .seealso: `PetscDS`, `PetscDSSetExactSolutionTimeDerivative()`, `PetscDSGetExactSolution()`
2553: @*/
2554: PetscErrorCode PetscDSGetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx)
2555: {
2558:   if (sol) {
2560:     *sol = prob->exactSol_t[f];
2561:   }
2562:   if (ctx) {
2564:     *ctx = prob->exactCtx_t[f];
2565:   }
2566:   return 0;
2567: }

2569: /*@C
2570:   PetscDSSetExactSolutionTimeDerivative - Set the pointwise time derivative of the exact solution function for a given test field

2572:   Not collective

2574:   Input Parameters:
2575: + prob - The `PetscDS`
2576: . f    - The test field number
2577: . sol  - time derivative of the solution function for the test fields
2578: - ctx  - time derivative of the solution context or NULL

2580:   Calling sequence for solution functions:
2581: .vb
2582:   sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2583: .ve
2584: + dim - the spatial dimension
2585: . t - current time
2586: . x - coordinates of the current point
2587: . Nc - the number of field components
2588: . u - the solution field evaluated at the current point
2589: - ctx - a user context

2591:   Level: intermediate

2593: .seealso: `PetscDS`, `PetscDSGetExactSolutionTimeDerivative()`, `PetscDSSetExactSolution()`
2594: @*/
2595: PetscErrorCode PetscDSSetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx)
2596: {
2599:   PetscDSEnlarge_Static(prob, f + 1);
2600:   if (sol) {
2602:     prob->exactSol_t[f] = sol;
2603:   }
2604:   if (ctx) {
2606:     prob->exactCtx_t[f] = ctx;
2607:   }
2608:   return 0;
2609: }

2611: /*@C
2612:   PetscDSGetConstants - Returns the array of constants passed to point functions

2614:   Not collective

2616:   Input Parameter:
2617: . prob - The `PetscDS` object

2619:   Output Parameters:
2620: + numConstants - The number of constants
2621: - constants    - The array of constants, NULL if there are none

2623:   Level: intermediate

2625: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSCreate()`
2626: @*/
2627: PetscErrorCode PetscDSGetConstants(PetscDS prob, PetscInt *numConstants, const PetscScalar *constants[])
2628: {
2630:   if (numConstants) {
2632:     *numConstants = prob->numConstants;
2633:   }
2634:   if (constants) {
2636:     *constants = prob->constants;
2637:   }
2638:   return 0;
2639: }

2641: /*@C
2642:   PetscDSSetConstants - Set the array of constants passed to point functions

2644:   Not collective

2646:   Input Parameters:
2647: + prob         - The `PetscDS` object
2648: . numConstants - The number of constants
2649: - constants    - The array of constants, NULL if there are none

2651:   Level: intermediate

2653: .seealso: `PetscDS`, `PetscDSGetConstants()`, `PetscDSCreate()`
2654: @*/
2655: PetscErrorCode PetscDSSetConstants(PetscDS prob, PetscInt numConstants, PetscScalar constants[])
2656: {
2658:   if (numConstants != prob->numConstants) {
2659:     PetscFree(prob->constants);
2660:     prob->numConstants = numConstants;
2661:     if (prob->numConstants) {
2662:       PetscMalloc1(prob->numConstants, &prob->constants);
2663:     } else {
2664:       prob->constants = NULL;
2665:     }
2666:   }
2667:   if (prob->numConstants) {
2669:     PetscArraycpy(prob->constants, constants, prob->numConstants);
2670:   }
2671:   return 0;
2672: }

2674: /*@
2675:   PetscDSGetFieldIndex - Returns the index of the given field

2677:   Not collective

2679:   Input Parameters:
2680: + prob - The `PetscDS` object
2681: - disc - The discretization object

2683:   Output Parameter:
2684: . f - The field number

2686:   Level: beginner

2688: .seealso: `PetscDS`, `PetscGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2689: @*/
2690: PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f)
2691: {
2692:   PetscInt g;

2696:   *f = -1;
2697:   for (g = 0; g < prob->Nf; ++g) {
2698:     if (disc == prob->disc[g]) break;
2699:   }
2701:   *f = g;
2702:   return 0;
2703: }

2705: /*@
2706:   PetscDSGetFieldSize - Returns the size of the given field in the full space basis

2708:   Not collective

2710:   Input Parameters:
2711: + prob - The `PetscDS` object
2712: - f - The field number

2714:   Output Parameter:
2715: . size - The size

2717:   Level: beginner

2719: .seealso: `PetscDS`, `PetscDSGetFieldOffset()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2720: @*/
2721: PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size)
2722: {
2726:   PetscDSSetUp(prob);
2727:   *size = prob->Nb[f];
2728:   return 0;
2729: }

2731: /*@
2732:   PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis

2734:   Not collective

2736:   Input Parameters:
2737: + prob - The `PetscDS` object
2738: - f - The field number

2740:   Output Parameter:
2741: . off - The offset

2743:   Level: beginner

2745: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2746: @*/
2747: PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off)
2748: {
2749:   PetscInt size, g;

2754:   *off = 0;
2755:   for (g = 0; g < f; ++g) {
2756:     PetscDSGetFieldSize(prob, g, &size);
2757:     *off += size;
2758:   }
2759:   return 0;
2760: }

2762: /*@
2763:   PetscDSGetFieldOffsetCohesive - Returns the offset of the given field in the full space basis on a cohesive cell

2765:   Not collective

2767:   Input Parameters:
2768: + prob - The `PetscDS` object
2769: - f - The field number

2771:   Output Parameter:
2772: . off - The offset

2774:   Level: beginner

2776: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2777: @*/
2778: PetscErrorCode PetscDSGetFieldOffsetCohesive(PetscDS ds, PetscInt f, PetscInt *off)
2779: {
2780:   PetscInt size, g;

2785:   *off = 0;
2786:   for (g = 0; g < f; ++g) {
2787:     PetscBool cohesive;

2789:     PetscDSGetCohesive(ds, g, &cohesive);
2790:     PetscDSGetFieldSize(ds, g, &size);
2791:     *off += cohesive ? size : size * 2;
2792:   }
2793:   return 0;
2794: }

2796: /*@
2797:   PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point

2799:   Not collective

2801:   Input Parameter:
2802: . prob - The `PetscDS` object

2804:   Output Parameter:
2805: . dimensions - The number of dimensions

2807:   Level: beginner

2809: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2810: @*/
2811: PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[])
2812: {
2814:   PetscDSSetUp(prob);
2816:   *dimensions = prob->Nb;
2817:   return 0;
2818: }

2820: /*@
2821:   PetscDSGetComponents - Returns the number of components for each field on an evaluation point

2823:   Not collective

2825:   Input Parameter:
2826: . prob - The `PetscDS` object

2828:   Output Parameter:
2829: . components - The number of components

2831:   Level: beginner

2833: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2834: @*/
2835: PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[])
2836: {
2838:   PetscDSSetUp(prob);
2840:   *components = prob->Nc;
2841:   return 0;
2842: }

2844: /*@
2845:   PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point

2847:   Not collective

2849:   Input Parameters:
2850: + prob - The `PetscDS` object
2851: - f - The field number

2853:   Output Parameter:
2854: . off - The offset

2856:   Level: beginner

2858: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2859: @*/
2860: PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off)
2861: {
2865:   PetscDSSetUp(prob);
2866:   *off = prob->off[f];
2867:   return 0;
2868: }

2870: /*@
2871:   PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point

2873:   Not collective

2875:   Input Parameter:
2876: . prob - The `PetscDS` object

2878:   Output Parameter:
2879: . offsets - The offsets

2881:   Level: beginner

2883: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2884: @*/
2885: PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[])
2886: {
2889:   PetscDSSetUp(prob);
2890:   *offsets = prob->off;
2891:   return 0;
2892: }

2894: /*@
2895:   PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point

2897:   Not collective

2899:   Input Parameter:
2900: . prob - The `PetscDS` object

2902:   Output Parameter:
2903: . offsets - The offsets

2905:   Level: beginner

2907: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2908: @*/
2909: PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[])
2910: {
2913:   PetscDSSetUp(prob);
2914:   *offsets = prob->offDer;
2915:   return 0;
2916: }

2918: /*@
2919:   PetscDSGetComponentOffsetsCohesive - Returns the offset of each field on an evaluation point

2921:   Not collective

2923:   Input Parameters:
2924: + ds - The `PetscDS` object
2925: - s  - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive

2927:   Output Parameter:
2928: . offsets - The offsets

2930:   Level: beginner

2932: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2933: @*/
2934: PetscErrorCode PetscDSGetComponentOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
2935: {
2940:   PetscDSSetUp(ds);
2941:   *offsets = ds->offCohesive[s];
2942:   return 0;
2943: }

2945: /*@
2946:   PetscDSGetComponentDerivativeOffsetsCohesive - Returns the offset of each field derivative on an evaluation point

2948:   Not collective

2950:   Input Parameters:
2951: + ds - The `PetscDS` object
2952: - s  - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive

2954:   Output Parameter:
2955: . offsets - The offsets

2957:   Level: beginner

2959: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2960: @*/
2961: PetscErrorCode PetscDSGetComponentDerivativeOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
2962: {
2967:   PetscDSSetUp(ds);
2968:   *offsets = ds->offDerCohesive[s];
2969:   return 0;
2970: }

2972: /*@C
2973:   PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization

2975:   Not collective

2977:   Input Parameter:
2978: . prob - The `PetscDS` object

2980:   Output Parameter:
2981: . T - The basis function and derivatives tabulation at quadrature points for each field

2983:   Level: intermediate

2985: .seealso: `PetscDS`, `PetscTabulation`, `PetscDSCreate()`
2986: @*/
2987: PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscTabulation *T[])
2988: {
2991:   PetscDSSetUp(prob);
2992:   *T = prob->T;
2993:   return 0;
2994: }

2996: /*@C
2997:   PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces

2999:   Not collective

3001:   Input Parameter:
3002: . prob - The `PetscDS` object

3004:   Output Parameter:
3005: . Tf - The basis function and derivative tabulation on each local face at quadrature points for each and field

3007:   Level: intermediate

3009: .seealso: `PetscTabulation`, `PetscDS`, `PetscDSGetTabulation()`, `PetscDSCreate()`
3010: @*/
3011: PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscTabulation *Tf[])
3012: {
3015:   PetscDSSetUp(prob);
3016:   *Tf = prob->Tf;
3017:   return 0;
3018: }

3020: PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x)
3021: {
3023:   PetscDSSetUp(prob);
3024:   if (u) {
3026:     *u = prob->u;
3027:   }
3028:   if (u_t) {
3030:     *u_t = prob->u_t;
3031:   }
3032:   if (u_x) {
3034:     *u_x = prob->u_x;
3035:   }
3036:   return 0;
3037: }

3039: PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3)
3040: {
3042:   PetscDSSetUp(prob);
3043:   if (f0) {
3045:     *f0 = prob->f0;
3046:   }
3047:   if (f1) {
3049:     *f1 = prob->f1;
3050:   }
3051:   if (g0) {
3053:     *g0 = prob->g0;
3054:   }
3055:   if (g1) {
3057:     *g1 = prob->g1;
3058:   }
3059:   if (g2) {
3061:     *g2 = prob->g2;
3062:   }
3063:   if (g3) {
3065:     *g3 = prob->g3;
3066:   }
3067:   return 0;
3068: }

3070: PetscErrorCode PetscDSGetWorkspace(PetscDS prob, PetscReal **x, PetscScalar **basisReal, PetscScalar **basisDerReal, PetscScalar **testReal, PetscScalar **testDerReal)
3071: {
3073:   PetscDSSetUp(prob);
3074:   if (x) {
3076:     *x = prob->x;
3077:   }
3078:   if (basisReal) {
3080:     *basisReal = prob->basisReal;
3081:   }
3082:   if (basisDerReal) {
3084:     *basisDerReal = prob->basisDerReal;
3085:   }
3086:   if (testReal) {
3088:     *testReal = prob->testReal;
3089:   }
3090:   if (testDerReal) {
3092:     *testDerReal = prob->testDerReal;
3093:   }
3094:   return 0;
3095: }

3097: /*@C
3098:   PetscDSAddBoundary - Add a boundary condition to the model. The pointwise functions are used to provide boundary values for essential boundary conditions.
3099:   In FEM, they are acting upon by dual basis functionals to generate FEM coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3100:   integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.

3102:   Collective on ds

3104:   Input Parameters:
3105: + ds       - The PetscDS object
3106: . type     - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3107: . name     - The BC name
3108: . label    - The label defining constrained points
3109: . Nv       - The number of `DMLabel` values for constrained points
3110: . values   - An array of label values for constrained points
3111: . field    - The field to constrain
3112: . Nc       - The number of constrained field components (0 will constrain all fields)
3113: . comps    - An array of constrained component numbers
3114: . bcFunc   - A pointwise function giving boundary values
3115: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3116: - ctx      - An optional user context for bcFunc

3118:   Output Parameters:
3119: - bd       - The boundary number

3121:   Options Database Keys:
3122: + -bc_<boundary name> <num> - Overrides the boundary ids
3123: - -bc_<boundary name>_comp <num> - Overrides the boundary components

3125:   Level: developer

3127:   Note:
3128:   Both bcFunc abd bcFunc_t will depend on the boundary condition type. If the type if `DM_BC_ESSENTIAL`, Then the calling sequence is:

3130: $ bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])

3132:   If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value, then the calling sequence is:
3133: .vb
3134:   bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3135:          const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3136:          const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3137:          PetscReal time, const PetscReal x[], PetscScalar bcval[])
3138: .ve
3139: + dim - the spatial dimension
3140: . Nf - the number of fields
3141: . uOff - the offset into u[] and u_t[] for each field
3142: . uOff_x - the offset into u_x[] for each field
3143: . u - each field evaluated at the current point
3144: . u_t - the time derivative of each field evaluated at the current point
3145: . u_x - the gradient of each field evaluated at the current point
3146: . aOff - the offset into a[] and a_t[] for each auxiliary field
3147: . aOff_x - the offset into a_x[] for each auxiliary field
3148: . a - each auxiliary field evaluated at the current point
3149: . a_t - the time derivative of each auxiliary field evaluated at the current point
3150: . a_x - the gradient of auxiliary each field evaluated at the current point
3151: . t - current time
3152: . x - coordinates of the current point
3153: . numConstants - number of constant parameters
3154: . constants - constant parameters
3155: - bcval - output values at the current point

3157: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundaryByName()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3158: @*/
3159: PetscErrorCode PetscDSAddBoundary(PetscDS ds, DMBoundaryConditionType type, const char name[], DMLabel label, PetscInt Nv, const PetscInt values[], PetscInt field, PetscInt Nc, const PetscInt comps[], void (*bcFunc)(void), void (*bcFunc_t)(void), void *ctx, PetscInt *bd)
3160: {
3161:   DSBoundary  head = ds->boundary, b;
3162:   PetscInt    n    = 0;
3163:   const char *lname;

3173:   if (Nc > 0) {
3174:     PetscInt *fcomps;
3175:     PetscInt  c;

3177:     PetscDSGetComponents(ds, &fcomps);
3179:     for (c = 0; c < Nc; ++c) {
3181:     }
3182:   }
3183:   PetscNew(&b);
3184:   PetscStrallocpy(name, (char **)&b->name);
3185:   PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf);
3186:   PetscWeakFormSetNumFields(b->wf, ds->Nf);
3187:   PetscMalloc1(Nv, &b->values);
3188:   if (Nv) PetscArraycpy(b->values, values, Nv);
3189:   PetscMalloc1(Nc, &b->comps);
3190:   if (Nc) PetscArraycpy(b->comps, comps, Nc);
3191:   PetscObjectGetName((PetscObject)label, &lname);
3192:   PetscStrallocpy(lname, (char **)&b->lname);
3193:   b->type   = type;
3194:   b->label  = label;
3195:   b->Nv     = Nv;
3196:   b->field  = field;
3197:   b->Nc     = Nc;
3198:   b->func   = bcFunc;
3199:   b->func_t = bcFunc_t;
3200:   b->ctx    = ctx;
3201:   b->next   = NULL;
3202:   /* Append to linked list so that we can preserve the order */
3203:   if (!head) ds->boundary = b;
3204:   while (head) {
3205:     if (!head->next) {
3206:       head->next = b;
3207:       head       = b;
3208:     }
3209:     head = head->next;
3210:     ++n;
3211:   }
3212:   if (bd) {
3214:     *bd = n;
3215:   }
3216:   return 0;
3217: }

3219: /*@C
3220:   PetscDSAddBoundaryByName - Add a boundary condition to the model. The pointwise functions are used to provide boundary values for essential boundary conditions.
3221:   In FEM, they are acting upon by dual basis functionals to generate FEM coefficients which are fixed. Natural boundary conditions signal to PETSc that
3222:   boundary integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.

3224:   Collective on ds

3226:   Input Parameters:
3227: + ds       - The `PetscDS` object
3228: . type     - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3229: . name     - The BC name
3230: . lname    - The naem of the label defining constrained points
3231: . Nv       - The number of `DMLabel` values for constrained points
3232: . values   - An array of label values for constrained points
3233: . field    - The field to constrain
3234: . Nc       - The number of constrained field components (0 will constrain all fields)
3235: . comps    - An array of constrained component numbers
3236: . bcFunc   - A pointwise function giving boundary values
3237: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3238: - ctx      - An optional user context for bcFunc

3240:   Output Parameters:
3241: - bd       - The boundary number

3243:   Options Database Keys:
3244: + -bc_<boundary name> <num> - Overrides the boundary ids
3245: - -bc_<boundary name>_comp <num> - Overrides the boundary components

3247:   Calling Sequence of bcFunc() and bcFunc_t():
3248:   If the type is `DM_BC_ESSENTIAL`
3249: .vb
3250:   bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])
3251: .ve
3252:   If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value,
3253: .vb
3254:   bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3255:          const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3256:          const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3257:          PetscReal time, const PetscReal x[], PetscScalar bcval[])
3258: .ve
3259: + dim - the spatial dimension
3260: . Nf - the number of fields
3261: . uOff - the offset into u[] and u_t[] for each field
3262: . uOff_x - the offset into u_x[] for each field
3263: . u - each field evaluated at the current point
3264: . u_t - the time derivative of each field evaluated at the current point
3265: . u_x - the gradient of each field evaluated at the current point
3266: . aOff - the offset into a[] and a_t[] for each auxiliary field
3267: . aOff_x - the offset into a_x[] for each auxiliary field
3268: . a - each auxiliary field evaluated at the current point
3269: . a_t - the time derivative of each auxiliary field evaluated at the current point
3270: . a_x - the gradient of auxiliary each field evaluated at the current point
3271: . t - current time
3272: . x - coordinates of the current point
3273: . numConstants - number of constant parameters
3274: . constants - constant parameters
3275: - bcval - output values at the current point

3277:   Level: developer

3279:   Note:
3280:   This function should only be used with `DMFOREST` currently, since labels cannot be defined before the underlying `DMPLEX` is built.

3282: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3283: @*/
3284: PetscErrorCode PetscDSAddBoundaryByName(PetscDS ds, DMBoundaryConditionType type, const char name[], const char lname[], PetscInt Nv, const PetscInt values[], PetscInt field, PetscInt Nc, const PetscInt comps[], void (*bcFunc)(void), void (*bcFunc_t)(void), void *ctx, PetscInt *bd)
3285: {
3286:   DSBoundary head = ds->boundary, b;
3287:   PetscInt   n    = 0;

3296:   PetscNew(&b);
3297:   PetscStrallocpy(name, (char **)&b->name);
3298:   PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf);
3299:   PetscWeakFormSetNumFields(b->wf, ds->Nf);
3300:   PetscMalloc1(Nv, &b->values);
3301:   if (Nv) PetscArraycpy(b->values, values, Nv);
3302:   PetscMalloc1(Nc, &b->comps);
3303:   if (Nc) PetscArraycpy(b->comps, comps, Nc);
3304:   PetscStrallocpy(lname, (char **)&b->lname);
3305:   b->type   = type;
3306:   b->label  = NULL;
3307:   b->Nv     = Nv;
3308:   b->field  = field;
3309:   b->Nc     = Nc;
3310:   b->func   = bcFunc;
3311:   b->func_t = bcFunc_t;
3312:   b->ctx    = ctx;
3313:   b->next   = NULL;
3314:   /* Append to linked list so that we can preserve the order */
3315:   if (!head) ds->boundary = b;
3316:   while (head) {
3317:     if (!head->next) {
3318:       head->next = b;
3319:       head       = b;
3320:     }
3321:     head = head->next;
3322:     ++n;
3323:   }
3324:   if (bd) {
3326:     *bd = n;
3327:   }
3328:   return 0;
3329: }

3331: /*@C
3332:   PetscDSUpdateBoundary - Change a boundary condition for the model. The pointwise functions are used to provide boundary values for essential boundary conditions.
3333:   In FEM, they are acting upon by dual basis functionals to generate FEM coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary integrals
3334:   should be performed, using the kernels from `PetscDSSetBdResidual()`.

3336:   Input Parameters:
3337: + ds       - The `PetscDS` object
3338: . bd       - The boundary condition number
3339: . type     - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3340: . name     - The BC name
3341: . label    - The label defining constrained points
3342: . Nv       - The number of `DMLabel` ids for constrained points
3343: . values   - An array of ids for constrained points
3344: . field    - The field to constrain
3345: . Nc       - The number of constrained field components
3346: . comps    - An array of constrained component numbers
3347: . bcFunc   - A pointwise function giving boundary values
3348: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3349: - ctx      - An optional user context for bcFunc

3351:   Level: developer

3353:   Note:
3354:   The boundary condition number is the order in which it was registered. The user can get the number of boundary conditions from `PetscDSGetNumBoundary()`.
3355:   See `PetscDSAddBoundary()` for a description of the calling sequences for the callbacks.

3357: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSGetNumBoundary()`, `DMLabel`
3358: @*/
3359: PetscErrorCode PetscDSUpdateBoundary(PetscDS ds, PetscInt bd, DMBoundaryConditionType type, const char name[], DMLabel label, PetscInt Nv, const PetscInt values[], PetscInt field, PetscInt Nc, const PetscInt comps[], void (*bcFunc)(void), void (*bcFunc_t)(void), void *ctx)
3360: {
3361:   DSBoundary b = ds->boundary;
3362:   PetscInt   n = 0;

3365:   while (b) {
3366:     if (n == bd) break;
3367:     b = b->next;
3368:     ++n;
3369:   }
3371:   if (name) {
3372:     PetscFree(b->name);
3373:     PetscStrallocpy(name, (char **)&b->name);
3374:   }
3375:   b->type = type;
3376:   if (label) {
3377:     const char *name;

3379:     b->label = label;
3380:     PetscFree(b->lname);
3381:     PetscObjectGetName((PetscObject)label, &name);
3382:     PetscStrallocpy(name, (char **)&b->lname);
3383:   }
3384:   if (Nv >= 0) {
3385:     b->Nv = Nv;
3386:     PetscFree(b->values);
3387:     PetscMalloc1(Nv, &b->values);
3388:     if (Nv) PetscArraycpy(b->values, values, Nv);
3389:   }
3390:   if (field >= 0) b->field = field;
3391:   if (Nc >= 0) {
3392:     b->Nc = Nc;
3393:     PetscFree(b->comps);
3394:     PetscMalloc1(Nc, &b->comps);
3395:     if (Nc) PetscArraycpy(b->comps, comps, Nc);
3396:   }
3397:   if (bcFunc) b->func = bcFunc;
3398:   if (bcFunc_t) b->func_t = bcFunc_t;
3399:   if (ctx) b->ctx = ctx;
3400:   return 0;
3401: }

3403: /*@
3404:   PetscDSGetNumBoundary - Get the number of registered BC

3406:   Input Parameters:
3407: . ds - The `PetscDS` object

3409:   Output Parameters:
3410: . numBd - The number of BC

3412:   Level: intermediate

3414: .seealso: `PetscDS`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`
3415: @*/
3416: PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd)
3417: {
3418:   DSBoundary b = ds->boundary;

3422:   *numBd = 0;
3423:   while (b) {
3424:     ++(*numBd);
3425:     b = b->next;
3426:   }
3427:   return 0;
3428: }

3430: /*@C
3431:   PetscDSGetBoundary - Gets a boundary condition to the model

3433:   Input Parameters:
3434: + ds          - The `PetscDS` object
3435: - bd          - The BC number

3437:   Output Parameters:
3438: + wf       - The `PetscWeakForm` holding the pointwise functions
3439: . type     - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3440: . name     - The BC name
3441: . label    - The label defining constrained points
3442: . Nv       - The number of `DMLabel` ids for constrained points
3443: . values   - An array of ids for constrained points
3444: . field    - The field to constrain
3445: . Nc       - The number of constrained field components
3446: . comps    - An array of constrained component numbers
3447: . bcFunc   - A pointwise function giving boundary values
3448: . bcFunc_t - A pointwise function giving the time derivative of the boundary values
3449: - ctx      - An optional user context for bcFunc

3451:   Options Database Keys:
3452: + -bc_<boundary name> <num> - Overrides the boundary ids
3453: - -bc_<boundary name>_comp <num> - Overrides the boundary components

3455:   Level: developer

3457: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `DMLabel`
3458: @*/
3459: PetscErrorCode PetscDSGetBoundary(PetscDS ds, PetscInt bd, PetscWeakForm *wf, DMBoundaryConditionType *type, const char *name[], DMLabel *label, PetscInt *Nv, const PetscInt *values[], PetscInt *field, PetscInt *Nc, const PetscInt *comps[], void (**func)(void), void (**func_t)(void), void **ctx)
3460: {
3461:   DSBoundary b = ds->boundary;
3462:   PetscInt   n = 0;

3465:   while (b) {
3466:     if (n == bd) break;
3467:     b = b->next;
3468:     ++n;
3469:   }
3471:   if (wf) {
3473:     *wf = b->wf;
3474:   }
3475:   if (type) {
3477:     *type = b->type;
3478:   }
3479:   if (name) {
3481:     *name = b->name;
3482:   }
3483:   if (label) {
3485:     *label = b->label;
3486:   }
3487:   if (Nv) {
3489:     *Nv = b->Nv;
3490:   }
3491:   if (values) {
3493:     *values = b->values;
3494:   }
3495:   if (field) {
3497:     *field = b->field;
3498:   }
3499:   if (Nc) {
3501:     *Nc = b->Nc;
3502:   }
3503:   if (comps) {
3505:     *comps = b->comps;
3506:   }
3507:   if (func) {
3509:     *func = b->func;
3510:   }
3511:   if (func_t) {
3513:     *func_t = b->func_t;
3514:   }
3515:   if (ctx) {
3517:     *ctx = b->ctx;
3518:   }
3519:   return 0;
3520: }

3522: static PetscErrorCode DSBoundaryDuplicate_Internal(DSBoundary b, DSBoundary *bNew)
3523: {
3524:   PetscNew(bNew);
3525:   PetscWeakFormCreate(PETSC_COMM_SELF, &(*bNew)->wf);
3526:   PetscWeakFormCopy(b->wf, (*bNew)->wf);
3527:   PetscStrallocpy(b->name, (char **)&((*bNew)->name));
3528:   PetscStrallocpy(b->lname, (char **)&((*bNew)->lname));
3529:   (*bNew)->type  = b->type;
3530:   (*bNew)->label = b->label;
3531:   (*bNew)->Nv    = b->Nv;
3532:   PetscMalloc1(b->Nv, &(*bNew)->values);
3533:   PetscArraycpy((*bNew)->values, b->values, b->Nv);
3534:   (*bNew)->field = b->field;
3535:   (*bNew)->Nc    = b->Nc;
3536:   PetscMalloc1(b->Nc, &(*bNew)->comps);
3537:   PetscArraycpy((*bNew)->comps, b->comps, b->Nc);
3538:   (*bNew)->func   = b->func;
3539:   (*bNew)->func_t = b->func_t;
3540:   (*bNew)->ctx    = b->ctx;
3541:   return 0;
3542: }

3544: /*@
3545:   PetscDSCopyBoundary - Copy all boundary condition objects to the new problem

3547:   Not collective

3549:   Input Parameters:
3550: + ds        - The source `PetscDS` object
3551: . numFields - The number of selected fields, or `PETSC_DEFAULT` for all fields
3552: - fields    - The selected fields, or NULL for all fields

3554:   Output Parameter:
3555: . newds     - The target `PetscDS`, now with a copy of the boundary conditions

3557:   Level: intermediate

3559: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3560: @*/
3561: PetscErrorCode PetscDSCopyBoundary(PetscDS ds, PetscInt numFields, const PetscInt fields[], PetscDS newds)
3562: {
3563:   DSBoundary b, *lastnext;

3567:   if (ds == newds) return 0;
3568:   PetscDSDestroyBoundary(newds);
3569:   lastnext = &(newds->boundary);
3570:   for (b = ds->boundary; b; b = b->next) {
3571:     DSBoundary bNew;
3572:     PetscInt   fieldNew = -1;

3574:     if (numFields > 0 && fields) {
3575:       PetscInt f;

3577:       for (f = 0; f < numFields; ++f)
3578:         if (b->field == fields[f]) break;
3579:       if (f == numFields) continue;
3580:       fieldNew = f;
3581:     }
3582:     DSBoundaryDuplicate_Internal(b, &bNew);
3583:     bNew->field = fieldNew < 0 ? b->field : fieldNew;
3584:     *lastnext   = bNew;
3585:     lastnext    = &(bNew->next);
3586:   }
3587:   return 0;
3588: }

3590: /*@
3591:   PetscDSDestroyBoundary - Remove all `DMBoundary` objects from the `PetscDS`

3593:   Not collective

3595:   Input Parameter:
3596: . ds - The `PetscDS` object

3598:   Level: intermediate

3600: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`
3601: @*/
3602: PetscErrorCode PetscDSDestroyBoundary(PetscDS ds)
3603: {
3604:   DSBoundary next = ds->boundary;

3606:   while (next) {
3607:     DSBoundary b = next;

3609:     next = b->next;
3610:     PetscWeakFormDestroy(&b->wf);
3611:     PetscFree(b->name);
3612:     PetscFree(b->lname);
3613:     PetscFree(b->values);
3614:     PetscFree(b->comps);
3615:     PetscFree(b);
3616:   }
3617:   return 0;
3618: }

3620: /*@
3621:   PetscDSSelectDiscretizations - Copy discretizations to the new problem with different field layout

3623:   Not collective

3625:   Input Parameters:
3626: + prob - The `PetscDS` object
3627: . numFields - Number of new fields
3628: - fields - Old field number for each new field

3630:   Output Parameter:
3631: . newprob - The `PetscDS` copy

3633:   Level: intermediate

3635: .seealso: `PetscDS`, `PetscDSSelectEquations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3636: @*/
3637: PetscErrorCode PetscDSSelectDiscretizations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
3638: {
3639:   PetscInt Nf, Nfn, fn;

3644:   PetscDSGetNumFields(prob, &Nf);
3645:   PetscDSGetNumFields(newprob, &Nfn);
3646:   numFields = numFields < 0 ? Nf : numFields;
3647:   for (fn = 0; fn < numFields; ++fn) {
3648:     const PetscInt f = fields ? fields[fn] : fn;
3649:     PetscObject    disc;

3651:     if (f >= Nf) continue;
3652:     PetscDSGetDiscretization(prob, f, &disc);
3653:     PetscDSSetDiscretization(newprob, fn, disc);
3654:   }
3655:   return 0;
3656: }

3658: /*@
3659:   PetscDSSelectEquations - Copy pointwise function pointers to the new problem with different field layout

3661:   Not collective

3663:   Input Parameters:
3664: + prob - The `PetscDS` object
3665: . numFields - Number of new fields
3666: - fields - Old field number for each new field

3668:   Output Parameter:
3669: . newprob - The `PetscDS` copy

3671:   Level: intermediate

3673: .seealso: `PetscDS`, `PetscDSSelectDiscretizations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3674: @*/
3675: PetscErrorCode PetscDSSelectEquations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
3676: {
3677:   PetscInt Nf, Nfn, fn, gn;

3682:   PetscDSGetNumFields(prob, &Nf);
3683:   PetscDSGetNumFields(newprob, &Nfn);
3685:   for (fn = 0; fn < numFields; ++fn) {
3686:     const PetscInt   f = fields ? fields[fn] : fn;
3687:     PetscPointFunc   obj;
3688:     PetscPointFunc   f0, f1;
3689:     PetscBdPointFunc f0Bd, f1Bd;
3690:     PetscRiemannFunc r;

3692:     if (f >= Nf) continue;
3693:     PetscDSGetObjective(prob, f, &obj);
3694:     PetscDSGetResidual(prob, f, &f0, &f1);
3695:     PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd);
3696:     PetscDSGetRiemannSolver(prob, f, &r);
3697:     PetscDSSetObjective(newprob, fn, obj);
3698:     PetscDSSetResidual(newprob, fn, f0, f1);
3699:     PetscDSSetBdResidual(newprob, fn, f0Bd, f1Bd);
3700:     PetscDSSetRiemannSolver(newprob, fn, r);
3701:     for (gn = 0; gn < numFields; ++gn) {
3702:       const PetscInt  g = fields ? fields[gn] : gn;
3703:       PetscPointJac   g0, g1, g2, g3;
3704:       PetscPointJac   g0p, g1p, g2p, g3p;
3705:       PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd;

3707:       if (g >= Nf) continue;
3708:       PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3);
3709:       PetscDSGetJacobianPreconditioner(prob, f, g, &g0p, &g1p, &g2p, &g3p);
3710:       PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd);
3711:       PetscDSSetJacobian(newprob, fn, gn, g0, g1, g2, g3);
3712:       PetscDSSetJacobianPreconditioner(newprob, fn, gn, g0p, g1p, g2p, g3p);
3713:       PetscDSSetBdJacobian(newprob, fn, gn, g0Bd, g1Bd, g2Bd, g3Bd);
3714:     }
3715:   }
3716:   return 0;
3717: }

3719: /*@
3720:   PetscDSCopyEquations - Copy all pointwise function pointers to another `PetscDS`

3722:   Not collective

3724:   Input Parameter:
3725: . prob - The `PetscDS` object

3727:   Output Parameter:
3728: . newprob - The `PetscDS` copy

3730:   Level: intermediate

3732: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3733: @*/
3734: PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob)
3735: {
3736:   PetscWeakForm wf, newwf;
3737:   PetscInt      Nf, Ng;

3741:   PetscDSGetNumFields(prob, &Nf);
3742:   PetscDSGetNumFields(newprob, &Ng);
3744:   PetscDSGetWeakForm(prob, &wf);
3745:   PetscDSGetWeakForm(newprob, &newwf);
3746:   PetscWeakFormCopy(wf, newwf);
3747:   return 0;
3748: }

3750: /*@
3751:   PetscDSCopyConstants - Copy all constants to another `PetscDS`

3753:   Not collective

3755:   Input Parameter:
3756: . prob - The `PetscDS` object

3758:   Output Parameter:
3759: . newprob - The `PetscDS` copy

3761:   Level: intermediate

3763: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3764: @*/
3765: PetscErrorCode PetscDSCopyConstants(PetscDS prob, PetscDS newprob)
3766: {
3767:   PetscInt           Nc;
3768:   const PetscScalar *constants;

3772:   PetscDSGetConstants(prob, &Nc, &constants);
3773:   PetscDSSetConstants(newprob, Nc, (PetscScalar *)constants);
3774:   return 0;
3775: }

3777: /*@
3778:   PetscDSCopyExactSolutions - Copy all exact solutions to another `PetscDS`

3780:   Not collective

3782:   Input Parameter:
3783: . ds - The `PetscDS` object

3785:   Output Parameter:
3786: . newds - The `PetscDS` copy

3788:   Level: intermediate

3790: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3791: @*/
3792: PetscErrorCode PetscDSCopyExactSolutions(PetscDS ds, PetscDS newds)
3793: {
3794:   PetscSimplePointFunc sol;
3795:   void                *ctx;
3796:   PetscInt             Nf, f;

3800:   PetscDSGetNumFields(ds, &Nf);
3801:   for (f = 0; f < Nf; ++f) {
3802:     PetscDSGetExactSolution(ds, f, &sol, &ctx);
3803:     PetscDSSetExactSolution(newds, f, sol, ctx);
3804:     PetscDSGetExactSolutionTimeDerivative(ds, f, &sol, &ctx);
3805:     PetscDSSetExactSolutionTimeDerivative(newds, f, sol, ctx);
3806:   }
3807:   return 0;
3808: }

3810: PetscErrorCode PetscDSGetHeightSubspace(PetscDS prob, PetscInt height, PetscDS *subprob)
3811: {
3812:   PetscInt dim, Nf, f;

3816:   if (height == 0) {
3817:     *subprob = prob;
3818:     return 0;
3819:   }
3820:   PetscDSGetNumFields(prob, &Nf);
3821:   PetscDSGetSpatialDimension(prob, &dim);
3823:   if (!prob->subprobs) PetscCalloc1(dim, &prob->subprobs);
3824:   if (!prob->subprobs[height - 1]) {
3825:     PetscInt cdim;

3827:     PetscDSCreate(PetscObjectComm((PetscObject)prob), &prob->subprobs[height - 1]);
3828:     PetscDSGetCoordinateDimension(prob, &cdim);
3829:     PetscDSSetCoordinateDimension(prob->subprobs[height - 1], cdim);
3830:     for (f = 0; f < Nf; ++f) {
3831:       PetscFE      subfe;
3832:       PetscObject  obj;
3833:       PetscClassId id;

3835:       PetscDSGetDiscretization(prob, f, &obj);
3836:       PetscObjectGetClassId(obj, &id);
3837:       if (id == PETSCFE_CLASSID) PetscFEGetHeightSubspace((PetscFE)obj, height, &subfe);
3838:       else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unsupported discretization type for field %" PetscInt_FMT, f);
3839:       PetscDSSetDiscretization(prob->subprobs[height - 1], f, (PetscObject)subfe);
3840:     }
3841:   }
3842:   *subprob = prob->subprobs[height - 1];
3843:   return 0;
3844: }

3846: PetscErrorCode PetscDSGetDiscType_Internal(PetscDS ds, PetscInt f, PetscDiscType *disctype)
3847: {
3848:   PetscObject  obj;
3849:   PetscClassId id;
3850:   PetscInt     Nf;

3854:   *disctype = PETSC_DISC_NONE;
3855:   PetscDSGetNumFields(ds, &Nf);
3857:   PetscDSGetDiscretization(ds, f, &obj);
3858:   if (obj) {
3859:     PetscObjectGetClassId(obj, &id);
3860:     if (id == PETSCFE_CLASSID) *disctype = PETSC_DISC_FE;
3861:     else *disctype = PETSC_DISC_FV;
3862:   }
3863:   return 0;
3864: }

3866: static PetscErrorCode PetscDSDestroy_Basic(PetscDS ds)
3867: {
3868:   PetscFree(ds->data);
3869:   return 0;
3870: }

3872: static PetscErrorCode PetscDSInitialize_Basic(PetscDS ds)
3873: {
3874:   ds->ops->setfromoptions = NULL;
3875:   ds->ops->setup          = NULL;
3876:   ds->ops->view           = NULL;
3877:   ds->ops->destroy        = PetscDSDestroy_Basic;
3878:   return 0;
3879: }

3881: /*MC
3882:   PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions

3884:   Level: intermediate

3886: .seealso: `PetscDSType`, `PetscDSCreate()`, `PetscDSSetType()`
3887: M*/

3889: PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS ds)
3890: {
3891:   PetscDS_Basic *b;

3894:   PetscNew(&b);
3895:   ds->data = b;

3897:   PetscDSInitialize_Basic(ds);
3898:   return 0;
3899: }