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; No Fortran Support

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

 30:   Example 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: .seealso: `PetscDSType`, `PetscDS`, `PetscDSRegisterAll()`, `PetscDSRegisterDestroy()`
 51: @*/
 52: PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS))
 53: {
 54:   PetscFunctionBegin;
 55:   PetscCall(PetscFunctionListAdd(&PetscDSList, sname, function));
 56:   PetscFunctionReturn(PETSC_SUCCESS);
 57: }

 59: /*@C
 60:   PetscDSSetType - Builds a particular `PetscDS`

 62:   Collective; No Fortran Support

 64:   Input Parameters:
 65: + prob - The `PetscDS` object
 66: - name - The `PetscDSType`

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

 71:   Level: intermediate

 73: .seealso: `PetscDSType`, `PetscDS`, `PetscDSGetType()`, `PetscDSCreate()`
 74: @*/
 75: PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name)
 76: {
 77:   PetscErrorCode (*r)(PetscDS);
 78:   PetscBool match;

 80:   PetscFunctionBegin;
 82:   PetscCall(PetscObjectTypeCompare((PetscObject)prob, name, &match));
 83:   if (match) PetscFunctionReturn(PETSC_SUCCESS);

 85:   PetscCall(PetscDSRegisterAll());
 86:   PetscCall(PetscFunctionListFind(PetscDSList, name, &r));
 87:   PetscCheck(r, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscDS type: %s", name);

 89:   PetscTryTypeMethod(prob, destroy);
 90:   prob->ops->destroy = NULL;

 92:   PetscCall((*r)(prob));
 93:   PetscCall(PetscObjectChangeTypeName((PetscObject)prob, name));
 94:   PetscFunctionReturn(PETSC_SUCCESS);
 95: }

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

100:   Not Collective; No Fortran Support

102:   Input Parameter:
103: . prob - The `PetscDS`

105:   Output Parameter:
106: . name - The `PetscDSType` name

108:   Level: intermediate

110: .seealso: `PetscDSType`, `PetscDS`, `PetscDSSetType()`, `PetscDSCreate()`
111: @*/
112: PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name)
113: {
114:   PetscFunctionBegin;
116:   PetscAssertPointer(name, 2);
117:   PetscCall(PetscDSRegisterAll());
118:   *name = ((PetscObject)prob)->type_name;
119:   PetscFunctionReturn(PETSC_SUCCESS);
120: }

122: static PetscErrorCode PetscDSView_Ascii(PetscDS ds, PetscViewer viewer)
123: {
124:   PetscViewerFormat  format;
125:   const PetscScalar *constants;
126:   PetscInt           Nf, numConstants, f;

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

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

173:     for (b = ds->boundary; b; b = b->next) {
174:       char    *name;
175:       PetscInt c, i;

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

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

237:   Collective

239:   Input Parameters:
240: + A    - the `PetscDS` object
241: . obj  - Optional object that provides the options prefix used in the search
242: - name - command line option

244:   Level: intermediate

246: .seealso: `PetscDSType`, `PetscDS`, `PetscDSView()`, `PetscObjectViewFromOptions()`, `PetscDSCreate()`
247: @*/
248: PetscErrorCode PetscDSViewFromOptions(PetscDS A, PetscObject obj, const char name[])
249: {
250:   PetscFunctionBegin;
252:   PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
253:   PetscFunctionReturn(PETSC_SUCCESS);
254: }

256: /*@C
257:   PetscDSView - Views a `PetscDS`

259:   Collective

261:   Input Parameters:
262: + prob - the `PetscDS` object to view
263: - v    - the viewer

265:   Level: developer

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

273:   PetscFunctionBegin;
275:   if (!v) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)prob), &v));
277:   PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &iascii));
278:   if (iascii) PetscCall(PetscDSView_Ascii(prob, v));
279:   PetscTryTypeMethod(prob, view, v);
280:   PetscFunctionReturn(PETSC_SUCCESS);
281: }

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

286:   Collective

288:   Input Parameter:
289: . prob - the `PetscDS` object to set options for

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

298:   Level: intermediate

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

309:   PetscFunctionBegin;
311:   if (!((PetscObject)prob)->type_name) {
312:     defaultType = PETSCDSBASIC;
313:   } else {
314:     defaultType = ((PetscObject)prob)->type_name;
315:   }
316:   PetscCall(PetscDSRegisterAll());

318:   PetscObjectOptionsBegin((PetscObject)prob);
319:   for (b = prob->boundary; b; b = b->next) {
320:     char      optname[1024];
321:     PetscInt  ids[1024], len = 1024;
322:     PetscBool flg;

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

361: /*@C
362:   PetscDSSetUp - Construct data structures for the `PetscDS`

364:   Collective

366:   Input Parameter:
367: . prob - the `PetscDS` object to setup

369:   Level: developer

371: .seealso: `PetscDS`, `PetscDSView()`, `PetscDSDestroy()`
372: @*/
373: PetscErrorCode PetscDSSetUp(PetscDS prob)
374: {
375:   const PetscInt Nf          = prob->Nf;
376:   PetscBool      hasH        = PETSC_FALSE;
377:   PetscInt       maxOrder[4] = {-1, -1, -1, -1};
378:   PetscInt       dim, dimEmbed, NbMax = 0, NcMax = 0, NqMax = 0, NsMax = 1, f;

380:   PetscFunctionBegin;
382:   if (prob->setup) PetscFunctionReturn(PETSC_SUCCESS);
383:   /* Calculate sizes */
384:   PetscCall(PetscDSGetSpatialDimension(prob, &dim));
385:   PetscCall(PetscDSGetCoordinateDimension(prob, &dimEmbed));
386:   prob->totDim = prob->totComp = 0;
387:   PetscCall(PetscMalloc2(Nf, &prob->Nc, Nf, &prob->Nb));
388:   PetscCall(PetscCalloc2(Nf + 1, &prob->off, Nf + 1, &prob->offDer));
389:   PetscCall(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]));
390:   PetscCall(PetscMalloc2(Nf, &prob->T, Nf, &prob->Tf));
391:   if (prob->forceQuad) {
392:     // Note: This assumes we have one kind of cell at each dimension.
393:     //       We can fix this by having quadrature hold the celltype
394:     PetscQuadrature maxQuad[4] = {NULL, NULL, NULL, NULL};

396:     for (f = 0; f < Nf; ++f) {
397:       PetscObject     obj;
398:       PetscClassId    id;
399:       PetscQuadrature q = NULL, fq = NULL;
400:       PetscInt        dim = -1, order = -1, forder = -1;

402:       PetscCall(PetscDSGetDiscretization(prob, f, &obj));
403:       if (!obj) continue;
404:       PetscCall(PetscObjectGetClassId(obj, &id));
405:       if (id == PETSCFE_CLASSID) {
406:         PetscFE fe = (PetscFE)obj;

408:         PetscCall(PetscFEGetQuadrature(fe, &q));
409:         PetscCall(PetscFEGetFaceQuadrature(fe, &fq));
410:       } else if (id == PETSCFV_CLASSID) {
411:         PetscFV fv = (PetscFV)obj;

413:         PetscCall(PetscFVGetQuadrature(fv, &q));
414:       }
415:       if (q) {
416:         PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
417:         PetscCall(PetscQuadratureGetOrder(q, &order));
418:         if (order > maxOrder[dim]) {
419:           maxOrder[dim] = order;
420:           maxQuad[dim]  = q;
421:         }
422:       }
423:       if (fq) {
424:         PetscCall(PetscQuadratureGetData(fq, &dim, NULL, NULL, NULL, NULL));
425:         PetscCall(PetscQuadratureGetOrder(fq, &forder));
426:         if (forder > maxOrder[dim]) {
427:           maxOrder[dim] = forder;
428:           maxQuad[dim]  = fq;
429:         }
430:       }
431:     }
432:     for (f = 0; f < Nf; ++f) {
433:       PetscObject     obj;
434:       PetscClassId    id;
435:       PetscQuadrature q;
436:       PetscInt        dim;

438:       PetscCall(PetscDSGetDiscretization(prob, f, &obj));
439:       if (!obj) continue;
440:       PetscCall(PetscObjectGetClassId(obj, &id));
441:       if (id == PETSCFE_CLASSID) {
442:         PetscFE fe = (PetscFE)obj;

444:         PetscCall(PetscFEGetQuadrature(fe, &q));
445:         PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
446:         PetscCall(PetscFESetQuadrature(fe, maxQuad[dim]));
447:         PetscCall(PetscFESetFaceQuadrature(fe, dim ? maxQuad[dim - 1] : NULL));
448:       } else if (id == PETSCFV_CLASSID) {
449:         PetscFV fv = (PetscFV)obj;

451:         PetscCall(PetscFVGetQuadrature(fv, &q));
452:         PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
453:         PetscCall(PetscFVSetQuadrature(fv, maxQuad[dim]));
454:       }
455:     }
456:   }
457:   for (f = 0; f < Nf; ++f) {
458:     PetscObject     obj;
459:     PetscClassId    id;
460:     PetscQuadrature q  = NULL;
461:     PetscInt        Nq = 0, Nb, Nc;

463:     PetscCall(PetscDSGetDiscretization(prob, f, &obj));
464:     if (prob->jetDegree[f] > 1) hasH = PETSC_TRUE;
465:     if (!obj) {
466:       /* Empty mesh */
467:       Nb = Nc    = 0;
468:       prob->T[f] = prob->Tf[f] = NULL;
469:     } else {
470:       PetscCall(PetscObjectGetClassId(obj, &id));
471:       if (id == PETSCFE_CLASSID) {
472:         PetscFE fe = (PetscFE)obj;

474:         PetscCall(PetscFEGetQuadrature(fe, &q));
475:         {
476:           PetscQuadrature fq;
477:           PetscInt        dim, order;

479:           PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
480:           PetscCall(PetscQuadratureGetOrder(q, &order));
481:           if (maxOrder[dim] < 0) maxOrder[dim] = order;
482:           PetscCheck(order == maxOrder[dim], PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Field %" PetscInt_FMT " cell quadrature order %" PetscInt_FMT " != %" PetscInt_FMT " DS cell quadrature order", f, order, maxOrder[dim]);
483:           PetscCall(PetscFEGetFaceQuadrature(fe, &fq));
484:           if (fq) {
485:             PetscCall(PetscQuadratureGetData(fq, &dim, NULL, NULL, NULL, NULL));
486:             PetscCall(PetscQuadratureGetOrder(fq, &order));
487:             if (maxOrder[dim] < 0) maxOrder[dim] = order;
488:             PetscCheck(order == maxOrder[dim], PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Field %" PetscInt_FMT " face quadrature order %" PetscInt_FMT " != %" PetscInt_FMT " DS face quadrature order", f, order, maxOrder[dim]);
489:           }
490:         }
491:         PetscCall(PetscFEGetDimension(fe, &Nb));
492:         PetscCall(PetscFEGetNumComponents(fe, &Nc));
493:         PetscCall(PetscFEGetCellTabulation(fe, prob->jetDegree[f], &prob->T[f]));
494:         PetscCall(PetscFEGetFaceTabulation(fe, prob->jetDegree[f], &prob->Tf[f]));
495:       } else if (id == PETSCFV_CLASSID) {
496:         PetscFV fv = (PetscFV)obj;

498:         PetscCall(PetscFVGetQuadrature(fv, &q));
499:         PetscCall(PetscFVGetNumComponents(fv, &Nc));
500:         Nb = Nc;
501:         PetscCall(PetscFVGetCellTabulation(fv, &prob->T[f]));
502:         /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */
503:       } else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f);
504:     }
505:     prob->Nc[f]                    = Nc;
506:     prob->Nb[f]                    = Nb;
507:     prob->off[f + 1]               = Nc + prob->off[f];
508:     prob->offDer[f + 1]            = Nc * dim + prob->offDer[f];
509:     prob->offCohesive[0][f + 1]    = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[0][f];
510:     prob->offDerCohesive[0][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[0][f];
511:     prob->offCohesive[1][f]        = (prob->cohesive[f] ? 0 : Nc) + prob->offCohesive[0][f];
512:     prob->offDerCohesive[1][f]     = (prob->cohesive[f] ? 0 : Nc) * dimEmbed + prob->offDerCohesive[0][f];
513:     prob->offCohesive[2][f + 1]    = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[2][f];
514:     prob->offDerCohesive[2][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[2][f];
515:     if (q) PetscCall(PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL));
516:     NqMax = PetscMax(NqMax, Nq);
517:     NbMax = PetscMax(NbMax, Nb);
518:     NcMax = PetscMax(NcMax, Nc);
519:     prob->totDim += Nb;
520:     prob->totComp += Nc;
521:     /* There are two faces for all fields on a cohesive cell, except for cohesive fields */
522:     if (prob->isCohesive && !prob->cohesive[f]) prob->totDim += Nb;
523:   }
524:   prob->offCohesive[1][Nf]    = prob->offCohesive[0][Nf];
525:   prob->offDerCohesive[1][Nf] = prob->offDerCohesive[0][Nf];
526:   /* Allocate works space */
527:   NsMax = 2; /* A non-cohesive discretizations can be used on a cohesive cell, so we need this extra workspace for all DS */
528:   PetscCall(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));
529:   PetscCall(PetscMalloc5(dimEmbed, &prob->x, NbMax * NcMax, &prob->basisReal, NbMax * NcMax * dimEmbed, &prob->basisDerReal, NbMax * NcMax, &prob->testReal, NbMax * NcMax * dimEmbed, &prob->testDerReal));
530:   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,
531:                          &prob->g2, NsMax * NsMax * NqMax * NcMax * NcMax * dimEmbed * dimEmbed, &prob->g3));
532:   PetscTryTypeMethod(prob, setup);
533:   prob->setup = PETSC_TRUE;
534:   PetscFunctionReturn(PETSC_SUCCESS);
535: }

537: static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob)
538: {
539:   PetscFunctionBegin;
540:   PetscCall(PetscFree2(prob->Nc, prob->Nb));
541:   PetscCall(PetscFree2(prob->off, prob->offDer));
542:   PetscCall(PetscFree6(prob->offCohesive[0], prob->offCohesive[1], prob->offCohesive[2], prob->offDerCohesive[0], prob->offDerCohesive[1], prob->offDerCohesive[2]));
543:   PetscCall(PetscFree2(prob->T, prob->Tf));
544:   PetscCall(PetscFree3(prob->u, prob->u_t, prob->u_x));
545:   PetscCall(PetscFree5(prob->x, prob->basisReal, prob->basisDerReal, prob->testReal, prob->testDerReal));
546:   PetscCall(PetscFree6(prob->f0, prob->f1, prob->g0, prob->g1, prob->g2, prob->g3));
547:   PetscFunctionReturn(PETSC_SUCCESS);
548: }

550: static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew)
551: {
552:   PetscObject          *tmpd;
553:   PetscBool            *tmpi;
554:   PetscInt             *tmpk;
555:   PetscBool            *tmpc;
556:   PetscPointFunc       *tmpup;
557:   PetscSimplePointFunc *tmpexactSol, *tmpexactSol_t;
558:   void                **tmpexactCtx, **tmpexactCtx_t;
559:   void                **tmpctx;
560:   PetscInt              Nf = prob->Nf, f;

562:   PetscFunctionBegin;
563:   if (Nf >= NfNew) PetscFunctionReturn(PETSC_SUCCESS);
564:   prob->setup = PETSC_FALSE;
565:   PetscCall(PetscDSDestroyStructs_Static(prob));
566:   PetscCall(PetscMalloc4(NfNew, &tmpd, NfNew, &tmpi, NfNew, &tmpc, NfNew, &tmpk));
567:   for (f = 0; f < Nf; ++f) {
568:     tmpd[f] = prob->disc[f];
569:     tmpi[f] = prob->implicit[f];
570:     tmpc[f] = prob->cohesive[f];
571:     tmpk[f] = prob->jetDegree[f];
572:   }
573:   for (f = Nf; f < NfNew; ++f) {
574:     tmpd[f] = NULL;
575:     tmpi[f] = PETSC_TRUE, tmpc[f] = PETSC_FALSE;
576:     tmpk[f] = 1;
577:   }
578:   PetscCall(PetscFree4(prob->disc, prob->implicit, prob->cohesive, prob->jetDegree));
579:   PetscCall(PetscWeakFormSetNumFields(prob->wf, NfNew));
580:   prob->Nf        = NfNew;
581:   prob->disc      = tmpd;
582:   prob->implicit  = tmpi;
583:   prob->cohesive  = tmpc;
584:   prob->jetDegree = tmpk;
585:   PetscCall(PetscCalloc2(NfNew, &tmpup, NfNew, &tmpctx));
586:   for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f];
587:   for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f];
588:   for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL;
589:   for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL;
590:   PetscCall(PetscFree2(prob->update, prob->ctx));
591:   prob->update = tmpup;
592:   prob->ctx    = tmpctx;
593:   PetscCall(PetscCalloc4(NfNew, &tmpexactSol, NfNew, &tmpexactCtx, NfNew, &tmpexactSol_t, NfNew, &tmpexactCtx_t));
594:   for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f];
595:   for (f = 0; f < Nf; ++f) tmpexactCtx[f] = prob->exactCtx[f];
596:   for (f = 0; f < Nf; ++f) tmpexactSol_t[f] = prob->exactSol_t[f];
597:   for (f = 0; f < Nf; ++f) tmpexactCtx_t[f] = prob->exactCtx_t[f];
598:   for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL;
599:   for (f = Nf; f < NfNew; ++f) tmpexactCtx[f] = NULL;
600:   for (f = Nf; f < NfNew; ++f) tmpexactSol_t[f] = NULL;
601:   for (f = Nf; f < NfNew; ++f) tmpexactCtx_t[f] = NULL;
602:   PetscCall(PetscFree4(prob->exactSol, prob->exactCtx, prob->exactSol_t, prob->exactCtx_t));
603:   prob->exactSol   = tmpexactSol;
604:   prob->exactCtx   = tmpexactCtx;
605:   prob->exactSol_t = tmpexactSol_t;
606:   prob->exactCtx_t = tmpexactCtx_t;
607:   PetscFunctionReturn(PETSC_SUCCESS);
608: }

610: /*@
611:   PetscDSDestroy - Destroys a `PetscDS` object

613:   Collective

615:   Input Parameter:
616: . ds - the `PetscDS` object to destroy

618:   Level: developer

620: .seealso: `PetscDSView()`
621: @*/
622: PetscErrorCode PetscDSDestroy(PetscDS *ds)
623: {
624:   PetscInt f;

626:   PetscFunctionBegin;
627:   if (!*ds) PetscFunctionReturn(PETSC_SUCCESS);

630:   if (--((PetscObject)(*ds))->refct > 0) {
631:     *ds = NULL;
632:     PetscFunctionReturn(PETSC_SUCCESS);
633:   }
634:   ((PetscObject)(*ds))->refct = 0;
635:   if ((*ds)->subprobs) {
636:     PetscInt dim, d;

638:     PetscCall(PetscDSGetSpatialDimension(*ds, &dim));
639:     for (d = 0; d < dim; ++d) PetscCall(PetscDSDestroy(&(*ds)->subprobs[d]));
640:   }
641:   PetscCall(PetscFree((*ds)->subprobs));
642:   PetscCall(PetscDSDestroyStructs_Static(*ds));
643:   for (f = 0; f < (*ds)->Nf; ++f) PetscCall(PetscObjectDereference((*ds)->disc[f]));
644:   PetscCall(PetscFree4((*ds)->disc, (*ds)->implicit, (*ds)->cohesive, (*ds)->jetDegree));
645:   PetscCall(PetscWeakFormDestroy(&(*ds)->wf));
646:   PetscCall(PetscFree2((*ds)->update, (*ds)->ctx));
647:   PetscCall(PetscFree4((*ds)->exactSol, (*ds)->exactCtx, (*ds)->exactSol_t, (*ds)->exactCtx_t));
648:   PetscTryTypeMethod((*ds), destroy);
649:   PetscCall(PetscDSDestroyBoundary(*ds));
650:   PetscCall(PetscFree((*ds)->constants));
651:   for (PetscInt c = 0; c < DM_NUM_POLYTOPES; ++c) {
652:     const PetscInt Na = DMPolytopeTypeGetNumArrangments((DMPolytopeType)c);
653:     if ((*ds)->quadPerm[c])
654:       for (PetscInt o = 0; o < Na; ++o) PetscCall(ISDestroy(&(*ds)->quadPerm[c][o]));
655:     PetscCall(PetscFree((*ds)->quadPerm[c]));
656:     (*ds)->quadPerm[c] = NULL;
657:   }
658:   PetscCall(PetscHeaderDestroy(ds));
659:   PetscFunctionReturn(PETSC_SUCCESS);
660: }

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

665:   Collective

667:   Input Parameter:
668: . comm - The communicator for the `PetscDS` object

670:   Output Parameter:
671: . ds - The `PetscDS` object

673:   Level: beginner

675: .seealso: `PetscDS`, `PetscDSSetType()`, `PETSCDSBASIC`, `PetscDSType`
676: @*/
677: PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *ds)
678: {
679:   PetscDS p;

681:   PetscFunctionBegin;
682:   PetscAssertPointer(ds, 2);
683:   *ds = NULL;
684:   PetscCall(PetscDSInitializePackage());

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

688:   p->Nf           = 0;
689:   p->setup        = PETSC_FALSE;
690:   p->numConstants = 0;
691:   p->constants    = NULL;
692:   p->dimEmbed     = -1;
693:   p->useJacPre    = PETSC_TRUE;
694:   p->forceQuad    = PETSC_TRUE;
695:   PetscCall(PetscWeakFormCreate(comm, &p->wf));
696:   PetscCall(PetscArrayzero(p->quadPerm, DM_NUM_POLYTOPES));

698:   *ds = p;
699:   PetscFunctionReturn(PETSC_SUCCESS);
700: }

702: /*@
703:   PetscDSGetNumFields - Returns the number of fields in the `PetscDS`

705:   Not Collective

707:   Input Parameter:
708: . prob - The `PetscDS` object

710:   Output Parameter:
711: . Nf - The number of fields

713:   Level: beginner

715: .seealso: `PetscDS`, `PetscDSGetSpatialDimension()`, `PetscDSCreate()`
716: @*/
717: PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf)
718: {
719:   PetscFunctionBegin;
721:   PetscAssertPointer(Nf, 2);
722:   *Nf = prob->Nf;
723:   PetscFunctionReturn(PETSC_SUCCESS);
724: }

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

729:   Not Collective

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

734:   Output Parameter:
735: . dim - The spatial dimension

737:   Level: beginner

739: .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
740: @*/
741: PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim)
742: {
743:   PetscFunctionBegin;
745:   PetscAssertPointer(dim, 2);
746:   *dim = 0;
747:   if (prob->Nf) {
748:     PetscObject  obj;
749:     PetscClassId id;

751:     PetscCall(PetscDSGetDiscretization(prob, 0, &obj));
752:     if (obj) {
753:       PetscCall(PetscObjectGetClassId(obj, &id));
754:       if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetSpatialDimension((PetscFE)obj, dim));
755:       else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetSpatialDimension((PetscFV)obj, dim));
756:       else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
757:     }
758:   }
759:   PetscFunctionReturn(PETSC_SUCCESS);
760: }

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

765:   Not Collective

767:   Input Parameter:
768: . prob - The `PetscDS` object

770:   Output Parameter:
771: . dimEmbed - The coordinate dimension

773:   Level: beginner

775: .seealso: `PetscDS`, `PetscDSSetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
776: @*/
777: PetscErrorCode PetscDSGetCoordinateDimension(PetscDS prob, PetscInt *dimEmbed)
778: {
779:   PetscFunctionBegin;
781:   PetscAssertPointer(dimEmbed, 2);
782:   PetscCheck(prob->dimEmbed >= 0, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONGSTATE, "No coordinate dimension set for this DS");
783:   *dimEmbed = prob->dimEmbed;
784:   PetscFunctionReturn(PETSC_SUCCESS);
785: }

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

790:   Logically Collective

792:   Input Parameters:
793: + prob     - The `PetscDS` object
794: - dimEmbed - The coordinate dimension

796:   Level: beginner

798: .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
799: @*/
800: PetscErrorCode PetscDSSetCoordinateDimension(PetscDS prob, PetscInt dimEmbed)
801: {
802:   PetscFunctionBegin;
804:   PetscCheck(dimEmbed >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Coordinate dimension must be non-negative, not %" PetscInt_FMT, dimEmbed);
805:   prob->dimEmbed = dimEmbed;
806:   PetscFunctionReturn(PETSC_SUCCESS);
807: }

809: /*@
810:   PetscDSGetForceQuad - Returns the flag to force matching quadratures among the field discretizations

812:   Not collective

814:   Input Parameter:
815: . ds - The `PetscDS` object

817:   Output Parameter:
818: . forceQuad - The flag

820:   Level: intermediate

822: .seealso: `PetscDS`, `PetscDSSetForceQuad()`, `PetscDSGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
823: @*/
824: PetscErrorCode PetscDSGetForceQuad(PetscDS ds, PetscBool *forceQuad)
825: {
826:   PetscFunctionBegin;
828:   PetscAssertPointer(forceQuad, 2);
829:   *forceQuad = ds->forceQuad;
830:   PetscFunctionReturn(PETSC_SUCCESS);
831: }

833: /*@
834:   PetscDSSetForceQuad - Set the flag to force matching quadratures among the field discretizations

836:   Logically collective on ds

838:   Input Parameters:
839: + ds        - The `PetscDS` object
840: - forceQuad - The flag

842:   Level: intermediate

844: .seealso: `PetscDS`, `PetscDSGetForceQuad()`, `PetscDSGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
845: @*/
846: PetscErrorCode PetscDSSetForceQuad(PetscDS ds, PetscBool forceQuad)
847: {
848:   PetscFunctionBegin;
850:   ds->forceQuad = forceQuad;
851:   PetscFunctionReturn(PETSC_SUCCESS);
852: }

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

857:   Not Collective

859:   Input Parameter:
860: . ds - The `PetscDS` object

862:   Output Parameter:
863: . isCohesive - The flag

865:   Level: developer

867: .seealso: `PetscDS`, `PetscDSGetNumCohesive()`, `PetscDSGetCohesive()`, `PetscDSSetCohesive()`, `PetscDSCreate()`
868: @*/
869: PetscErrorCode PetscDSIsCohesive(PetscDS ds, PetscBool *isCohesive)
870: {
871:   PetscFunctionBegin;
873:   PetscAssertPointer(isCohesive, 2);
874:   *isCohesive = ds->isCohesive;
875:   PetscFunctionReturn(PETSC_SUCCESS);
876: }

878: /*@
879:   PetscDSGetNumCohesive - Returns the number of cohesive fields, meaning those defined on the interior of a cohesive cell

881:   Not Collective

883:   Input Parameter:
884: . ds - The `PetscDS` object

886:   Output Parameter:
887: . numCohesive - The number of cohesive fields

889:   Level: developer

891: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSCreate()`
892: @*/
893: PetscErrorCode PetscDSGetNumCohesive(PetscDS ds, PetscInt *numCohesive)
894: {
895:   PetscInt f;

897:   PetscFunctionBegin;
899:   PetscAssertPointer(numCohesive, 2);
900:   *numCohesive = 0;
901:   for (f = 0; f < ds->Nf; ++f) *numCohesive += ds->cohesive[f] ? 1 : 0;
902:   PetscFunctionReturn(PETSC_SUCCESS);
903: }

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

908:   Not Collective

910:   Input Parameters:
911: + ds - The `PetscDS` object
912: - f  - The field index

914:   Output Parameter:
915: . isCohesive - The flag

917:   Level: developer

919: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
920: @*/
921: PetscErrorCode PetscDSGetCohesive(PetscDS ds, PetscInt f, PetscBool *isCohesive)
922: {
923:   PetscFunctionBegin;
925:   PetscAssertPointer(isCohesive, 3);
926:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
927:   *isCohesive = ds->cohesive[f];
928:   PetscFunctionReturn(PETSC_SUCCESS);
929: }

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

934:   Not Collective

936:   Input Parameters:
937: + ds         - The `PetscDS` object
938: . f          - The field index
939: - isCohesive - The flag for a cohesive field

941:   Level: developer

943: .seealso: `PetscDS`, `PetscDSGetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
944: @*/
945: PetscErrorCode PetscDSSetCohesive(PetscDS ds, PetscInt f, PetscBool isCohesive)
946: {
947:   PetscInt i;

949:   PetscFunctionBegin;
951:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
952:   ds->cohesive[f] = isCohesive;
953:   ds->isCohesive  = PETSC_FALSE;
954:   for (i = 0; i < ds->Nf; ++i) ds->isCohesive = ds->isCohesive || ds->cohesive[f] ? PETSC_TRUE : PETSC_FALSE;
955:   PetscFunctionReturn(PETSC_SUCCESS);
956: }

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

961:   Not Collective

963:   Input Parameter:
964: . prob - The `PetscDS` object

966:   Output Parameter:
967: . dim - The total problem dimension

969:   Level: beginner

971: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
972: @*/
973: PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim)
974: {
975:   PetscFunctionBegin;
977:   PetscCall(PetscDSSetUp(prob));
978:   PetscAssertPointer(dim, 2);
979:   *dim = prob->totDim;
980:   PetscFunctionReturn(PETSC_SUCCESS);
981: }

983: /*@
984:   PetscDSGetTotalComponents - Returns the total number of components in this system

986:   Not Collective

988:   Input Parameter:
989: . prob - The `PetscDS` object

991:   Output Parameter:
992: . Nc - The total number of components

994:   Level: beginner

996: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
997: @*/
998: PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc)
999: {
1000:   PetscFunctionBegin;
1002:   PetscCall(PetscDSSetUp(prob));
1003:   PetscAssertPointer(Nc, 2);
1004:   *Nc = prob->totComp;
1005:   PetscFunctionReturn(PETSC_SUCCESS);
1006: }

1008: /*@
1009:   PetscDSGetDiscretization - Returns the discretization object for the given field

1011:   Not Collective

1013:   Input Parameters:
1014: + prob - The `PetscDS` object
1015: - f    - The field number

1017:   Output Parameter:
1018: . disc - The discretization object

1020:   Level: beginner

1022: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1023: @*/
1024: PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc)
1025: {
1026:   PetscFunctionBeginHot;
1028:   PetscAssertPointer(disc, 3);
1029:   PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
1030:   *disc = prob->disc[f];
1031:   PetscFunctionReturn(PETSC_SUCCESS);
1032: }

1034: /*@
1035:   PetscDSSetDiscretization - Sets the discretization object for the given field

1037:   Not Collective

1039:   Input Parameters:
1040: + prob - The `PetscDS` object
1041: . f    - The field number
1042: - disc - The discretization object

1044:   Level: beginner

1046: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSGetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1047: @*/
1048: PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc)
1049: {
1050:   PetscFunctionBegin;
1052:   if (disc) PetscAssertPointer(disc, 3);
1053:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1054:   PetscCall(PetscDSEnlarge_Static(prob, f + 1));
1055:   PetscCall(PetscObjectDereference(prob->disc[f]));
1056:   prob->disc[f] = disc;
1057:   PetscCall(PetscObjectReference(disc));
1058:   if (disc) {
1059:     PetscClassId id;

1061:     PetscCall(PetscObjectGetClassId(disc, &id));
1062:     if (id == PETSCFE_CLASSID) {
1063:       PetscCall(PetscDSSetImplicit(prob, f, PETSC_TRUE));
1064:     } else if (id == PETSCFV_CLASSID) {
1065:       PetscCall(PetscDSSetImplicit(prob, f, PETSC_FALSE));
1066:     }
1067:     PetscCall(PetscDSSetJetDegree(prob, f, 1));
1068:   }
1069:   PetscFunctionReturn(PETSC_SUCCESS);
1070: }

1072: /*@
1073:   PetscDSGetWeakForm - Returns the weak form object

1075:   Not Collective

1077:   Input Parameter:
1078: . ds - The `PetscDS` object

1080:   Output Parameter:
1081: . wf - The weak form object

1083:   Level: beginner

1085: .seealso: `PetscWeakForm`, `PetscDSSetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1086: @*/
1087: PetscErrorCode PetscDSGetWeakForm(PetscDS ds, PetscWeakForm *wf)
1088: {
1089:   PetscFunctionBegin;
1091:   PetscAssertPointer(wf, 2);
1092:   *wf = ds->wf;
1093:   PetscFunctionReturn(PETSC_SUCCESS);
1094: }

1096: /*@
1097:   PetscDSSetWeakForm - Sets the weak form object

1099:   Not Collective

1101:   Input Parameters:
1102: + ds - The `PetscDS` object
1103: - wf - The weak form object

1105:   Level: beginner

1107: .seealso: `PetscWeakForm`, `PetscDSGetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1108: @*/
1109: PetscErrorCode PetscDSSetWeakForm(PetscDS ds, PetscWeakForm wf)
1110: {
1111:   PetscFunctionBegin;
1114:   PetscCall(PetscObjectDereference((PetscObject)ds->wf));
1115:   ds->wf = wf;
1116:   PetscCall(PetscObjectReference((PetscObject)wf));
1117:   PetscCall(PetscWeakFormSetNumFields(wf, ds->Nf));
1118:   PetscFunctionReturn(PETSC_SUCCESS);
1119: }

1121: /*@
1122:   PetscDSAddDiscretization - Adds a discretization object

1124:   Not Collective

1126:   Input Parameters:
1127: + prob - The `PetscDS` object
1128: - disc - The boundary discretization object

1130:   Level: beginner

1132: .seealso: `PetscWeakForm`, `PetscDSGetDiscretization()`, `PetscDSSetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1133: @*/
1134: PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc)
1135: {
1136:   PetscFunctionBegin;
1137:   PetscCall(PetscDSSetDiscretization(prob, prob->Nf, disc));
1138:   PetscFunctionReturn(PETSC_SUCCESS);
1139: }

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

1144:   Not Collective

1146:   Input Parameter:
1147: . prob - The `PetscDS` object

1149:   Output Parameter:
1150: . q - The quadrature object

1152:   Level: intermediate

1154: .seealso: `PetscDS`, `PetscQuadrature`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1155: @*/
1156: PetscErrorCode PetscDSGetQuadrature(PetscDS prob, PetscQuadrature *q)
1157: {
1158:   PetscObject  obj;
1159:   PetscClassId id;

1161:   PetscFunctionBegin;
1162:   *q = NULL;
1163:   if (!prob->Nf) PetscFunctionReturn(PETSC_SUCCESS);
1164:   PetscCall(PetscDSGetDiscretization(prob, 0, &obj));
1165:   PetscCall(PetscObjectGetClassId(obj, &id));
1166:   if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetQuadrature((PetscFE)obj, q));
1167:   else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetQuadrature((PetscFV)obj, q));
1168:   else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
1169:   PetscFunctionReturn(PETSC_SUCCESS);
1170: }

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

1175:   Not Collective

1177:   Input Parameters:
1178: + prob - The `PetscDS` object
1179: - f    - The field number

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

1184:   Level: developer

1186: .seealso: `TSIMEX`, `PetscDS`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1187: @*/
1188: PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit)
1189: {
1190:   PetscFunctionBegin;
1192:   PetscAssertPointer(implicit, 3);
1193:   PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
1194:   *implicit = prob->implicit[f];
1195:   PetscFunctionReturn(PETSC_SUCCESS);
1196: }

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

1201:   Not Collective

1203:   Input Parameters:
1204: + prob     - The `PetscDS` object
1205: . f        - The field number
1206: - implicit - The flag indicating what kind of solve to use for this field

1208:   Level: developer

1210: .seealso: `TSIMEX`, `PetscDSGetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1211: @*/
1212: PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit)
1213: {
1214:   PetscFunctionBegin;
1216:   PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
1217:   prob->implicit[f] = implicit;
1218:   PetscFunctionReturn(PETSC_SUCCESS);
1219: }

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

1224:   Not Collective

1226:   Input Parameters:
1227: + ds - The `PetscDS` object
1228: - f  - The field number

1230:   Output Parameter:
1231: . k - The highest derivative we need to tabulate

1233:   Level: developer

1235: .seealso: `PetscDS`, `PetscDSSetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1236: @*/
1237: PetscErrorCode PetscDSGetJetDegree(PetscDS ds, PetscInt f, PetscInt *k)
1238: {
1239:   PetscFunctionBegin;
1241:   PetscAssertPointer(k, 3);
1242:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1243:   *k = ds->jetDegree[f];
1244:   PetscFunctionReturn(PETSC_SUCCESS);
1245: }

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

1250:   Not Collective

1252:   Input Parameters:
1253: + ds - The `PetscDS` object
1254: . f  - The field number
1255: - k  - The highest derivative we need to tabulate

1257:   Level: developer

1259: .seealso: ``PetscDS`, `PetscDSGetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1260: @*/
1261: PetscErrorCode PetscDSSetJetDegree(PetscDS ds, PetscInt f, PetscInt k)
1262: {
1263:   PetscFunctionBegin;
1265:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1266:   ds->jetDegree[f] = k;
1267:   PetscFunctionReturn(PETSC_SUCCESS);
1268: }

1270: /*@C
1271:   PetscDSGetObjective - Get the pointwise objective function for a given test field

1273:   Not Collective

1275:   Input Parameters:
1276: + ds - The `PetscDS`
1277: - f  - The test field number

1279:   Output Parameter:
1280: . obj - integrand for the test function term

1282:   Calling sequence of `obj`:
1283: + dim          - the spatial dimension
1284: . Nf           - the number of fields
1285: . NfAux        - the number of auxiliary fields
1286: . uOff         - the offset into u[] and u_t[] for each field
1287: . uOff_x       - the offset into u_x[] for each field
1288: . u            - each field evaluated at the current point
1289: . u_t          - the time derivative of each field evaluated at the current point
1290: . u_x          - the gradient of each field evaluated at the current point
1291: . aOff         - the offset into a[] and a_t[] for each auxiliary field
1292: . aOff_x       - the offset into a_x[] for each auxiliary field
1293: . a            - each auxiliary field evaluated at the current point
1294: . a_t          - the time derivative of each auxiliary field evaluated at the current point
1295: . a_x          - the gradient of auxiliary each field evaluated at the current point
1296: . t            - current time
1297: . x            - coordinates of the current point
1298: . numConstants - number of constant parameters
1299: . constants    - constant parameters
1300: - obj          - output values at the current point

1302:   Level: intermediate

1304:   Note:
1305:   We are using a first order FEM model for the weak form:  \int_\Omega \phi obj(u, u_t, \nabla u, x, t)

1307: .seealso: `PetscDS`, `PetscDSSetObjective()`, `PetscDSGetResidual()`
1308: @*/
1309: 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[]))
1310: {
1311:   PetscPointFunc *tmp;
1312:   PetscInt        n;

1314:   PetscFunctionBegin;
1316:   PetscAssertPointer(obj, 3);
1317:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1318:   PetscCall(PetscWeakFormGetObjective(ds->wf, NULL, 0, f, 0, &n, &tmp));
1319:   *obj = tmp ? tmp[0] : NULL;
1320:   PetscFunctionReturn(PETSC_SUCCESS);
1321: }

1323: /*@C
1324:   PetscDSSetObjective - Set the pointwise objective function for a given test field

1326:   Not Collective

1328:   Input Parameters:
1329: + ds  - The `PetscDS`
1330: . f   - The test field number
1331: - obj - integrand for the test function term

1333:   Calling sequence of `obj`:
1334: + dim          - the spatial dimension
1335: . Nf           - the number of fields
1336: . NfAux        - the number of auxiliary fields
1337: . uOff         - the offset into u[] and u_t[] for each field
1338: . uOff_x       - the offset into u_x[] for each field
1339: . u            - each field evaluated at the current point
1340: . u_t          - the time derivative of each field evaluated at the current point
1341: . u_x          - the gradient of each field evaluated at the current point
1342: . aOff         - the offset into a[] and a_t[] for each auxiliary field
1343: . aOff_x       - the offset into a_x[] for each auxiliary field
1344: . a            - each auxiliary field evaluated at the current point
1345: . a_t          - the time derivative of each auxiliary field evaluated at the current point
1346: . a_x          - the gradient of auxiliary each field evaluated at the current point
1347: . t            - current time
1348: . x            - coordinates of the current point
1349: . numConstants - number of constant parameters
1350: . constants    - constant parameters
1351: - obj          - output values at the current point

1353:   Level: intermediate

1355:   Note:
1356:   We are using a first order FEM model for the weak form:  \int_\Omega \phi obj(u, u_t, \nabla u, x, t)

1358: .seealso: `PetscDS`, `PetscDSGetObjective()`, `PetscDSSetResidual()`
1359: @*/
1360: 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[]))
1361: {
1362:   PetscFunctionBegin;
1365:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1366:   PetscCall(PetscWeakFormSetIndexObjective(ds->wf, NULL, 0, f, 0, 0, obj));
1367:   PetscFunctionReturn(PETSC_SUCCESS);
1368: }

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

1373:   Not Collective

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

1379:   Output Parameters:
1380: + f0 - integrand for the test function term
1381: - f1 - integrand for the test function gradient term

1383:   Calling sequence of `f0`:
1384: + dim          - the spatial dimension
1385: . Nf           - the number of fields
1386: . NfAux        - the number of auxiliary fields
1387: . uOff         - the offset into u[] and u_t[] for each field
1388: . uOff_x       - the offset into u_x[] for each field
1389: . u            - each field evaluated at the current point
1390: . u_t          - the time derivative of each field evaluated at the current point
1391: . u_x          - the gradient of each field evaluated at the current point
1392: . aOff         - the offset into a[] and a_t[] for each auxiliary field
1393: . aOff_x       - the offset into a_x[] for each auxiliary field
1394: . a            - each auxiliary field evaluated at the current point
1395: . a_t          - the time derivative of each auxiliary field evaluated at the current point
1396: . a_x          - the gradient of auxiliary each field evaluated at the current point
1397: . t            - current time
1398: . x            - coordinates of the current point
1399: . numConstants - number of constant parameters
1400: . constants    - constant parameters
1401: - f0           - output values at the current point

1403:   Level: intermediate

1405:   Note:
1406:   `f1` has an identical form and is omitted for brevity.

1408:   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)

1410: .seealso: `PetscDS`, `PetscDSSetResidual()`
1411: @*/
1412: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1413: {
1414:   PetscPointFunc *tmp0, *tmp1;
1415:   PetscInt        n0, n1;

1417:   PetscFunctionBegin;
1419:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1420:   PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1));
1421:   *f0 = tmp0 ? tmp0[0] : NULL;
1422:   *f1 = tmp1 ? tmp1[0] : NULL;
1423:   PetscFunctionReturn(PETSC_SUCCESS);
1424: }

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

1429:   Not Collective

1431:   Input Parameters:
1432: + ds - The `PetscDS`
1433: . f  - The test field number
1434: . f0 - integrand for the test function term
1435: - f1 - integrand for the test function gradient term

1437:   Calling sequence of `f0`:
1438: + dim          - the spatial dimension
1439: . Nf           - the number of fields
1440: . NfAux        - the number of auxiliary fields
1441: . uOff         - the offset into u[] and u_t[] for each field
1442: . uOff_x       - the offset into u_x[] for each field
1443: . u            - each field evaluated at the current point
1444: . u_t          - the time derivative of each field evaluated at the current point
1445: . u_x          - the gradient of each field evaluated at the current point
1446: . aOff         - the offset into a[] and a_t[] for each auxiliary field
1447: . aOff_x       - the offset into a_x[] for each auxiliary field
1448: . a            - each auxiliary field evaluated at the current point
1449: . a_t          - the time derivative of each auxiliary field evaluated at the current point
1450: . a_x          - the gradient of auxiliary each field evaluated at the current point
1451: . t            - current time
1452: . x            - coordinates of the current point
1453: . numConstants - number of constant parameters
1454: . constants    - constant parameters
1455: - f0           - output values at the current point

1457:   Level: intermediate

1459:   Note:
1460:   `f1` has an identical form and is omitted for brevity.

1462:   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)

1464: .seealso: `PetscDS`, `PetscDSGetResidual()`
1465: @*/
1466: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1467: {
1468:   PetscFunctionBegin;
1472:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1473:   PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1));
1474:   PetscFunctionReturn(PETSC_SUCCESS);
1475: }

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

1480:   Not Collective

1482:   Input Parameters:
1483: + ds - The `PetscDS`
1484: - f  - The test field number

1486:   Output Parameters:
1487: + f0 - integrand for the test function term
1488: - f1 - integrand for the test function gradient term

1490:   Calling sequence of `f0`:
1491: + dim          - the spatial dimension
1492: . Nf           - the number of fields
1493: . NfAux        - the number of auxiliary fields
1494: . uOff         - the offset into u[] and u_t[] for each field
1495: . uOff_x       - the offset into u_x[] for each field
1496: . u            - each field evaluated at the current point
1497: . u_t          - the time derivative of each field evaluated at the current point
1498: . u_x          - the gradient of each field evaluated at the current point
1499: . aOff         - the offset into a[] and a_t[] for each auxiliary field
1500: . aOff_x       - the offset into a_x[] for each auxiliary field
1501: . a            - each auxiliary field evaluated at the current point
1502: . a_t          - the time derivative of each auxiliary field evaluated at the current point
1503: . a_x          - the gradient of auxiliary each field evaluated at the current point
1504: . t            - current time
1505: . x            - coordinates of the current point
1506: . numConstants - number of constant parameters
1507: . constants    - constant parameters
1508: - f0           - output values at the current point

1510:   Level: intermediate

1512:   Note:
1513:   `f1` has an identical form and is omitted for brevity.

1515:   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)

1517: .seealso: `PetscDS`, `PetscDSSetRHSResidual()`
1518: @*/
1519: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1520: {
1521:   PetscPointFunc *tmp0, *tmp1;
1522:   PetscInt        n0, n1;

1524:   PetscFunctionBegin;
1526:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1527:   PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 100, &n0, &tmp0, &n1, &tmp1));
1528:   *f0 = tmp0 ? tmp0[0] : NULL;
1529:   *f1 = tmp1 ? tmp1[0] : NULL;
1530:   PetscFunctionReturn(PETSC_SUCCESS);
1531: }

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

1536:   Not Collective

1538:   Input Parameters:
1539: + ds - The `PetscDS`
1540: . f  - The test field number
1541: . f0 - integrand for the test function term
1542: - f1 - integrand for the test function gradient term

1544:   Calling sequence for the callbacks `f0`:
1545: + dim          - the spatial dimension
1546: . Nf           - the number of fields
1547: . NfAux        - the number of auxiliary fields
1548: . uOff         - the offset into u[] and u_t[] for each field
1549: . uOff_x       - the offset into u_x[] for each field
1550: . u            - each field evaluated at the current point
1551: . u_t          - the time derivative of each field evaluated at the current point
1552: . u_x          - the gradient of each field evaluated at the current point
1553: . aOff         - the offset into a[] and a_t[] for each auxiliary field
1554: . aOff_x       - the offset into a_x[] for each auxiliary field
1555: . a            - each auxiliary field evaluated at the current point
1556: . a_t          - the time derivative of each auxiliary field evaluated at the current point
1557: . a_x          - the gradient of auxiliary each field evaluated at the current point
1558: . t            - current time
1559: . x            - coordinates of the current point
1560: . numConstants - number of constant parameters
1561: . constants    - constant parameters
1562: - f0           - output values at the current point

1564:   Level: intermediate

1566:   Note:
1567:   `f1` has an identical form and is omitted for brevity.

1569:   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)

1571: .seealso: `PetscDS`, `PetscDSGetResidual()`
1572: @*/
1573: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1574: {
1575:   PetscFunctionBegin;
1579:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1580:   PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 100, 0, f0, 0, f1));
1581:   PetscFunctionReturn(PETSC_SUCCESS);
1582: }

1584: /*@C
1585:   PetscDSHasJacobian - Checks that the Jacobian functions have been set

1587:   Not Collective

1589:   Input Parameter:
1590: . ds - The `PetscDS`

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

1595:   Level: intermediate

1597: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1598: @*/
1599: PetscErrorCode PetscDSHasJacobian(PetscDS ds, PetscBool *hasJac)
1600: {
1601:   PetscFunctionBegin;
1603:   PetscCall(PetscWeakFormHasJacobian(ds->wf, hasJac));
1604:   PetscFunctionReturn(PETSC_SUCCESS);
1605: }

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

1610:   Not Collective

1612:   Input Parameters:
1613: + ds - The `PetscDS`
1614: . f  - The test field number
1615: - g  - The field number

1617:   Output Parameters:
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 of `g0`:
1624: + dim          - the spatial dimension
1625: . Nf           - the number of fields
1626: . NfAux        - the number of auxiliary fields
1627: . uOff         - the offset into u[] and u_t[] for each field
1628: . uOff_x       - the offset into u_x[] for each field
1629: . u            - each field evaluated at the current point
1630: . u_t          - the time derivative of each field evaluated at the current point
1631: . u_x          - the gradient of each field evaluated at the current point
1632: . aOff         - the offset into a[] and a_t[] for each auxiliary field
1633: . aOff_x       - the offset into a_x[] for each auxiliary field
1634: . a            - each auxiliary field evaluated at the current point
1635: . a_t          - the time derivative of each auxiliary field evaluated at the current point
1636: . a_x          - the gradient of auxiliary each field evaluated at the current point
1637: . t            - current time
1638: . u_tShift     - the multiplier a for dF/dU_t
1639: . x            - coordinates of the current point
1640: . numConstants - number of constant parameters
1641: . constants    - constant parameters
1642: - g0           - output values at the current point

1644:   Level: intermediate

1646:   Note:
1647:   `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.

1649:   We are using a first order FEM model for the weak form\:
1650:   \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

1652: .seealso: `PetscDS`, `PetscDSSetJacobian()`
1653: @*/
1654: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1655: {
1656:   PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1657:   PetscInt       n0, n1, n2, n3;

1659:   PetscFunctionBegin;
1661:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1662:   PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf);
1663:   PetscCall(PetscWeakFormGetJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1664:   *g0 = tmp0 ? tmp0[0] : NULL;
1665:   *g1 = tmp1 ? tmp1[0] : NULL;
1666:   *g2 = tmp2 ? tmp2[0] : NULL;
1667:   *g3 = tmp3 ? tmp3[0] : NULL;
1668:   PetscFunctionReturn(PETSC_SUCCESS);
1669: }

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

1674:   Not Collective

1676:   Input Parameters:
1677: + ds - The `PetscDS`
1678: . f  - The test field number
1679: . g  - The field number
1680: . g0 - integrand for the test and basis function term
1681: . g1 - integrand for the test function and basis function gradient term
1682: . g2 - integrand for the test function gradient and basis function term
1683: - g3 - integrand for the test function gradient and basis function gradient term

1685:   Calling sequence of `g0`:
1686: + dim          - the spatial dimension
1687: . Nf           - the number of fields
1688: . NfAux        - the number of auxiliary fields
1689: . uOff         - the offset into u[] and u_t[] for each field
1690: . uOff_x       - the offset into u_x[] for each field
1691: . u            - each field evaluated at the current point
1692: . u_t          - the time derivative of each field evaluated at the current point
1693: . u_x          - the gradient of each field evaluated at the current point
1694: . aOff         - the offset into a[] and a_t[] for each auxiliary field
1695: . aOff_x       - the offset into a_x[] for each auxiliary field
1696: . a            - each auxiliary field evaluated at the current point
1697: . a_t          - the time derivative of each auxiliary field evaluated at the current point
1698: . a_x          - the gradient of auxiliary each field evaluated at the current point
1699: . t            - current time
1700: . u_tShift     - the multiplier a for dF/dU_t
1701: . x            - coordinates of the current point
1702: . numConstants - number of constant parameters
1703: . constants    - constant parameters
1704: - g0           - output values at the current point

1706:   Level: intermediate

1708:   Note:
1709:   `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.

1711:   We are using a first order FEM model for the weak form\:
1712:   \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

1714: .seealso: `PetscDS`, `PetscDSGetJacobian()`
1715: @*/
1716: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1717: {
1718:   PetscFunctionBegin;
1724:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1725:   PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1726:   PetscCall(PetscWeakFormSetIndexJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1727:   PetscFunctionReturn(PETSC_SUCCESS);
1728: }

1730: /*@C
1731:   PetscDSUseJacobianPreconditioner - Set whether to construct a Jacobian preconditioner

1733:   Not Collective

1735:   Input Parameters:
1736: + prob      - The `PetscDS`
1737: - useJacPre - flag that enables construction of a Jacobian preconditioner

1739:   Level: intermediate

1741: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1742: @*/
1743: PetscErrorCode PetscDSUseJacobianPreconditioner(PetscDS prob, PetscBool useJacPre)
1744: {
1745:   PetscFunctionBegin;
1747:   prob->useJacPre = useJacPre;
1748:   PetscFunctionReturn(PETSC_SUCCESS);
1749: }

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

1754:   Not Collective

1756:   Input Parameter:
1757: . ds - The `PetscDS`

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

1762:   Level: intermediate

1764: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1765: @*/
1766: PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS ds, PetscBool *hasJacPre)
1767: {
1768:   PetscFunctionBegin;
1770:   *hasJacPre = PETSC_FALSE;
1771:   if (!ds->useJacPre) PetscFunctionReturn(PETSC_SUCCESS);
1772:   PetscCall(PetscWeakFormHasJacobianPreconditioner(ds->wf, hasJacPre));
1773:   PetscFunctionReturn(PETSC_SUCCESS);
1774: }

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

1780:   Not Collective

1782:   Input Parameters:
1783: + ds - The `PetscDS`
1784: . f  - The test field number
1785: - g  - The field number

1787:   Output Parameters:
1788: + g0 - integrand for the test and basis function term
1789: . g1 - integrand for the test function and basis function gradient term
1790: . g2 - integrand for the test function gradient and basis function term
1791: - g3 - integrand for the test function gradient and basis function gradient term

1793:   Calling sequence of `g0`:
1794: + dim          - the spatial dimension
1795: . Nf           - the number of fields
1796: . NfAux        - the number of auxiliary fields
1797: . uOff         - the offset into u[] and u_t[] for each field
1798: . uOff_x       - the offset into u_x[] for each field
1799: . u            - each field evaluated at the current point
1800: . u_t          - the time derivative of each field evaluated at the current point
1801: . u_x          - the gradient of each field evaluated at the current point
1802: . aOff         - the offset into a[] and a_t[] for each auxiliary field
1803: . aOff_x       - the offset into a_x[] for each auxiliary field
1804: . a            - each auxiliary field evaluated at the current point
1805: . a_t          - the time derivative of each auxiliary field evaluated at the current point
1806: . a_x          - the gradient of auxiliary each field evaluated at the current point
1807: . t            - current time
1808: . u_tShift     - the multiplier a for dF/dU_t
1809: . x            - coordinates of the current point
1810: . numConstants - number of constant parameters
1811: . constants    - constant parameters
1812: - g0           - output values at the current point

1814:   Level: intermediate

1816:   Note:
1817:   `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1818:   We are using a first order FEM model for the weak form\:
1819:   \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

1821: .seealso: `PetscDS`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1822: @*/
1823: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1824: {
1825:   PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1826:   PetscInt       n0, n1, n2, n3;

1828:   PetscFunctionBegin;
1830:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1831:   PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf);
1832:   PetscCall(PetscWeakFormGetJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1833:   *g0 = tmp0 ? tmp0[0] : NULL;
1834:   *g1 = tmp1 ? tmp1[0] : NULL;
1835:   *g2 = tmp2 ? tmp2[0] : NULL;
1836:   *g3 = tmp3 ? tmp3[0] : NULL;
1837:   PetscFunctionReturn(PETSC_SUCCESS);
1838: }

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

1844:   Not Collective

1846:   Input Parameters:
1847: + ds - The `PetscDS`
1848: . f  - The test field number
1849: . g  - The field number
1850: . g0 - integrand for the test and basis function term
1851: . g1 - integrand for the test function and basis function gradient term
1852: . g2 - integrand for the test function gradient and basis function term
1853: - g3 - integrand for the test function gradient and basis function gradient term

1855:   Calling sequence of `g0`:
1856: + dim          - the spatial dimension
1857: . Nf           - the number of fields
1858: . NfAux        - the number of auxiliary fields
1859: . uOff         - the offset into u[] and u_t[] for each field
1860: . uOff_x       - the offset into u_x[] for each field
1861: . u            - each field evaluated at the current point
1862: . u_t          - the time derivative of each field evaluated at the current point
1863: . u_x          - the gradient of each field evaluated at the current point
1864: . aOff         - the offset into a[] and a_t[] for each auxiliary field
1865: . aOff_x       - the offset into a_x[] for each auxiliary field
1866: . a            - each auxiliary field evaluated at the current point
1867: . a_t          - the time derivative of each auxiliary field evaluated at the current point
1868: . a_x          - the gradient of auxiliary each field evaluated at the current point
1869: . t            - current time
1870: . u_tShift     - the multiplier a for dF/dU_t
1871: . x            - coordinates of the current point
1872: . numConstants - number of constant parameters
1873: . constants    - constant parameters
1874: - g0           - output values at the current point

1876:   Level: intermediate

1878:   Note:
1879:   `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.

1881:   We are using a first order FEM model for the weak form\:
1882:   \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

1884: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobian()`
1885: @*/
1886: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1887: {
1888:   PetscFunctionBegin;
1894:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1895:   PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1896:   PetscCall(PetscWeakFormSetIndexJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1897:   PetscFunctionReturn(PETSC_SUCCESS);
1898: }

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

1903:   Not Collective

1905:   Input Parameter:
1906: . ds - The `PetscDS`

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

1911:   Level: intermediate

1913: .seealso: `PetscDS`, `PetscDSGetDynamicJacobian()`, `PetscDSSetDynamicJacobian()`, `PetscDSGetJacobian()`
1914: @*/
1915: PetscErrorCode PetscDSHasDynamicJacobian(PetscDS ds, PetscBool *hasDynJac)
1916: {
1917:   PetscFunctionBegin;
1919:   PetscCall(PetscWeakFormHasDynamicJacobian(ds->wf, hasDynJac));
1920:   PetscFunctionReturn(PETSC_SUCCESS);
1921: }

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

1926:   Not Collective

1928:   Input Parameters:
1929: + ds - The `PetscDS`
1930: . f  - The test field number
1931: - g  - The field number

1933:   Output Parameters:
1934: + g0 - integrand for the test and basis function term
1935: . g1 - integrand for the test function and basis function gradient term
1936: . g2 - integrand for the test function gradient and basis function term
1937: - g3 - integrand for the test function gradient and basis function gradient term

1939:   Calling sequence of `g0`:
1940: + dim          - the spatial dimension
1941: . Nf           - the number of fields
1942: . NfAux        - the number of auxiliary fields
1943: . uOff         - the offset into u[] and u_t[] for each field
1944: . uOff_x       - the offset into u_x[] for each field
1945: . u            - each field evaluated at the current point
1946: . u_t          - the time derivative of each field evaluated at the current point
1947: . u_x          - the gradient of each field evaluated at the current point
1948: . aOff         - the offset into a[] and a_t[] for each auxiliary field
1949: . aOff_x       - the offset into a_x[] for each auxiliary field
1950: . a            - each auxiliary field evaluated at the current point
1951: . a_t          - the time derivative of each auxiliary field evaluated at the current point
1952: . a_x          - the gradient of auxiliary each field evaluated at the current point
1953: . t            - current time
1954: . u_tShift     - the multiplier a for dF/dU_t
1955: . x            - coordinates of the current point
1956: . numConstants - number of constant parameters
1957: . constants    - constant parameters
1958: - g0           - output values at the current point

1960:   Level: intermediate

1962:   Note:
1963:   `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.

1965:   We are using a first order FEM model for the weak form\:
1966:   \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

1968: .seealso: `PetscDS`, `PetscDSSetJacobian()`
1969: @*/
1970: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1971: {
1972:   PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1973:   PetscInt       n0, n1, n2, n3;

1975:   PetscFunctionBegin;
1977:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1978:   PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf);
1979:   PetscCall(PetscWeakFormGetDynamicJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1980:   *g0 = tmp0 ? tmp0[0] : NULL;
1981:   *g1 = tmp1 ? tmp1[0] : NULL;
1982:   *g2 = tmp2 ? tmp2[0] : NULL;
1983:   *g3 = tmp3 ? tmp3[0] : NULL;
1984:   PetscFunctionReturn(PETSC_SUCCESS);
1985: }

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

1990:   Not Collective

1992:   Input Parameters:
1993: + ds - The `PetscDS`
1994: . f  - The test field number
1995: . g  - The field number
1996: . g0 - integrand for the test and basis function term
1997: . g1 - integrand for the test function and basis function gradient term
1998: . g2 - integrand for the test function gradient and basis function term
1999: - g3 - integrand for the test function gradient and basis function gradient term

2001:   Calling sequence of `g0`:
2002: + dim          - the spatial dimension
2003: . Nf           - the number of fields
2004: . NfAux        - the number of auxiliary fields
2005: . uOff         - the offset into u[] and u_t[] for each field
2006: . uOff_x       - the offset into u_x[] for each field
2007: . u            - each field evaluated at the current point
2008: . u_t          - the time derivative of each field evaluated at the current point
2009: . u_x          - the gradient of each field evaluated at the current point
2010: . aOff         - the offset into a[] and a_t[] for each auxiliary field
2011: . aOff_x       - the offset into a_x[] for each auxiliary field
2012: . a            - each auxiliary field evaluated at the current point
2013: . a_t          - the time derivative of each auxiliary field evaluated at the current point
2014: . a_x          - the gradient of auxiliary each field evaluated at the current point
2015: . t            - current time
2016: . u_tShift     - the multiplier a for dF/dU_t
2017: . x            - coordinates of the current point
2018: . numConstants - number of constant parameters
2019: . constants    - constant parameters
2020: - g0           - output values at the current point

2022:   Level: intermediate

2024:   Note:
2025:   `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.

2027:   We are using a first order FEM model for the weak form\:
2028:   \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

2030: .seealso: `PetscDS`, `PetscDSGetJacobian()`
2031: @*/
2032: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
2033: {
2034:   PetscFunctionBegin;
2040:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2041:   PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2042:   PetscCall(PetscWeakFormSetIndexDynamicJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2043:   PetscFunctionReturn(PETSC_SUCCESS);
2044: }

2046: /*@C
2047:   PetscDSGetRiemannSolver - Returns the Riemann solver for the given field

2049:   Not Collective

2051:   Input Parameters:
2052: + ds - The `PetscDS` object
2053: - f  - The field number

2055:   Output Parameter:
2056: . r - Riemann solver

2058:   Calling sequence of `r`:
2059: + dim          - The spatial dimension
2060: . Nf           - The number of fields
2061: . x            - The coordinates at a point on the interface
2062: . n            - The normal vector to the interface
2063: . uL           - The state vector to the left of the interface
2064: . uR           - The state vector to the right of the interface
2065: . flux         - output array of flux through the interface
2066: . numConstants - number of constant parameters
2067: . constants    - constant parameters
2068: - ctx          - optional user context

2070:   Level: intermediate

2072: .seealso: `PetscDS`, `PetscDSSetRiemannSolver()`
2073: @*/
2074: 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))
2075: {
2076:   PetscRiemannFunc *tmp;
2077:   PetscInt          n;

2079:   PetscFunctionBegin;
2081:   PetscAssertPointer(r, 3);
2082:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
2083:   PetscCall(PetscWeakFormGetRiemannSolver(ds->wf, NULL, 0, f, 0, &n, &tmp));
2084:   *r = tmp ? tmp[0] : NULL;
2085:   PetscFunctionReturn(PETSC_SUCCESS);
2086: }

2088: /*@C
2089:   PetscDSSetRiemannSolver - Sets the Riemann solver for the given field

2091:   Not Collective

2093:   Input Parameters:
2094: + ds - The `PetscDS` object
2095: . f  - The field number
2096: - r  - Riemann solver

2098:   Calling sequence of `r`:
2099: + dim          - The spatial dimension
2100: . Nf           - The number of fields
2101: . x            - The coordinates at a point on the interface
2102: . n            - The normal vector to the interface
2103: . uL           - The state vector to the left of the interface
2104: . uR           - The state vector to the right of the interface
2105: . flux         - output array of flux through the interface
2106: . numConstants - number of constant parameters
2107: . constants    - constant parameters
2108: - ctx          - optional user context

2110:   Level: intermediate

2112: .seealso: `PetscDS`, `PetscDSGetRiemannSolver()`
2113: @*/
2114: 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))
2115: {
2116:   PetscFunctionBegin;
2119:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2120:   PetscCall(PetscWeakFormSetIndexRiemannSolver(ds->wf, NULL, 0, f, 0, 0, r));
2121:   PetscFunctionReturn(PETSC_SUCCESS);
2122: }

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

2127:   Not Collective

2129:   Input Parameters:
2130: + ds - The `PetscDS`
2131: - f  - The field number

2133:   Output Parameter:
2134: . update - update function

2136:   Calling sequence of `update`:
2137: + dim          - the spatial dimension
2138: . Nf           - the number of fields
2139: . NfAux        - the number of auxiliary fields
2140: . uOff         - the offset into u[] and u_t[] for each field
2141: . uOff_x       - the offset into u_x[] for each field
2142: . u            - each field evaluated at the current point
2143: . u_t          - the time derivative of each field evaluated at the current point
2144: . u_x          - the gradient of each field evaluated at the current point
2145: . aOff         - the offset into a[] and a_t[] for each auxiliary field
2146: . aOff_x       - the offset into a_x[] for each auxiliary field
2147: . a            - each auxiliary field evaluated at the current point
2148: . a_t          - the time derivative of each auxiliary field evaluated at the current point
2149: . a_x          - the gradient of auxiliary each field evaluated at the current point
2150: . t            - current time
2151: . x            - coordinates of the current point
2152: . numConstants - number of constant parameters
2153: . constants    - constant parameters
2154: - uNew         - new value for field at the current point

2156:   Level: intermediate

2158: .seealso: `PetscDS`, `PetscDSSetUpdate()`, `PetscDSSetResidual()`
2159: @*/
2160: 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[]))
2161: {
2162:   PetscFunctionBegin;
2164:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
2165:   if (update) {
2166:     PetscAssertPointer(update, 3);
2167:     *update = ds->update[f];
2168:   }
2169:   PetscFunctionReturn(PETSC_SUCCESS);
2170: }

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

2175:   Not Collective

2177:   Input Parameters:
2178: + ds     - The `PetscDS`
2179: . f      - The field number
2180: - update - update function

2182:   Calling sequence of `update`:
2183: + dim          - the spatial dimension
2184: . Nf           - the number of fields
2185: . NfAux        - the number of auxiliary fields
2186: . uOff         - the offset into u[] and u_t[] for each field
2187: . uOff_x       - the offset into u_x[] for each field
2188: . u            - each field evaluated at the current point
2189: . u_t          - the time derivative of each field evaluated at the current point
2190: . u_x          - the gradient of each field evaluated at the current point
2191: . aOff         - the offset into a[] and a_t[] for each auxiliary field
2192: . aOff_x       - the offset into a_x[] for each auxiliary field
2193: . a            - each auxiliary field evaluated at the current point
2194: . a_t          - the time derivative of each auxiliary field evaluated at the current point
2195: . a_x          - the gradient of auxiliary each field evaluated at the current point
2196: . t            - current time
2197: . x            - coordinates of the current point
2198: . numConstants - number of constant parameters
2199: . constants    - constant parameters
2200: - uNew         - new field values at the current point

2202:   Level: intermediate

2204: .seealso: `PetscDS`, `PetscDSGetResidual()`
2205: @*/
2206: 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[]))
2207: {
2208:   PetscFunctionBegin;
2211:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2212:   PetscCall(PetscDSEnlarge_Static(ds, f + 1));
2213:   ds->update[f] = update;
2214:   PetscFunctionReturn(PETSC_SUCCESS);
2215: }

2217: PetscErrorCode PetscDSGetContext(PetscDS ds, PetscInt f, void *ctx)
2218: {
2219:   PetscFunctionBegin;
2221:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
2222:   PetscAssertPointer(ctx, 3);
2223:   *(void **)ctx = ds->ctx[f];
2224:   PetscFunctionReturn(PETSC_SUCCESS);
2225: }

2227: PetscErrorCode PetscDSSetContext(PetscDS ds, PetscInt f, void *ctx)
2228: {
2229:   PetscFunctionBegin;
2231:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2232:   PetscCall(PetscDSEnlarge_Static(ds, f + 1));
2233:   ds->ctx[f] = ctx;
2234:   PetscFunctionReturn(PETSC_SUCCESS);
2235: }

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

2240:   Not Collective

2242:   Input Parameters:
2243: + ds - The PetscDS
2244: - f  - The test field number

2246:   Output Parameters:
2247: + f0 - boundary integrand for the test function term
2248: - f1 - boundary integrand for the test function gradient term

2250:   Calling sequence of `f0`:
2251: + dim          - the spatial dimension
2252: . Nf           - the number of fields
2253: . NfAux        - the number of auxiliary fields
2254: . uOff         - the offset into u[] and u_t[] for each field
2255: . uOff_x       - the offset into u_x[] for each field
2256: . u            - each field evaluated at the current point
2257: . u_t          - the time derivative of each field evaluated at the current point
2258: . u_x          - the gradient of each field evaluated at the current point
2259: . aOff         - the offset into a[] and a_t[] for each auxiliary field
2260: . aOff_x       - the offset into a_x[] for each auxiliary field
2261: . a            - each auxiliary field evaluated at the current point
2262: . a_t          - the time derivative of each auxiliary field evaluated at the current point
2263: . a_x          - the gradient of auxiliary each field evaluated at the current point
2264: . t            - current time
2265: . x            - coordinates of the current point
2266: . n            - unit normal at the current point
2267: . numConstants - number of constant parameters
2268: . constants    - constant parameters
2269: - f0           - output values at the current point

2271:   Level: intermediate

2273:   Note:
2274:   The calling sequence of `f1` is identical, and therefore omitted for brevity.

2276:   We are using a first order FEM model for the weak form\:
2277:   \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

2279: .seealso: `PetscDS`, `PetscDSSetBdResidual()`
2280: @*/
2281: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
2282: {
2283:   PetscBdPointFunc *tmp0, *tmp1;
2284:   PetscInt          n0, n1;

2286:   PetscFunctionBegin;
2288:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
2289:   PetscCall(PetscWeakFormGetBdResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1));
2290:   *f0 = tmp0 ? tmp0[0] : NULL;
2291:   *f1 = tmp1 ? tmp1[0] : NULL;
2292:   PetscFunctionReturn(PETSC_SUCCESS);
2293: }

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

2298:   Not Collective

2300:   Input Parameters:
2301: + ds - The `PetscDS`
2302: . f  - The test field number
2303: . f0 - boundary integrand for the test function term
2304: - f1 - boundary integrand for the test function gradient term

2306:   Calling sequence of `f0`:
2307: + dim          - the spatial dimension
2308: . Nf           - the number of fields
2309: . NfAux        - the number of auxiliary fields
2310: . uOff         - the offset into u[] and u_t[] for each field
2311: . uOff_x       - the offset into u_x[] for each field
2312: . u            - each field evaluated at the current point
2313: . u_t          - the time derivative of each field evaluated at the current point
2314: . u_x          - the gradient of each field evaluated at the current point
2315: . aOff         - the offset into a[] and a_t[] for each auxiliary field
2316: . aOff_x       - the offset into a_x[] for each auxiliary field
2317: . a            - each auxiliary field evaluated at the current point
2318: . a_t          - the time derivative of each auxiliary field evaluated at the current point
2319: . a_x          - the gradient of auxiliary each field evaluated at the current point
2320: . t            - current time
2321: . x            - coordinates of the current point
2322: . n            - unit normal at the current point
2323: . numConstants - number of constant parameters
2324: . constants    - constant parameters
2325: - f0           - output values at the current point

2327:   Level: intermediate

2329:   Note:
2330:   The calling sequence of `f1` is identical, and therefore omitted for brevity.

2332:   We are using a first order FEM model for the weak form\:
2333:   \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

2335: .seealso: `PetscDS`, `PetscDSGetBdResidual()`
2336: @*/
2337: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
2338: {
2339:   PetscFunctionBegin;
2341:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2342:   PetscCall(PetscWeakFormSetIndexBdResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1));
2343:   PetscFunctionReturn(PETSC_SUCCESS);
2344: }

2346: /*@
2347:   PetscDSHasBdJacobian - Indicates that boundary Jacobian functions have been set

2349:   Not Collective

2351:   Input Parameter:
2352: . ds - The `PetscDS`

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

2357:   Level: intermediate

2359: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2360: @*/
2361: PetscErrorCode PetscDSHasBdJacobian(PetscDS ds, PetscBool *hasBdJac)
2362: {
2363:   PetscFunctionBegin;
2365:   PetscAssertPointer(hasBdJac, 2);
2366:   PetscCall(PetscWeakFormHasBdJacobian(ds->wf, hasBdJac));
2367:   PetscFunctionReturn(PETSC_SUCCESS);
2368: }

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

2373:   Not Collective

2375:   Input Parameters:
2376: + ds - The `PetscDS`
2377: . f  - The test field number
2378: - g  - The field number

2380:   Output Parameters:
2381: + g0 - integrand for the test and basis function term
2382: . g1 - integrand for the test function and basis function gradient term
2383: . g2 - integrand for the test function gradient and basis function term
2384: - g3 - integrand for the test function gradient and basis function gradient term

2386:   Calling sequence of `g0`:
2387: + dim          - the spatial dimension
2388: . Nf           - the number of fields
2389: . NfAux        - the number of auxiliary fields
2390: . uOff         - the offset into u[] and u_t[] for each field
2391: . uOff_x       - the offset into u_x[] for each field
2392: . u            - each field evaluated at the current point
2393: . u_t          - the time derivative of each field evaluated at the current point
2394: . u_x          - the gradient of each field evaluated at the current point
2395: . aOff         - the offset into a[] and a_t[] for each auxiliary field
2396: . aOff_x       - the offset into a_x[] for each auxiliary field
2397: . a            - each auxiliary field evaluated at the current point
2398: . a_t          - the time derivative of each auxiliary field evaluated at the current point
2399: . a_x          - the gradient of auxiliary each field evaluated at the current point
2400: . t            - current time
2401: . u_tShift     - the multiplier a for dF/dU_t
2402: . x            - coordinates of the current point
2403: . n            - normal at the current point
2404: . numConstants - number of constant parameters
2405: . constants    - constant parameters
2406: - g0           - output values at the current point

2408:   Level: intermediate

2410:   Note:
2411:   `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.

2413:   We are using a first order FEM model for the weak form\:
2414:   \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

2416: .seealso: `PetscDS`, `PetscDSSetBdJacobian()`
2417: @*/
2418: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
2419: {
2420:   PetscBdPointJac *tmp0, *tmp1, *tmp2, *tmp3;
2421:   PetscInt         n0, n1, n2, n3;

2423:   PetscFunctionBegin;
2425:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
2426:   PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf);
2427:   PetscCall(PetscWeakFormGetBdJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2428:   *g0 = tmp0 ? tmp0[0] : NULL;
2429:   *g1 = tmp1 ? tmp1[0] : NULL;
2430:   *g2 = tmp2 ? tmp2[0] : NULL;
2431:   *g3 = tmp3 ? tmp3[0] : NULL;
2432:   PetscFunctionReturn(PETSC_SUCCESS);
2433: }

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

2438:   Not Collective

2440:   Input Parameters:
2441: + ds - The PetscDS
2442: . f  - The test field number
2443: . g  - The field number
2444: . g0 - integrand for the test and basis function term
2445: . g1 - integrand for the test function and basis function gradient term
2446: . g2 - integrand for the test function gradient and basis function term
2447: - g3 - integrand for the test function gradient and basis function gradient term

2449:   Calling sequence of `g0`:
2450: + dim          - the spatial dimension
2451: . Nf           - the number of fields
2452: . NfAux        - the number of auxiliary fields
2453: . uOff         - the offset into u[] and u_t[] for each field
2454: . uOff_x       - the offset into u_x[] for each field
2455: . u            - each field evaluated at the current point
2456: . u_t          - the time derivative of each field evaluated at the current point
2457: . u_x          - the gradient of each field evaluated at the current point
2458: . aOff         - the offset into a[] and a_t[] for each auxiliary field
2459: . aOff_x       - the offset into a_x[] for each auxiliary field
2460: . a            - each auxiliary field evaluated at the current point
2461: . a_t          - the time derivative of each auxiliary field evaluated at the current point
2462: . a_x          - the gradient of auxiliary each field evaluated at the current point
2463: . t            - current time
2464: . u_tShift     - the multiplier a for dF/dU_t
2465: . x            - coordinates of the current point
2466: . n            - normal at the current point
2467: . numConstants - number of constant parameters
2468: . constants    - constant parameters
2469: - g0           - output values at the current point

2471:   Level: intermediate

2473:   Note:
2474:   `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.

2476:   We are using a first order FEM model for the weak form\:
2477:   \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

2479: .seealso: `PetscDS`, `PetscDSGetBdJacobian()`
2480: @*/
2481: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
2482: {
2483:   PetscFunctionBegin;
2489:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2490:   PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2491:   PetscCall(PetscWeakFormSetIndexBdJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2492:   PetscFunctionReturn(PETSC_SUCCESS);
2493: }

2495: /*@
2496:   PetscDSHasBdJacobianPreconditioner - Signals that boundary Jacobian preconditioner functions have been set

2498:   Not Collective

2500:   Input Parameter:
2501: . ds - The `PetscDS`

2503:   Output Parameter:
2504: . hasBdJacPre - flag that pointwise function for the boundary Jacobian preconditioner has been set

2506:   Level: intermediate

2508: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2509: @*/
2510: PetscErrorCode PetscDSHasBdJacobianPreconditioner(PetscDS ds, PetscBool *hasBdJacPre)
2511: {
2512:   PetscFunctionBegin;
2514:   PetscAssertPointer(hasBdJacPre, 2);
2515:   PetscCall(PetscWeakFormHasBdJacobianPreconditioner(ds->wf, hasBdJacPre));
2516:   PetscFunctionReturn(PETSC_SUCCESS);
2517: }

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

2522:   Not Collective; No Fortran Support

2524:   Input Parameters:
2525: + ds - The `PetscDS`
2526: . f  - The test field number
2527: - g  - The field number

2529:   Output Parameters:
2530: + g0 - integrand for the test and basis function term
2531: . g1 - integrand for the test function and basis function gradient term
2532: . g2 - integrand for the test function gradient and basis function term
2533: - g3 - integrand for the test function gradient and basis function gradient term

2535:   Calling sequence of `g0`:
2536: + dim          - the spatial dimension
2537: . Nf           - the number of fields
2538: . NfAux        - the number of auxiliary fields
2539: . uOff         - the offset into u[] and u_t[] for each field
2540: . uOff_x       - the offset into u_x[] for each field
2541: . u            - each field evaluated at the current point
2542: . u_t          - the time derivative of each field evaluated at the current point
2543: . u_x          - the gradient of each field evaluated at the current point
2544: . aOff         - the offset into a[] and a_t[] for each auxiliary field
2545: . aOff_x       - the offset into a_x[] for each auxiliary field
2546: . a            - each auxiliary field evaluated at the current point
2547: . a_t          - the time derivative of each auxiliary field evaluated at the current point
2548: . a_x          - the gradient of auxiliary each field evaluated at the current point
2549: . t            - current time
2550: . u_tShift     - the multiplier a for dF/dU_t
2551: . x            - coordinates of the current point
2552: . n            - normal at the current point
2553: . numConstants - number of constant parameters
2554: . constants    - constant parameters
2555: - g0           - output values at the current point

2557:   Level: intermediate

2559:   Note:
2560:   `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.

2562:   We are using a first order FEM model for the weak form\:
2563:   \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

2565: .seealso: `PetscDS`, `PetscDSSetBdJacobianPreconditioner()`
2566: @*/
2567: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
2568: {
2569:   PetscBdPointJac *tmp0, *tmp1, *tmp2, *tmp3;
2570:   PetscInt         n0, n1, n2, n3;

2572:   PetscFunctionBegin;
2574:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
2575:   PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf);
2576:   PetscCall(PetscWeakFormGetBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2577:   *g0 = tmp0 ? tmp0[0] : NULL;
2578:   *g1 = tmp1 ? tmp1[0] : NULL;
2579:   *g2 = tmp2 ? tmp2[0] : NULL;
2580:   *g3 = tmp3 ? tmp3[0] : NULL;
2581:   PetscFunctionReturn(PETSC_SUCCESS);
2582: }

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

2587:   Not Collective; No Fortran Support

2589:   Input Parameters:
2590: + ds - The `PetscDS`
2591: . f  - The test field number
2592: . g  - The field number
2593: . g0 - integrand for the test and basis function term
2594: . g1 - integrand for the test function and basis function gradient term
2595: . g2 - integrand for the test function gradient and basis function term
2596: - g3 - integrand for the test function gradient and basis function gradient term

2598:   Calling sequence of `g0':
2599: + dim          - the spatial dimension
2600: . Nf           - the number of fields
2601: . NfAux        - the number of auxiliary fields
2602: . uOff         - the offset into u[] and u_t[] for each field
2603: . uOff_x       - the offset into u_x[] for each field
2604: . u            - each field evaluated at the current point
2605: . u_t          - the time derivative of each field evaluated at the current point
2606: . u_x          - the gradient of each field evaluated at the current point
2607: . aOff         - the offset into a[] and a_t[] for each auxiliary field
2608: . aOff_x       - the offset into a_x[] for each auxiliary field
2609: . a            - each auxiliary field evaluated at the current point
2610: . a_t          - the time derivative of each auxiliary field evaluated at the current point
2611: . a_x          - the gradient of auxiliary each field evaluated at the current point
2612: . t            - current time
2613: . u_tShift     - the multiplier a for dF/dU_t
2614: . x            - coordinates of the current point
2615: . n            - normal at the current point
2616: . numConstants - number of constant parameters
2617: . constants    - constant parameters
2618: - g0           - output values at the current point

2620:   Level: intermediate

2622:   Note:
2623:   `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.

2625:   We are using a first order FEM model for the weak form\:
2626:   \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

2628: .seealso: `PetscDS`, `PetscDSGetBdJacobianPreconditioner()`
2629: @*/
2630: 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, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
2631: {
2632:   PetscFunctionBegin;
2638:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2639:   PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2640:   PetscCall(PetscWeakFormSetIndexBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2641:   PetscFunctionReturn(PETSC_SUCCESS);
2642: }

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

2647:   Not Collective

2649:   Input Parameters:
2650: + prob - The PetscDS
2651: - f    - The test field number

2653:   Output Parameters:
2654: + sol - exact solution for the test field
2655: - ctx - exact solution context

2657:   Calling sequence of `exactSol`:
2658: + dim - the spatial dimension
2659: . t   - current time
2660: . x   - coordinates of the current point
2661: . Nc  - the number of field components
2662: . u   - the solution field evaluated at the current point
2663: - ctx - a user context

2665:   Level: intermediate

2667: .seealso: `PetscDS`, `PetscDSSetExactSolution()`, `PetscDSGetExactSolutionTimeDerivative()`
2668: @*/
2669: PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx)
2670: {
2671:   PetscFunctionBegin;
2673:   PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
2674:   if (sol) {
2675:     PetscAssertPointer(sol, 3);
2676:     *sol = prob->exactSol[f];
2677:   }
2678:   if (ctx) {
2679:     PetscAssertPointer(ctx, 4);
2680:     *ctx = prob->exactCtx[f];
2681:   }
2682:   PetscFunctionReturn(PETSC_SUCCESS);
2683: }

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

2688:   Not Collective

2690:   Input Parameters:
2691: + prob - The `PetscDS`
2692: . f    - The test field number
2693: . sol  - solution function for the test fields
2694: - ctx  - solution context or `NULL`

2696:   Calling sequence of `sol`:
2697: .vb
2698:   PetscErrorCode sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2699: .ve
2700: + dim - the spatial dimension
2701: . t   - current time
2702: . x   - coordinates of the current point
2703: . Nc  - the number of field components
2704: . u   - the solution field evaluated at the current point
2705: - ctx - a user context

2707:   Level: intermediate

2709: .seealso: `PetscDS`, `PetscDSGetExactSolution()`
2710: @*/
2711: PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx)
2712: {
2713:   PetscFunctionBegin;
2715:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2716:   PetscCall(PetscDSEnlarge_Static(prob, f + 1));
2717:   if (sol) {
2719:     prob->exactSol[f] = sol;
2720:   }
2721:   if (ctx) {
2723:     prob->exactCtx[f] = ctx;
2724:   }
2725:   PetscFunctionReturn(PETSC_SUCCESS);
2726: }

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

2731:   Not Collective

2733:   Input Parameters:
2734: + prob - The `PetscDS`
2735: - f    - The test field number

2737:   Output Parameters:
2738: + sol - time derivative of the exact solution for the test field
2739: - ctx - time derivative of the exact solution context

2741:   Calling sequence of `exactSol`:
2742: + dim - the spatial dimension
2743: . t   - current time
2744: . x   - coordinates of the current point
2745: . Nc  - the number of field components
2746: . u   - the solution field evaluated at the current point
2747: - ctx - a user context

2749:   Level: intermediate

2751: .seealso: `PetscDS`, `PetscDSSetExactSolutionTimeDerivative()`, `PetscDSGetExactSolution()`
2752: @*/
2753: PetscErrorCode PetscDSGetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx)
2754: {
2755:   PetscFunctionBegin;
2757:   PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
2758:   if (sol) {
2759:     PetscAssertPointer(sol, 3);
2760:     *sol = prob->exactSol_t[f];
2761:   }
2762:   if (ctx) {
2763:     PetscAssertPointer(ctx, 4);
2764:     *ctx = prob->exactCtx_t[f];
2765:   }
2766:   PetscFunctionReturn(PETSC_SUCCESS);
2767: }

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

2772:   Not Collective

2774:   Input Parameters:
2775: + prob - The `PetscDS`
2776: . f    - The test field number
2777: . sol  - time derivative of the solution function for the test fields
2778: - ctx  - time derivative of the solution context or `NULL`

2780:   Calling sequence of `sol`:
2781: .vb
2782:   PetscErrorCode sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2783: .ve
2784: + dim - the spatial dimension
2785: . t   - current time
2786: . x   - coordinates of the current point
2787: . Nc  - the number of field components
2788: . u   - the solution field evaluated at the current point
2789: - ctx - a user context

2791:   Level: intermediate

2793: .seealso: `PetscDS`, `PetscDSGetExactSolutionTimeDerivative()`, `PetscDSSetExactSolution()`
2794: @*/
2795: PetscErrorCode PetscDSSetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx)
2796: {
2797:   PetscFunctionBegin;
2799:   PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2800:   PetscCall(PetscDSEnlarge_Static(prob, f + 1));
2801:   if (sol) {
2803:     prob->exactSol_t[f] = sol;
2804:   }
2805:   if (ctx) {
2807:     prob->exactCtx_t[f] = ctx;
2808:   }
2809:   PetscFunctionReturn(PETSC_SUCCESS);
2810: }

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

2815:   Not Collective

2817:   Input Parameter:
2818: . prob - The `PetscDS` object

2820:   Output Parameters:
2821: + numConstants - The number of constants
2822: - constants    - The array of constants, NULL if there are none

2824:   Level: intermediate

2826: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSCreate()`
2827: @*/
2828: PetscErrorCode PetscDSGetConstants(PetscDS prob, PetscInt *numConstants, const PetscScalar *constants[])
2829: {
2830:   PetscFunctionBegin;
2832:   if (numConstants) {
2833:     PetscAssertPointer(numConstants, 2);
2834:     *numConstants = prob->numConstants;
2835:   }
2836:   if (constants) {
2837:     PetscAssertPointer(constants, 3);
2838:     *constants = prob->constants;
2839:   }
2840:   PetscFunctionReturn(PETSC_SUCCESS);
2841: }

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

2846:   Not Collective

2848:   Input Parameters:
2849: + prob         - The `PetscDS` object
2850: . numConstants - The number of constants
2851: - constants    - The array of constants, NULL if there are none

2853:   Level: intermediate

2855: .seealso: `PetscDS`, `PetscDSGetConstants()`, `PetscDSCreate()`
2856: @*/
2857: PetscErrorCode PetscDSSetConstants(PetscDS prob, PetscInt numConstants, PetscScalar constants[])
2858: {
2859:   PetscFunctionBegin;
2861:   if (numConstants != prob->numConstants) {
2862:     PetscCall(PetscFree(prob->constants));
2863:     prob->numConstants = numConstants;
2864:     if (prob->numConstants) {
2865:       PetscCall(PetscMalloc1(prob->numConstants, &prob->constants));
2866:     } else {
2867:       prob->constants = NULL;
2868:     }
2869:   }
2870:   if (prob->numConstants) {
2871:     PetscAssertPointer(constants, 3);
2872:     PetscCall(PetscArraycpy(prob->constants, constants, prob->numConstants));
2873:   }
2874:   PetscFunctionReturn(PETSC_SUCCESS);
2875: }

2877: /*@
2878:   PetscDSGetFieldIndex - Returns the index of the given field

2880:   Not Collective

2882:   Input Parameters:
2883: + prob - The `PetscDS` object
2884: - disc - The discretization object

2886:   Output Parameter:
2887: . f - The field number

2889:   Level: beginner

2891: .seealso: `PetscDS`, `PetscGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2892: @*/
2893: PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f)
2894: {
2895:   PetscInt g;

2897:   PetscFunctionBegin;
2899:   PetscAssertPointer(f, 3);
2900:   *f = -1;
2901:   for (g = 0; g < prob->Nf; ++g) {
2902:     if (disc == prob->disc[g]) break;
2903:   }
2904:   PetscCheck(g != prob->Nf, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS.");
2905:   *f = g;
2906:   PetscFunctionReturn(PETSC_SUCCESS);
2907: }

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

2912:   Not Collective

2914:   Input Parameters:
2915: + prob - The `PetscDS` object
2916: - f    - The field number

2918:   Output Parameter:
2919: . size - The size

2921:   Level: beginner

2923: .seealso: `PetscDS`, `PetscDSGetFieldOffset()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2924: @*/
2925: PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size)
2926: {
2927:   PetscFunctionBegin;
2929:   PetscAssertPointer(size, 3);
2930:   PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
2931:   PetscCall(PetscDSSetUp(prob));
2932:   *size = prob->Nb[f];
2933:   PetscFunctionReturn(PETSC_SUCCESS);
2934: }

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

2939:   Not Collective

2941:   Input Parameters:
2942: + prob - The `PetscDS` object
2943: - f    - The field number

2945:   Output Parameter:
2946: . off - The offset

2948:   Level: beginner

2950: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2951: @*/
2952: PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off)
2953: {
2954:   PetscInt size, g;

2956:   PetscFunctionBegin;
2958:   PetscAssertPointer(off, 3);
2959:   PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
2960:   *off = 0;
2961:   for (g = 0; g < f; ++g) {
2962:     PetscCall(PetscDSGetFieldSize(prob, g, &size));
2963:     *off += size;
2964:   }
2965:   PetscFunctionReturn(PETSC_SUCCESS);
2966: }

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

2971:   Not Collective

2973:   Input Parameters:
2974: + ds - The `PetscDS` object
2975: - f  - The field number

2977:   Output Parameter:
2978: . off - The offset

2980:   Level: beginner

2982: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2983: @*/
2984: PetscErrorCode PetscDSGetFieldOffsetCohesive(PetscDS ds, PetscInt f, PetscInt *off)
2985: {
2986:   PetscInt size, g;

2988:   PetscFunctionBegin;
2990:   PetscAssertPointer(off, 3);
2991:   PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
2992:   *off = 0;
2993:   for (g = 0; g < f; ++g) {
2994:     PetscBool cohesive;

2996:     PetscCall(PetscDSGetCohesive(ds, g, &cohesive));
2997:     PetscCall(PetscDSGetFieldSize(ds, g, &size));
2998:     *off += cohesive ? size : size * 2;
2999:   }
3000:   PetscFunctionReturn(PETSC_SUCCESS);
3001: }

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

3006:   Not Collective

3008:   Input Parameter:
3009: . prob - The `PetscDS` object

3011:   Output Parameter:
3012: . dimensions - The number of dimensions

3014:   Level: beginner

3016: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3017: @*/
3018: PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[])
3019: {
3020:   PetscFunctionBegin;
3022:   PetscCall(PetscDSSetUp(prob));
3023:   PetscAssertPointer(dimensions, 2);
3024:   *dimensions = prob->Nb;
3025:   PetscFunctionReturn(PETSC_SUCCESS);
3026: }

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

3031:   Not Collective

3033:   Input Parameter:
3034: . prob - The `PetscDS` object

3036:   Output Parameter:
3037: . components - The number of components

3039:   Level: beginner

3041: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3042: @*/
3043: PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[])
3044: {
3045:   PetscFunctionBegin;
3047:   PetscCall(PetscDSSetUp(prob));
3048:   PetscAssertPointer(components, 2);
3049:   *components = prob->Nc;
3050:   PetscFunctionReturn(PETSC_SUCCESS);
3051: }

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

3056:   Not Collective

3058:   Input Parameters:
3059: + prob - The `PetscDS` object
3060: - f    - The field number

3062:   Output Parameter:
3063: . off - The offset

3065:   Level: beginner

3067: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3068: @*/
3069: PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off)
3070: {
3071:   PetscFunctionBegin;
3073:   PetscAssertPointer(off, 3);
3074:   PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
3075:   PetscCall(PetscDSSetUp(prob));
3076:   *off = prob->off[f];
3077:   PetscFunctionReturn(PETSC_SUCCESS);
3078: }

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

3083:   Not Collective

3085:   Input Parameter:
3086: . prob - The `PetscDS` object

3088:   Output Parameter:
3089: . offsets - The offsets

3091:   Level: beginner

3093: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3094: @*/
3095: PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[])
3096: {
3097:   PetscFunctionBegin;
3099:   PetscAssertPointer(offsets, 2);
3100:   PetscCall(PetscDSSetUp(prob));
3101:   *offsets = prob->off;
3102:   PetscFunctionReturn(PETSC_SUCCESS);
3103: }

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

3108:   Not Collective

3110:   Input Parameter:
3111: . prob - The `PetscDS` object

3113:   Output Parameter:
3114: . offsets - The offsets

3116:   Level: beginner

3118: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3119: @*/
3120: PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[])
3121: {
3122:   PetscFunctionBegin;
3124:   PetscAssertPointer(offsets, 2);
3125:   PetscCall(PetscDSSetUp(prob));
3126:   *offsets = prob->offDer;
3127:   PetscFunctionReturn(PETSC_SUCCESS);
3128: }

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

3133:   Not Collective

3135:   Input Parameters:
3136: + ds - The `PetscDS` object
3137: - s  - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive

3139:   Output Parameter:
3140: . offsets - The offsets

3142:   Level: beginner

3144: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3145: @*/
3146: PetscErrorCode PetscDSGetComponentOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
3147: {
3148:   PetscFunctionBegin;
3150:   PetscAssertPointer(offsets, 3);
3151:   PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS");
3152:   PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s);
3153:   PetscCall(PetscDSSetUp(ds));
3154:   *offsets = ds->offCohesive[s];
3155:   PetscFunctionReturn(PETSC_SUCCESS);
3156: }

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

3161:   Not Collective

3163:   Input Parameters:
3164: + ds - The `PetscDS` object
3165: - s  - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive

3167:   Output Parameter:
3168: . offsets - The offsets

3170:   Level: beginner

3172: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3173: @*/
3174: PetscErrorCode PetscDSGetComponentDerivativeOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
3175: {
3176:   PetscFunctionBegin;
3178:   PetscAssertPointer(offsets, 3);
3179:   PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS");
3180:   PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s);
3181:   PetscCall(PetscDSSetUp(ds));
3182:   *offsets = ds->offDerCohesive[s];
3183:   PetscFunctionReturn(PETSC_SUCCESS);
3184: }

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

3189:   Not Collective

3191:   Input Parameter:
3192: . prob - The `PetscDS` object

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

3197:   Level: intermediate

3199: .seealso: `PetscDS`, `PetscTabulation`, `PetscDSCreate()`
3200: @*/
3201: PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscTabulation *T[])
3202: {
3203:   PetscFunctionBegin;
3205:   PetscAssertPointer(T, 2);
3206:   PetscCall(PetscDSSetUp(prob));
3207:   *T = prob->T;
3208:   PetscFunctionReturn(PETSC_SUCCESS);
3209: }

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

3214:   Not Collective

3216:   Input Parameter:
3217: . prob - The `PetscDS` object

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

3222:   Level: intermediate

3224: .seealso: `PetscTabulation`, `PetscDS`, `PetscDSGetTabulation()`, `PetscDSCreate()`
3225: @*/
3226: PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscTabulation *Tf[])
3227: {
3228:   PetscFunctionBegin;
3230:   PetscAssertPointer(Tf, 2);
3231:   PetscCall(PetscDSSetUp(prob));
3232:   *Tf = prob->Tf;
3233:   PetscFunctionReturn(PETSC_SUCCESS);
3234: }

3236: PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x)
3237: {
3238:   PetscFunctionBegin;
3240:   PetscCall(PetscDSSetUp(prob));
3241:   if (u) {
3242:     PetscAssertPointer(u, 2);
3243:     *u = prob->u;
3244:   }
3245:   if (u_t) {
3246:     PetscAssertPointer(u_t, 3);
3247:     *u_t = prob->u_t;
3248:   }
3249:   if (u_x) {
3250:     PetscAssertPointer(u_x, 4);
3251:     *u_x = prob->u_x;
3252:   }
3253:   PetscFunctionReturn(PETSC_SUCCESS);
3254: }

3256: PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3)
3257: {
3258:   PetscFunctionBegin;
3260:   PetscCall(PetscDSSetUp(prob));
3261:   if (f0) {
3262:     PetscAssertPointer(f0, 2);
3263:     *f0 = prob->f0;
3264:   }
3265:   if (f1) {
3266:     PetscAssertPointer(f1, 3);
3267:     *f1 = prob->f1;
3268:   }
3269:   if (g0) {
3270:     PetscAssertPointer(g0, 4);
3271:     *g0 = prob->g0;
3272:   }
3273:   if (g1) {
3274:     PetscAssertPointer(g1, 5);
3275:     *g1 = prob->g1;
3276:   }
3277:   if (g2) {
3278:     PetscAssertPointer(g2, 6);
3279:     *g2 = prob->g2;
3280:   }
3281:   if (g3) {
3282:     PetscAssertPointer(g3, 7);
3283:     *g3 = prob->g3;
3284:   }
3285:   PetscFunctionReturn(PETSC_SUCCESS);
3286: }

3288: PetscErrorCode PetscDSGetWorkspace(PetscDS prob, PetscReal **x, PetscScalar **basisReal, PetscScalar **basisDerReal, PetscScalar **testReal, PetscScalar **testDerReal)
3289: {
3290:   PetscFunctionBegin;
3292:   PetscCall(PetscDSSetUp(prob));
3293:   if (x) {
3294:     PetscAssertPointer(x, 2);
3295:     *x = prob->x;
3296:   }
3297:   if (basisReal) {
3298:     PetscAssertPointer(basisReal, 3);
3299:     *basisReal = prob->basisReal;
3300:   }
3301:   if (basisDerReal) {
3302:     PetscAssertPointer(basisDerReal, 4);
3303:     *basisDerReal = prob->basisDerReal;
3304:   }
3305:   if (testReal) {
3306:     PetscAssertPointer(testReal, 5);
3307:     *testReal = prob->testReal;
3308:   }
3309:   if (testDerReal) {
3310:     PetscAssertPointer(testDerReal, 6);
3311:     *testDerReal = prob->testDerReal;
3312:   }
3313:   PetscFunctionReturn(PETSC_SUCCESS);
3314: }

3316: /*@C
3317:   PetscDSAddBoundary - Add a boundary condition to the model.

3319:   Collective

3321:   Input Parameters:
3322: + ds       - The PetscDS object
3323: . type     - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3324: . name     - The BC name
3325: . label    - The label defining constrained points
3326: . Nv       - The number of `DMLabel` values for constrained points
3327: . values   - An array of label values for constrained points
3328: . field    - The field to constrain
3329: . Nc       - The number of constrained field components (0 will constrain all fields)
3330: . comps    - An array of constrained component numbers
3331: . bcFunc   - A pointwise function giving boundary values
3332: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3333: - ctx      - An optional user context for bcFunc

3335:   Output Parameter:
3336: . bd - The boundary number

3338:   Options Database Keys:
3339: + -bc_<boundary name> <num>      - Overrides the boundary ids
3340: - -bc_<boundary name>_comp <num> - Overrides the boundary components

3342:   Level: developer

3344:   Note:
3345:   Both `bcFunc` and `bcFunc_t` will depend on the boundary condition type. If the type if `DM_BC_ESSENTIAL`, then the calling sequence is\:

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

3349:   If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value, then the calling sequence is\:
3350: .vb
3351:   void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3352:               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3353:               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3354:               PetscReal time, const PetscReal x[], PetscScalar bcval[])
3355: .ve
3356: + dim - the spatial dimension
3357: . Nf - the number of fields
3358: . uOff - the offset into u[] and u_t[] for each field
3359: . uOff_x - the offset into u_x[] for each field
3360: . u - each field evaluated at the current point
3361: . u_t - the time derivative of each field evaluated at the current point
3362: . u_x - the gradient of each field evaluated at the current point
3363: . aOff - the offset into a[] and a_t[] for each auxiliary field
3364: . aOff_x - the offset into a_x[] for each auxiliary field
3365: . a - each auxiliary field evaluated at the current point
3366: . a_t - the time derivative of each auxiliary field evaluated at the current point
3367: . a_x - the gradient of auxiliary each field evaluated at the current point
3368: . t - current time
3369: . x - coordinates of the current point
3370: . numConstants - number of constant parameters
3371: . constants - constant parameters
3372: - bcval - output values at the current point

3374:   Notes:
3375:   The pointwise functions are used to provide boundary values for essential boundary
3376:   conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3377:   coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3378:   integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.

3380: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundaryByName()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3381: @*/
3382: 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)
3383: {
3384:   DSBoundary  head = ds->boundary, b;
3385:   PetscInt    n    = 0;
3386:   const char *lname;

3388:   PetscFunctionBegin;
3391:   PetscAssertPointer(name, 3);
3396:   PetscCheck(field >= 0 && field < ds->Nf, PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_OUTOFRANGE, "Field %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", field, ds->Nf);
3397:   if (Nc > 0) {
3398:     PetscInt *fcomps;
3399:     PetscInt  c;

3401:     PetscCall(PetscDSGetComponents(ds, &fcomps));
3402:     PetscCheck(Nc <= fcomps[field], PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_OUTOFRANGE, "Number of constrained components %" PetscInt_FMT " > %" PetscInt_FMT " components for field %" PetscInt_FMT, Nc, fcomps[field], field);
3403:     for (c = 0; c < Nc; ++c) {
3404:       PetscCheck(comps[c] >= 0 && comps[c] < fcomps[field], PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_OUTOFRANGE, "Constrained component[%" PetscInt_FMT "] %" PetscInt_FMT " not in [0, %" PetscInt_FMT ") components for field %" PetscInt_FMT, c, comps[c], fcomps[field], field);
3405:     }
3406:   }
3407:   PetscCall(PetscNew(&b));
3408:   PetscCall(PetscStrallocpy(name, (char **)&b->name));
3409:   PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf));
3410:   PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf));
3411:   PetscCall(PetscMalloc1(Nv, &b->values));
3412:   if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3413:   PetscCall(PetscMalloc1(Nc, &b->comps));
3414:   if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3415:   PetscCall(PetscObjectGetName((PetscObject)label, &lname));
3416:   PetscCall(PetscStrallocpy(lname, (char **)&b->lname));
3417:   b->type   = type;
3418:   b->label  = label;
3419:   b->Nv     = Nv;
3420:   b->field  = field;
3421:   b->Nc     = Nc;
3422:   b->func   = bcFunc;
3423:   b->func_t = bcFunc_t;
3424:   b->ctx    = ctx;
3425:   b->next   = NULL;
3426:   /* Append to linked list so that we can preserve the order */
3427:   if (!head) ds->boundary = b;
3428:   while (head) {
3429:     if (!head->next) {
3430:       head->next = b;
3431:       head       = b;
3432:     }
3433:     head = head->next;
3434:     ++n;
3435:   }
3436:   if (bd) {
3437:     PetscAssertPointer(bd, 13);
3438:     *bd = n;
3439:   }
3440:   PetscFunctionReturn(PETSC_SUCCESS);
3441: }

3443: // PetscClangLinter pragma ignore: -fdoc-section-header-unknown
3444: /*@C
3445:   PetscDSAddBoundaryByName - Add a boundary condition to the model.

3447:   Collective

3449:   Input Parameters:
3450: + ds       - The `PetscDS` object
3451: . type     - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3452: . name     - The BC name
3453: . lname    - The naem of the label defining constrained points
3454: . Nv       - The number of `DMLabel` values for constrained points
3455: . values   - An array of label values for constrained points
3456: . field    - The field to constrain
3457: . Nc       - The number of constrained field components (0 will constrain all fields)
3458: . comps    - An array of constrained component numbers
3459: . bcFunc   - A pointwise function giving boundary values
3460: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3461: - ctx      - An optional user context for bcFunc

3463:   Output Parameter:
3464: . bd - The boundary number

3466:   Options Database Keys:
3467: + -bc_<boundary name> <num>      - Overrides the boundary ids
3468: - -bc_<boundary name>_comp <num> - Overrides the boundary components

3470:   Calling Sequence of `bcFunc` and `bcFunc_t`:
3471:   If the type is `DM_BC_ESSENTIAL`
3472: .vb
3473:   void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])
3474: .ve
3475:   If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value,
3476: .vb
3477:   void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3478:               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3479:               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3480:               PetscReal time, const PetscReal x[], PetscScalar bcval[])
3481: .ve
3482: + dim - the spatial dimension
3483: . Nf - the number of fields
3484: . uOff - the offset into u[] and u_t[] for each field
3485: . uOff_x - the offset into u_x[] for each field
3486: . u - each field evaluated at the current point
3487: . u_t - the time derivative of each field evaluated at the current point
3488: . u_x - the gradient of each field evaluated at the current point
3489: . aOff - the offset into a[] and a_t[] for each auxiliary field
3490: . aOff_x - the offset into a_x[] for each auxiliary field
3491: . a - each auxiliary field evaluated at the current point
3492: . a_t - the time derivative of each auxiliary field evaluated at the current point
3493: . a_x - the gradient of auxiliary each field evaluated at the current point
3494: . t - current time
3495: . x - coordinates of the current point
3496: . numConstants - number of constant parameters
3497: . constants - constant parameters
3498: - bcval - output values at the current point

3500:   Level: developer

3502:   Notes:
3503:   The pointwise functions are used to provide boundary values for essential boundary
3504:   conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3505:   coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3506:   integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.

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

3510: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3511: @*/
3512: 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)
3513: {
3514:   DSBoundary head = ds->boundary, b;
3515:   PetscInt   n    = 0;

3517:   PetscFunctionBegin;
3520:   PetscAssertPointer(name, 3);
3521:   PetscAssertPointer(lname, 4);
3525:   PetscCall(PetscNew(&b));
3526:   PetscCall(PetscStrallocpy(name, (char **)&b->name));
3527:   PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf));
3528:   PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf));
3529:   PetscCall(PetscMalloc1(Nv, &b->values));
3530:   if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3531:   PetscCall(PetscMalloc1(Nc, &b->comps));
3532:   if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3533:   PetscCall(PetscStrallocpy(lname, (char **)&b->lname));
3534:   b->type   = type;
3535:   b->label  = NULL;
3536:   b->Nv     = Nv;
3537:   b->field  = field;
3538:   b->Nc     = Nc;
3539:   b->func   = bcFunc;
3540:   b->func_t = bcFunc_t;
3541:   b->ctx    = ctx;
3542:   b->next   = NULL;
3543:   /* Append to linked list so that we can preserve the order */
3544:   if (!head) ds->boundary = b;
3545:   while (head) {
3546:     if (!head->next) {
3547:       head->next = b;
3548:       head       = b;
3549:     }
3550:     head = head->next;
3551:     ++n;
3552:   }
3553:   if (bd) {
3554:     PetscAssertPointer(bd, 13);
3555:     *bd = n;
3556:   }
3557:   PetscFunctionReturn(PETSC_SUCCESS);
3558: }

3560: /*@C
3561:   PetscDSUpdateBoundary - Change a boundary condition for the model.

3563:   Input Parameters:
3564: + ds       - The `PetscDS` object
3565: . bd       - The boundary condition number
3566: . type     - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3567: . name     - The BC name
3568: . label    - The label defining constrained points
3569: . Nv       - The number of `DMLabel` ids for constrained points
3570: . values   - An array of ids for constrained points
3571: . field    - The field to constrain
3572: . Nc       - The number of constrained field components
3573: . comps    - An array of constrained component numbers
3574: . bcFunc   - A pointwise function giving boundary values
3575: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3576: - ctx      - An optional user context for bcFunc

3578:   Level: developer

3580:   Notes:
3581:   The pointwise functions are used to provide boundary values for essential boundary
3582:   conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3583:   coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3584:   integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.

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

3589: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSGetNumBoundary()`, `DMLabel`
3590: @*/
3591: 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)
3592: {
3593:   DSBoundary b = ds->boundary;
3594:   PetscInt   n = 0;

3596:   PetscFunctionBegin;
3598:   while (b) {
3599:     if (n == bd) break;
3600:     b = b->next;
3601:     ++n;
3602:   }
3603:   PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n);
3604:   if (name) {
3605:     PetscCall(PetscFree(b->name));
3606:     PetscCall(PetscStrallocpy(name, (char **)&b->name));
3607:   }
3608:   b->type = type;
3609:   if (label) {
3610:     const char *name;

3612:     b->label = label;
3613:     PetscCall(PetscFree(b->lname));
3614:     PetscCall(PetscObjectGetName((PetscObject)label, &name));
3615:     PetscCall(PetscStrallocpy(name, (char **)&b->lname));
3616:   }
3617:   if (Nv >= 0) {
3618:     b->Nv = Nv;
3619:     PetscCall(PetscFree(b->values));
3620:     PetscCall(PetscMalloc1(Nv, &b->values));
3621:     if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3622:   }
3623:   if (field >= 0) b->field = field;
3624:   if (Nc >= 0) {
3625:     b->Nc = Nc;
3626:     PetscCall(PetscFree(b->comps));
3627:     PetscCall(PetscMalloc1(Nc, &b->comps));
3628:     if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3629:   }
3630:   if (bcFunc) b->func = bcFunc;
3631:   if (bcFunc_t) b->func_t = bcFunc_t;
3632:   if (ctx) b->ctx = ctx;
3633:   PetscFunctionReturn(PETSC_SUCCESS);
3634: }

3636: /*@
3637:   PetscDSGetNumBoundary - Get the number of registered BC

3639:   Input Parameter:
3640: . ds - The `PetscDS` object

3642:   Output Parameter:
3643: . numBd - The number of BC

3645:   Level: intermediate

3647: .seealso: `PetscDS`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`
3648: @*/
3649: PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd)
3650: {
3651:   DSBoundary b = ds->boundary;

3653:   PetscFunctionBegin;
3655:   PetscAssertPointer(numBd, 2);
3656:   *numBd = 0;
3657:   while (b) {
3658:     ++(*numBd);
3659:     b = b->next;
3660:   }
3661:   PetscFunctionReturn(PETSC_SUCCESS);
3662: }

3664: /*@C
3665:   PetscDSGetBoundary - Gets a boundary condition to the model

3667:   Input Parameters:
3668: + ds - The `PetscDS` object
3669: - bd - The BC number

3671:   Output Parameters:
3672: + wf     - The `PetscWeakForm` holding the pointwise functions
3673: . type   - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3674: . name   - The BC name
3675: . label  - The label defining constrained points
3676: . Nv     - The number of `DMLabel` ids for constrained points
3677: . values - An array of ids for constrained points
3678: . field  - The field to constrain
3679: . Nc     - The number of constrained field components
3680: . comps  - An array of constrained component numbers
3681: . func   - A pointwise function giving boundary values
3682: . func_t - A pointwise function giving the time derivative of the boundary values
3683: - ctx    - An optional user context for bcFunc

3685:   Options Database Keys:
3686: + -bc_<boundary name> <num>      - Overrides the boundary ids
3687: - -bc_<boundary name>_comp <num> - Overrides the boundary components

3689:   Level: developer

3691: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `DMLabel`
3692: @*/
3693: 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)
3694: {
3695:   DSBoundary b = ds->boundary;
3696:   PetscInt   n = 0;

3698:   PetscFunctionBegin;
3700:   while (b) {
3701:     if (n == bd) break;
3702:     b = b->next;
3703:     ++n;
3704:   }
3705:   PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n);
3706:   if (wf) {
3707:     PetscAssertPointer(wf, 3);
3708:     *wf = b->wf;
3709:   }
3710:   if (type) {
3711:     PetscAssertPointer(type, 4);
3712:     *type = b->type;
3713:   }
3714:   if (name) {
3715:     PetscAssertPointer(name, 5);
3716:     *name = b->name;
3717:   }
3718:   if (label) {
3719:     PetscAssertPointer(label, 6);
3720:     *label = b->label;
3721:   }
3722:   if (Nv) {
3723:     PetscAssertPointer(Nv, 7);
3724:     *Nv = b->Nv;
3725:   }
3726:   if (values) {
3727:     PetscAssertPointer(values, 8);
3728:     *values = b->values;
3729:   }
3730:   if (field) {
3731:     PetscAssertPointer(field, 9);
3732:     *field = b->field;
3733:   }
3734:   if (Nc) {
3735:     PetscAssertPointer(Nc, 10);
3736:     *Nc = b->Nc;
3737:   }
3738:   if (comps) {
3739:     PetscAssertPointer(comps, 11);
3740:     *comps = b->comps;
3741:   }
3742:   if (func) {
3743:     PetscAssertPointer(func, 12);
3744:     *func = b->func;
3745:   }
3746:   if (func_t) {
3747:     PetscAssertPointer(func_t, 13);
3748:     *func_t = b->func_t;
3749:   }
3750:   if (ctx) {
3751:     PetscAssertPointer(ctx, 14);
3752:     *ctx = b->ctx;
3753:   }
3754:   PetscFunctionReturn(PETSC_SUCCESS);
3755: }

3757: static PetscErrorCode DSBoundaryDuplicate_Internal(DSBoundary b, DSBoundary *bNew)
3758: {
3759:   PetscFunctionBegin;
3760:   PetscCall(PetscNew(bNew));
3761:   PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &(*bNew)->wf));
3762:   PetscCall(PetscWeakFormCopy(b->wf, (*bNew)->wf));
3763:   PetscCall(PetscStrallocpy(b->name, (char **)&((*bNew)->name)));
3764:   PetscCall(PetscStrallocpy(b->lname, (char **)&((*bNew)->lname)));
3765:   (*bNew)->type  = b->type;
3766:   (*bNew)->label = b->label;
3767:   (*bNew)->Nv    = b->Nv;
3768:   PetscCall(PetscMalloc1(b->Nv, &(*bNew)->values));
3769:   PetscCall(PetscArraycpy((*bNew)->values, b->values, b->Nv));
3770:   (*bNew)->field = b->field;
3771:   (*bNew)->Nc    = b->Nc;
3772:   PetscCall(PetscMalloc1(b->Nc, &(*bNew)->comps));
3773:   PetscCall(PetscArraycpy((*bNew)->comps, b->comps, b->Nc));
3774:   (*bNew)->func   = b->func;
3775:   (*bNew)->func_t = b->func_t;
3776:   (*bNew)->ctx    = b->ctx;
3777:   PetscFunctionReturn(PETSC_SUCCESS);
3778: }

3780: /*@
3781:   PetscDSCopyBoundary - Copy all boundary condition objects to the new problem

3783:   Not Collective

3785:   Input Parameters:
3786: + ds        - The source `PetscDS` object
3787: . numFields - The number of selected fields, or `PETSC_DEFAULT` for all fields
3788: - fields    - The selected fields, or NULL for all fields

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

3793:   Level: intermediate

3795: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3796: @*/
3797: PetscErrorCode PetscDSCopyBoundary(PetscDS ds, PetscInt numFields, const PetscInt fields[], PetscDS newds)
3798: {
3799:   DSBoundary b, *lastnext;

3801:   PetscFunctionBegin;
3804:   if (ds == newds) PetscFunctionReturn(PETSC_SUCCESS);
3805:   PetscCall(PetscDSDestroyBoundary(newds));
3806:   lastnext = &(newds->boundary);
3807:   for (b = ds->boundary; b; b = b->next) {
3808:     DSBoundary bNew;
3809:     PetscInt   fieldNew = -1;

3811:     if (numFields > 0 && fields) {
3812:       PetscInt f;

3814:       for (f = 0; f < numFields; ++f)
3815:         if (b->field == fields[f]) break;
3816:       if (f == numFields) continue;
3817:       fieldNew = f;
3818:     }
3819:     PetscCall(DSBoundaryDuplicate_Internal(b, &bNew));
3820:     bNew->field = fieldNew < 0 ? b->field : fieldNew;
3821:     *lastnext   = bNew;
3822:     lastnext    = &(bNew->next);
3823:   }
3824:   PetscFunctionReturn(PETSC_SUCCESS);
3825: }

3827: /*@
3828:   PetscDSDestroyBoundary - Remove all `DMBoundary` objects from the `PetscDS`

3830:   Not Collective

3832:   Input Parameter:
3833: . ds - The `PetscDS` object

3835:   Level: intermediate

3837: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`
3838: @*/
3839: PetscErrorCode PetscDSDestroyBoundary(PetscDS ds)
3840: {
3841:   DSBoundary next = ds->boundary;

3843:   PetscFunctionBegin;
3844:   while (next) {
3845:     DSBoundary b = next;

3847:     next = b->next;
3848:     PetscCall(PetscWeakFormDestroy(&b->wf));
3849:     PetscCall(PetscFree(b->name));
3850:     PetscCall(PetscFree(b->lname));
3851:     PetscCall(PetscFree(b->values));
3852:     PetscCall(PetscFree(b->comps));
3853:     PetscCall(PetscFree(b));
3854:   }
3855:   PetscFunctionReturn(PETSC_SUCCESS);
3856: }

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

3861:   Not Collective

3863:   Input Parameters:
3864: + prob      - The `PetscDS` object
3865: . numFields - Number of new fields
3866: - fields    - Old field number for each new field

3868:   Output Parameter:
3869: . newprob - The `PetscDS` copy

3871:   Level: intermediate

3873: .seealso: `PetscDS`, `PetscDSSelectEquations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3874: @*/
3875: PetscErrorCode PetscDSSelectDiscretizations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
3876: {
3877:   PetscInt Nf, Nfn, fn;

3879:   PetscFunctionBegin;
3881:   if (fields) PetscAssertPointer(fields, 3);
3883:   PetscCall(PetscDSGetNumFields(prob, &Nf));
3884:   PetscCall(PetscDSGetNumFields(newprob, &Nfn));
3885:   numFields = numFields < 0 ? Nf : numFields;
3886:   for (fn = 0; fn < numFields; ++fn) {
3887:     const PetscInt f = fields ? fields[fn] : fn;
3888:     PetscObject    disc;

3890:     if (f >= Nf) continue;
3891:     PetscCall(PetscDSGetDiscretization(prob, f, &disc));
3892:     PetscCall(PetscDSSetDiscretization(newprob, fn, disc));
3893:   }
3894:   PetscFunctionReturn(PETSC_SUCCESS);
3895: }

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

3900:   Not Collective

3902:   Input Parameters:
3903: + prob      - The `PetscDS` object
3904: . numFields - Number of new fields
3905: - fields    - Old field number for each new field

3907:   Output Parameter:
3908: . newprob - The `PetscDS` copy

3910:   Level: intermediate

3912: .seealso: `PetscDS`, `PetscDSSelectDiscretizations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3913: @*/
3914: PetscErrorCode PetscDSSelectEquations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
3915: {
3916:   PetscInt Nf, Nfn, fn, gn;

3918:   PetscFunctionBegin;
3920:   if (fields) PetscAssertPointer(fields, 3);
3922:   PetscCall(PetscDSGetNumFields(prob, &Nf));
3923:   PetscCall(PetscDSGetNumFields(newprob, &Nfn));
3924:   PetscCheck(numFields <= Nfn, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_SIZ, "Number of fields %" PetscInt_FMT " to transfer must not be greater then the total number of fields %" PetscInt_FMT, numFields, Nfn);
3925:   for (fn = 0; fn < numFields; ++fn) {
3926:     const PetscInt   f = fields ? fields[fn] : fn;
3927:     PetscPointFunc   obj;
3928:     PetscPointFunc   f0, f1;
3929:     PetscBdPointFunc f0Bd, f1Bd;
3930:     PetscRiemannFunc r;

3932:     if (f >= Nf) continue;
3933:     PetscCall(PetscDSGetObjective(prob, f, &obj));
3934:     PetscCall(PetscDSGetResidual(prob, f, &f0, &f1));
3935:     PetscCall(PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd));
3936:     PetscCall(PetscDSGetRiemannSolver(prob, f, &r));
3937:     PetscCall(PetscDSSetObjective(newprob, fn, obj));
3938:     PetscCall(PetscDSSetResidual(newprob, fn, f0, f1));
3939:     PetscCall(PetscDSSetBdResidual(newprob, fn, f0Bd, f1Bd));
3940:     PetscCall(PetscDSSetRiemannSolver(newprob, fn, r));
3941:     for (gn = 0; gn < numFields; ++gn) {
3942:       const PetscInt  g = fields ? fields[gn] : gn;
3943:       PetscPointJac   g0, g1, g2, g3;
3944:       PetscPointJac   g0p, g1p, g2p, g3p;
3945:       PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd;

3947:       if (g >= Nf) continue;
3948:       PetscCall(PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3));
3949:       PetscCall(PetscDSGetJacobianPreconditioner(prob, f, g, &g0p, &g1p, &g2p, &g3p));
3950:       PetscCall(PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd));
3951:       PetscCall(PetscDSSetJacobian(newprob, fn, gn, g0, g1, g2, g3));
3952:       PetscCall(PetscDSSetJacobianPreconditioner(newprob, fn, gn, g0p, g1p, g2p, g3p));
3953:       PetscCall(PetscDSSetBdJacobian(newprob, fn, gn, g0Bd, g1Bd, g2Bd, g3Bd));
3954:     }
3955:   }
3956:   PetscFunctionReturn(PETSC_SUCCESS);
3957: }

3959: /*@
3960:   PetscDSCopyEquations - Copy all pointwise function pointers to another `PetscDS`

3962:   Not Collective

3964:   Input Parameter:
3965: . prob - The `PetscDS` object

3967:   Output Parameter:
3968: . newprob - The `PetscDS` copy

3970:   Level: intermediate

3972: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3973: @*/
3974: PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob)
3975: {
3976:   PetscWeakForm wf, newwf;
3977:   PetscInt      Nf, Ng;

3979:   PetscFunctionBegin;
3982:   PetscCall(PetscDSGetNumFields(prob, &Nf));
3983:   PetscCall(PetscDSGetNumFields(newprob, &Ng));
3984:   PetscCheck(Nf == Ng, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %" PetscInt_FMT " != %" PetscInt_FMT, Nf, Ng);
3985:   PetscCall(PetscDSGetWeakForm(prob, &wf));
3986:   PetscCall(PetscDSGetWeakForm(newprob, &newwf));
3987:   PetscCall(PetscWeakFormCopy(wf, newwf));
3988:   PetscFunctionReturn(PETSC_SUCCESS);
3989: }

3991: /*@
3992:   PetscDSCopyConstants - Copy all constants to another `PetscDS`

3994:   Not Collective

3996:   Input Parameter:
3997: . prob - The `PetscDS` object

3999:   Output Parameter:
4000: . newprob - The `PetscDS` copy

4002:   Level: intermediate

4004: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
4005: @*/
4006: PetscErrorCode PetscDSCopyConstants(PetscDS prob, PetscDS newprob)
4007: {
4008:   PetscInt           Nc;
4009:   const PetscScalar *constants;

4011:   PetscFunctionBegin;
4014:   PetscCall(PetscDSGetConstants(prob, &Nc, &constants));
4015:   PetscCall(PetscDSSetConstants(newprob, Nc, (PetscScalar *)constants));
4016:   PetscFunctionReturn(PETSC_SUCCESS);
4017: }

4019: /*@
4020:   PetscDSCopyExactSolutions - Copy all exact solutions to another `PetscDS`

4022:   Not Collective

4024:   Input Parameter:
4025: . ds - The `PetscDS` object

4027:   Output Parameter:
4028: . newds - The `PetscDS` copy

4030:   Level: intermediate

4032: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
4033: @*/
4034: PetscErrorCode PetscDSCopyExactSolutions(PetscDS ds, PetscDS newds)
4035: {
4036:   PetscSimplePointFunc sol;
4037:   void                *ctx;
4038:   PetscInt             Nf, f;

4040:   PetscFunctionBegin;
4043:   PetscCall(PetscDSGetNumFields(ds, &Nf));
4044:   for (f = 0; f < Nf; ++f) {
4045:     PetscCall(PetscDSGetExactSolution(ds, f, &sol, &ctx));
4046:     PetscCall(PetscDSSetExactSolution(newds, f, sol, ctx));
4047:     PetscCall(PetscDSGetExactSolutionTimeDerivative(ds, f, &sol, &ctx));
4048:     PetscCall(PetscDSSetExactSolutionTimeDerivative(newds, f, sol, ctx));
4049:   }
4050:   PetscFunctionReturn(PETSC_SUCCESS);
4051: }

4053: PetscErrorCode PetscDSCopy(PetscDS ds, DM dmNew, PetscDS dsNew)
4054: {
4055:   DSBoundary b;
4056:   PetscInt   cdim, Nf, f, d;
4057:   PetscBool  isCohesive;
4058:   void      *ctx;

4060:   PetscFunctionBegin;
4061:   PetscCall(PetscDSCopyConstants(ds, dsNew));
4062:   PetscCall(PetscDSCopyExactSolutions(ds, dsNew));
4063:   PetscCall(PetscDSSelectDiscretizations(ds, PETSC_DETERMINE, NULL, dsNew));
4064:   PetscCall(PetscDSCopyEquations(ds, dsNew));
4065:   PetscCall(PetscDSGetNumFields(ds, &Nf));
4066:   for (f = 0; f < Nf; ++f) {
4067:     PetscCall(PetscDSGetContext(ds, f, &ctx));
4068:     PetscCall(PetscDSSetContext(dsNew, f, ctx));
4069:     PetscCall(PetscDSGetCohesive(ds, f, &isCohesive));
4070:     PetscCall(PetscDSSetCohesive(dsNew, f, isCohesive));
4071:     PetscCall(PetscDSGetJetDegree(ds, f, &d));
4072:     PetscCall(PetscDSSetJetDegree(dsNew, f, d));
4073:   }
4074:   if (Nf) {
4075:     PetscCall(PetscDSGetCoordinateDimension(ds, &cdim));
4076:     PetscCall(PetscDSSetCoordinateDimension(dsNew, cdim));
4077:   }
4078:   PetscCall(PetscDSCopyBoundary(ds, PETSC_DETERMINE, NULL, dsNew));
4079:   for (b = dsNew->boundary; b; b = b->next) {
4080:     PetscCall(DMGetLabel(dmNew, b->lname, &b->label));
4081:     /* Do not check if label exists here, since p4est calls this for the reference tree which does not have the labels */
4082:     //PetscCheck(b->label,PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Label %s missing in new DM", name);
4083:   }
4084:   PetscFunctionReturn(PETSC_SUCCESS);
4085: }

4087: PetscErrorCode PetscDSGetHeightSubspace(PetscDS prob, PetscInt height, PetscDS *subprob)
4088: {
4089:   PetscInt dim, Nf, f;

4091:   PetscFunctionBegin;
4093:   PetscAssertPointer(subprob, 3);
4094:   if (height == 0) {
4095:     *subprob = prob;
4096:     PetscFunctionReturn(PETSC_SUCCESS);
4097:   }
4098:   PetscCall(PetscDSGetNumFields(prob, &Nf));
4099:   PetscCall(PetscDSGetSpatialDimension(prob, &dim));
4100:   PetscCheck(height <= dim, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_OUTOFRANGE, "DS can only handle height in [0, %" PetscInt_FMT "], not %" PetscInt_FMT, dim, height);
4101:   if (!prob->subprobs) PetscCall(PetscCalloc1(dim, &prob->subprobs));
4102:   if (!prob->subprobs[height - 1]) {
4103:     PetscInt cdim;

4105:     PetscCall(PetscDSCreate(PetscObjectComm((PetscObject)prob), &prob->subprobs[height - 1]));
4106:     PetscCall(PetscDSGetCoordinateDimension(prob, &cdim));
4107:     PetscCall(PetscDSSetCoordinateDimension(prob->subprobs[height - 1], cdim));
4108:     for (f = 0; f < Nf; ++f) {
4109:       PetscFE      subfe;
4110:       PetscObject  obj;
4111:       PetscClassId id;

4113:       PetscCall(PetscDSGetDiscretization(prob, f, &obj));
4114:       PetscCall(PetscObjectGetClassId(obj, &id));
4115:       if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetHeightSubspace((PetscFE)obj, height, &subfe));
4116:       else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unsupported discretization type for field %" PetscInt_FMT, f);
4117:       PetscCall(PetscDSSetDiscretization(prob->subprobs[height - 1], f, (PetscObject)subfe));
4118:     }
4119:   }
4120:   *subprob = prob->subprobs[height - 1];
4121:   PetscFunctionReturn(PETSC_SUCCESS);
4122: }

4124: PetscErrorCode PetscDSPermuteQuadPoint(PetscDS ds, PetscInt ornt, PetscInt field, PetscInt q, PetscInt *qperm)
4125: {
4126:   IS              permIS;
4127:   PetscQuadrature quad;
4128:   DMPolytopeType  ct;
4129:   const PetscInt *perm;
4130:   PetscInt        Na, Nq;

4132:   PetscFunctionBeginHot;
4133:   PetscCall(PetscFEGetQuadrature((PetscFE)ds->disc[field], &quad));
4134:   PetscCall(PetscQuadratureGetData(quad, NULL, NULL, &Nq, NULL, NULL));
4135:   PetscCall(PetscQuadratureGetCellType(quad, &ct));
4136:   PetscCheck(q >= 0 && q < Nq, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Quadrature point %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", q, Nq);
4137:   Na = DMPolytopeTypeGetNumArrangments(ct) / 2;
4138:   PetscCheck(ornt >= -Na && ornt < Na, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Orientation %" PetscInt_FMT " of %s is not in [%" PetscInt_FMT ", %" PetscInt_FMT ")", ornt, DMPolytopeTypes[ct], -Na, Na);
4139:   if (!ds->quadPerm[(PetscInt)ct]) PetscCall(PetscQuadratureComputePermutations(quad, NULL, &ds->quadPerm[(PetscInt)ct]));
4140:   permIS = ds->quadPerm[(PetscInt)ct][ornt + Na];
4141:   PetscCall(ISGetIndices(permIS, &perm));
4142:   *qperm = perm[q];
4143:   PetscCall(ISRestoreIndices(permIS, &perm));
4144:   PetscFunctionReturn(PETSC_SUCCESS);
4145: }

4147: PetscErrorCode PetscDSGetDiscType_Internal(PetscDS ds, PetscInt f, PetscDiscType *disctype)
4148: {
4149:   PetscObject  obj;
4150:   PetscClassId id;
4151:   PetscInt     Nf;

4153:   PetscFunctionBegin;
4155:   PetscAssertPointer(disctype, 3);
4156:   *disctype = PETSC_DISC_NONE;
4157:   PetscCall(PetscDSGetNumFields(ds, &Nf));
4158:   PetscCheck(f < Nf, PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_SIZ, "Field %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, Nf);
4159:   PetscCall(PetscDSGetDiscretization(ds, f, &obj));
4160:   if (obj) {
4161:     PetscCall(PetscObjectGetClassId(obj, &id));
4162:     if (id == PETSCFE_CLASSID) *disctype = PETSC_DISC_FE;
4163:     else *disctype = PETSC_DISC_FV;
4164:   }
4165:   PetscFunctionReturn(PETSC_SUCCESS);
4166: }

4168: static PetscErrorCode PetscDSDestroy_Basic(PetscDS ds)
4169: {
4170:   PetscFunctionBegin;
4171:   PetscCall(PetscFree(ds->data));
4172:   PetscFunctionReturn(PETSC_SUCCESS);
4173: }

4175: static PetscErrorCode PetscDSInitialize_Basic(PetscDS ds)
4176: {
4177:   PetscFunctionBegin;
4178:   ds->ops->setfromoptions = NULL;
4179:   ds->ops->setup          = NULL;
4180:   ds->ops->view           = NULL;
4181:   ds->ops->destroy        = PetscDSDestroy_Basic;
4182:   PetscFunctionReturn(PETSC_SUCCESS);
4183: }

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

4188:   Level: intermediate

4190: .seealso: `PetscDSType`, `PetscDSCreate()`, `PetscDSSetType()`
4191: M*/

4193: PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS ds)
4194: {
4195:   PetscDS_Basic *b;

4197:   PetscFunctionBegin;
4199:   PetscCall(PetscNew(&b));
4200:   ds->data = b;

4202:   PetscCall(PetscDSInitialize_Basic(ds));
4203:   PetscFunctionReturn(PETSC_SUCCESS);
4204: }