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: /*@
9: PetscDSRegister - Adds a new `PetscDS` implementation
11: Not Collective; No Fortran Support
13: Input Parameters:
14: + sname - The name of a new user-defined creation routine
15: - function - The creation routine itself
17: Example Usage:
18: .vb
19: PetscDSRegister("my_ds", MyPetscDSCreate);
20: .ve
22: Then, your PetscDS type can be chosen with the procedural interface via
23: .vb
24: PetscDSCreate(MPI_Comm, PetscDS *);
25: PetscDSSetType(PetscDS, "my_ds");
26: .ve
27: or at runtime via the option
28: .vb
29: -petscds_type my_ds
30: .ve
32: Level: advanced
34: Note:
35: `PetscDSRegister()` may be called multiple times to add several user-defined `PetscDSs`
37: .seealso: `PetscDSType`, `PetscDS`, `PetscDSRegisterAll()`
38: @*/
39: PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS))
40: {
41: PetscFunctionBegin;
42: PetscCall(PetscFunctionListAdd(&PetscDSList, sname, function));
43: PetscFunctionReturn(PETSC_SUCCESS);
44: }
46: /*@
47: PetscDSSetType - Builds a particular `PetscDS`
49: Collective; No Fortran Support
51: Input Parameters:
52: + prob - The `PetscDS` object
53: - name - The `PetscDSType`
55: Options Database Key:
56: . -petscds_type type - Sets the PetscDS type; use -help for a list of available types
58: Level: intermediate
60: .seealso: `PetscDSType`, `PetscDS`, `PetscDSGetType()`, `PetscDSCreate()`
61: @*/
62: PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name)
63: {
64: PetscErrorCode (*r)(PetscDS);
65: PetscBool match;
67: PetscFunctionBegin;
69: PetscCall(PetscObjectTypeCompare((PetscObject)prob, name, &match));
70: if (match) PetscFunctionReturn(PETSC_SUCCESS);
72: PetscCall(PetscDSRegisterAll());
73: PetscCall(PetscFunctionListFind(PetscDSList, name, &r));
74: PetscCheck(r, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscDS type: %s", name);
76: PetscTryTypeMethod(prob, destroy);
77: prob->ops->destroy = NULL;
79: PetscCall((*r)(prob));
80: PetscCall(PetscObjectChangeTypeName((PetscObject)prob, name));
81: PetscFunctionReturn(PETSC_SUCCESS);
82: }
84: /*@
85: PetscDSGetType - Gets the `PetscDSType` name (as a string) from the `PetscDS`
87: Not Collective; No Fortran Support
89: Input Parameter:
90: . prob - The `PetscDS`
92: Output Parameter:
93: . name - The `PetscDSType` name
95: Level: intermediate
97: .seealso: `PetscDSType`, `PetscDS`, `PetscDSSetType()`, `PetscDSCreate()`
98: @*/
99: PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name)
100: {
101: PetscFunctionBegin;
103: PetscAssertPointer(name, 2);
104: PetscCall(PetscDSRegisterAll());
105: *name = ((PetscObject)prob)->type_name;
106: PetscFunctionReturn(PETSC_SUCCESS);
107: }
109: static PetscErrorCode PetscDSView_Ascii(PetscDS ds, PetscViewer viewer)
110: {
111: PetscViewerFormat format;
112: const PetscScalar *constants;
113: PetscInt Nf, numConstants, f;
115: PetscFunctionBegin;
116: PetscCall(PetscDSGetNumFields(ds, &Nf));
117: PetscCall(PetscViewerGetFormat(viewer, &format));
118: PetscCall(PetscViewerASCIIPrintf(viewer, "Discrete System with %" PetscInt_FMT " fields\n", Nf));
119: PetscCall(PetscViewerASCIIPushTab(viewer));
120: PetscCall(PetscViewerASCIIPrintf(viewer, " cell total dim %" PetscInt_FMT " total comp %" PetscInt_FMT "\n", ds->totDim, ds->totComp));
121: if (ds->isCohesive) PetscCall(PetscViewerASCIIPrintf(viewer, " cohesive cell\n"));
122: for (f = 0; f < Nf; ++f) {
123: DSBoundary b;
124: PetscObject obj;
125: PetscClassId id;
126: PetscQuadrature q;
127: const char *name;
128: PetscInt Nc, Nq, Nqc;
130: PetscCall(PetscDSGetDiscretization(ds, f, &obj));
131: PetscCall(PetscObjectGetClassId(obj, &id));
132: PetscCall(PetscObjectGetName(obj, &name));
133: PetscCall(PetscViewerASCIIPrintf(viewer, "Field %s", name ? name : "<unknown>"));
134: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
135: if (id == PETSCFE_CLASSID) {
136: PetscCall(PetscFEGetNumComponents((PetscFE)obj, &Nc));
137: PetscCall(PetscFEGetQuadrature((PetscFE)obj, &q));
138: PetscCall(PetscViewerASCIIPrintf(viewer, " FEM"));
139: } else if (id == PETSCFV_CLASSID) {
140: PetscCall(PetscFVGetNumComponents((PetscFV)obj, &Nc));
141: PetscCall(PetscFVGetQuadrature((PetscFV)obj, &q));
142: PetscCall(PetscViewerASCIIPrintf(viewer, " FVM"));
143: } else SETERRQ(PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f);
144: if (Nc > 1) PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " components", Nc));
145: else PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " component ", Nc));
146: if (ds->implicit[f]) PetscCall(PetscViewerASCIIPrintf(viewer, " (implicit)"));
147: else PetscCall(PetscViewerASCIIPrintf(viewer, " (explicit)"));
148: if (q) {
149: PetscCall(PetscQuadratureGetData(q, NULL, &Nqc, &Nq, NULL, NULL));
150: PetscCall(PetscViewerASCIIPrintf(viewer, " (Nq %" PetscInt_FMT " Nqc %" PetscInt_FMT ")", Nq, Nqc));
151: }
152: PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT "-jet", ds->jetDegree[f]));
153: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
154: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
155: PetscCall(PetscViewerASCIIPushTab(viewer));
156: if (id == PETSCFE_CLASSID) PetscCall(PetscFEView((PetscFE)obj, viewer));
157: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVView((PetscFV)obj, viewer));
158: PetscCall(PetscViewerASCIIPopTab(viewer));
160: for (b = ds->boundary; b; b = b->next) {
161: char *name;
163: if (b->field != f) continue;
164: PetscCall(PetscViewerASCIIPushTab(viewer));
165: PetscCall(PetscViewerASCIIPrintf(viewer, "Boundary %s (%s) %s\n", b->name, b->lname, DMBoundaryConditionTypes[b->type]));
166: if (!b->Nc) {
167: PetscCall(PetscViewerASCIIPrintf(viewer, " all components\n"));
168: } else {
169: PetscCall(PetscViewerASCIIPrintf(viewer, " components: "));
170: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
171: for (PetscInt c = 0; c < b->Nc; ++c) {
172: if (c > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", "));
173: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->comps[c]));
174: }
175: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
176: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
177: }
178: PetscCall(PetscViewerASCIIPrintf(viewer, " values: "));
179: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
180: for (PetscInt i = 0; i < b->Nv; ++i) {
181: if (i > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", "));
182: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->values[i]));
183: }
184: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
185: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
186: #if defined(__clang__)
187: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat-pedantic")
188: #elif defined(__GNUC__) || defined(__GNUG__)
189: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat")
190: #endif
191: if (b->func) {
192: PetscCall(PetscDLAddr(b->func, &name));
193: if (name) PetscCall(PetscViewerASCIIPrintf(viewer, " func: %s\n", name));
194: else PetscCall(PetscViewerASCIIPrintf(viewer, " func: %p\n", b->func));
195: PetscCall(PetscFree(name));
196: }
197: if (b->func_t) {
198: PetscCall(PetscDLAddr(b->func_t, &name));
199: if (name) PetscCall(PetscViewerASCIIPrintf(viewer, " func_t: %s\n", name));
200: else PetscCall(PetscViewerASCIIPrintf(viewer, " func_t: %p\n", b->func_t));
201: PetscCall(PetscFree(name));
202: }
203: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_END()
204: PetscCall(PetscWeakFormView(b->wf, viewer));
205: PetscCall(PetscViewerASCIIPopTab(viewer));
206: }
207: }
208: PetscCall(PetscDSGetConstants(ds, &numConstants, &constants));
209: if (numConstants) {
210: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " constants\n", numConstants));
211: PetscCall(PetscViewerASCIIPushTab(viewer));
212: for (f = 0; f < numConstants; ++f) PetscCall(PetscViewerASCIIPrintf(viewer, "%g\n", (double)PetscRealPart(constants[f])));
213: PetscCall(PetscViewerASCIIPopTab(viewer));
214: }
215: PetscCall(PetscWeakFormView(ds->wf, viewer));
216: PetscCall(PetscViewerASCIIPopTab(viewer));
217: PetscFunctionReturn(PETSC_SUCCESS);
218: }
220: /*@
221: PetscDSViewFromOptions - View a `PetscDS` based on values in the options database
223: Collective
225: Input Parameters:
226: + A - the `PetscDS` object
227: . obj - optional object that provides the options prefix used in the search of the options database, pass `NULL` to use the options prefix of `A`
228: - name - command line option
230: Options Database Key:
231: . -name viewer_specification - See `PetscOptionsCreateViewer()` for the values of `viewer_specification`
233: Level: intermediate
235: Note:
236: This checks the options database, creates the viewer on-the-fly, uses it and then destroys it. Hence it should not be called in heavily used routines,
237: rather `PetscOptionsCreateViewer()` should be used to construct the viewer once which can then be utilized in the heavily used routine.
239: .seealso: `PetscDSType`, `PetscDS`, `PetscDSView()`, `PetscObjectViewFromOptions()`, `PetscDSCreate()`, `PetscOptionsCreateViewer()`
240: @*/
241: PetscErrorCode PetscDSViewFromOptions(PetscDS A, PetscObject obj, const char name[])
242: {
243: PetscFunctionBegin;
245: PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
246: PetscFunctionReturn(PETSC_SUCCESS);
247: }
249: /*@
250: PetscDSView - Views a `PetscDS`
252: Collective
254: Input Parameters:
255: + prob - the `PetscDS` object to view
256: - v - the viewer
258: Level: developer
260: .seealso: `PetscDSType`, `PetscDS`, `PetscViewer`, `PetscDSDestroy()`, `PetscDSViewFromOptions()`
261: @*/
262: PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v)
263: {
264: PetscBool isascii;
266: PetscFunctionBegin;
268: if (!v) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)prob), &v));
270: PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &isascii));
271: if (isascii) PetscCall(PetscDSView_Ascii(prob, v));
272: PetscTryTypeMethod(prob, view, v);
273: PetscFunctionReturn(PETSC_SUCCESS);
274: }
276: /*@
277: PetscDSSetFromOptions - sets parameters in a `PetscDS` from the options database
279: Collective
281: Input Parameter:
282: . prob - the `PetscDS` object to set options for
284: Options Database Keys:
285: + -petscds_type type - set the `PetscDS` type
286: . -petscds_view viewer_specification - view the `PetscDS` at the end of this call, see `PetscOptionsCreateViewer()` for the format of `viewer_specification`
287: . -petscds_jac_pre (true|false) - turn formation of a separate Jacobian preconditioner on or off
288: . -bc_NAME ids - comma separated list of label ids for the boundary condition NAME
289: - -bc_NAME_comp comps - comma separated list of field components to constrain for the boundary condition NAME
291: Level: intermediate
293: .seealso: `PetscDS`, `PetscDSView()`, `PetscOptionsCreateViewer()`
294: @*/
295: PetscErrorCode PetscDSSetFromOptions(PetscDS prob)
296: {
297: DSBoundary b;
298: const char *defaultType;
299: char name[256];
300: PetscBool flg;
302: PetscFunctionBegin;
304: if (!((PetscObject)prob)->type_name) {
305: defaultType = PETSCDSBASIC;
306: } else {
307: defaultType = ((PetscObject)prob)->type_name;
308: }
309: PetscCall(PetscDSRegisterAll());
311: PetscObjectOptionsBegin((PetscObject)prob);
312: for (b = prob->boundary; b; b = b->next) {
313: char optname[1024];
314: PetscInt ids[1024], len = 1024;
315: PetscBool flg;
317: PetscCall(PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name));
318: PetscCall(PetscMemzero(ids, sizeof(ids)));
319: PetscCall(PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg));
320: if (flg) {
321: b->Nv = len;
322: PetscCall(PetscFree(b->values));
323: PetscCall(PetscMalloc1(len, &b->values));
324: PetscCall(PetscArraycpy(b->values, ids, len));
325: PetscCall(PetscWeakFormRewriteKeys(b->wf, b->label, len, b->values));
326: }
327: len = 1024;
328: PetscCall(PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name));
329: PetscCall(PetscMemzero(ids, sizeof(ids)));
330: PetscCall(PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg));
331: if (flg) {
332: b->Nc = len;
333: PetscCall(PetscFree(b->comps));
334: PetscCall(PetscMalloc1(len, &b->comps));
335: PetscCall(PetscArraycpy(b->comps, ids, len));
336: }
337: }
338: PetscCall(PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, sizeof(name), &flg));
339: if (flg) {
340: PetscCall(PetscDSSetType(prob, name));
341: } else if (!((PetscObject)prob)->type_name) {
342: PetscCall(PetscDSSetType(prob, defaultType));
343: }
344: PetscCall(PetscOptionsBool("-petscds_jac_pre", "Discrete System", "PetscDSUseJacobianPreconditioner", prob->useJacPre, &prob->useJacPre, &flg));
345: PetscCall(PetscOptionsBool("-petscds_force_quad", "Discrete System", "PetscDSSetForceQuad", prob->forceQuad, &prob->forceQuad, &flg));
346: PetscCall(PetscOptionsInt("-petscds_print_integrate", "Discrete System", "", prob->printIntegrate, &prob->printIntegrate, NULL));
347: PetscTryTypeMethod(prob, setfromoptions);
348: /* process any options handlers added with PetscObjectAddOptionsHandler() */
349: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)prob, PetscOptionsObject));
350: PetscOptionsEnd();
351: if (prob->Nf) PetscCall(PetscDSViewFromOptions(prob, NULL, "-petscds_view"));
352: PetscFunctionReturn(PETSC_SUCCESS);
353: }
355: /*@
356: PetscDSSetUp - Construct data structures for the `PetscDS`
358: Collective
360: Input Parameter:
361: . prob - the `PetscDS` object to setup
363: Level: developer
365: .seealso: `PetscDS`, `PetscDSView()`, `PetscDSDestroy()`
366: @*/
367: PetscErrorCode PetscDSSetUp(PetscDS prob)
368: {
369: const PetscInt Nf = prob->Nf;
370: PetscBool hasH = PETSC_FALSE;
371: PetscInt maxOrder[4] = {-2, -2, -2, -2};
372: PetscInt dim, dimEmbed, NbMax = 0, NcMax = 0, NqMax = 0, NsMax = 1, f;
374: PetscFunctionBegin;
376: if (prob->setup) PetscFunctionReturn(PETSC_SUCCESS);
377: /* Calculate sizes */
378: PetscCall(PetscDSGetSpatialDimension(prob, &dim));
379: PetscCall(PetscDSGetCoordinateDimension(prob, &dimEmbed));
380: prob->totDim = prob->totComp = 0;
381: PetscCall(PetscMalloc2(Nf, &prob->Nc, Nf, &prob->Nb));
382: PetscCall(PetscCalloc2(Nf + 1, &prob->off, Nf + 1, &prob->offDer));
383: 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]));
384: PetscCall(PetscMalloc2(Nf, &prob->T, Nf, &prob->Tf));
385: if (prob->forceQuad) {
386: // Note: This assumes we have one kind of cell at each dimension.
387: // We can fix this by having quadrature hold the celltype
388: PetscQuadrature maxQuad[4] = {NULL, NULL, NULL, NULL};
390: for (f = 0; f < Nf; ++f) {
391: PetscObject obj;
392: PetscClassId id;
393: PetscQuadrature q = NULL, fq = NULL;
394: PetscInt dim = -1, order = -1, forder = -1;
396: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
397: if (!obj) continue;
398: PetscCall(PetscObjectGetClassId(obj, &id));
399: if (id == PETSCFE_CLASSID) {
400: PetscFE fe = (PetscFE)obj;
402: PetscCall(PetscFEGetQuadrature(fe, &q));
403: PetscCall(PetscFEGetFaceQuadrature(fe, &fq));
404: } else if (id == PETSCFV_CLASSID) {
405: PetscFV fv = (PetscFV)obj;
407: PetscCall(PetscFVGetQuadrature(fv, &q));
408: }
409: if (q) {
410: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
411: PetscCall(PetscQuadratureGetOrder(q, &order));
412: if (order > maxOrder[dim]) {
413: maxOrder[dim] = order;
414: maxQuad[dim] = q;
415: }
416: }
417: if (fq) {
418: PetscCall(PetscQuadratureGetData(fq, &dim, NULL, NULL, NULL, NULL));
419: PetscCall(PetscQuadratureGetOrder(fq, &forder));
420: if (forder > maxOrder[dim]) {
421: maxOrder[dim] = forder;
422: maxQuad[dim] = fq;
423: }
424: }
425: }
426: for (f = 0; f < Nf; ++f) {
427: PetscObject obj;
428: PetscClassId id;
429: PetscQuadrature q;
430: PetscInt dim;
432: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
433: if (!obj) continue;
434: PetscCall(PetscObjectGetClassId(obj, &id));
435: if (id == PETSCFE_CLASSID) {
436: PetscFE fe = (PetscFE)obj;
438: PetscCall(PetscFEGetQuadrature(fe, &q));
439: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
440: PetscCall(PetscFESetQuadrature(fe, maxQuad[dim]));
441: PetscCall(PetscFESetFaceQuadrature(fe, dim ? maxQuad[dim - 1] : NULL));
442: } else if (id == PETSCFV_CLASSID) {
443: PetscFV fv = (PetscFV)obj;
445: PetscCall(PetscFVGetQuadrature(fv, &q));
446: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
447: PetscCall(PetscFVSetQuadrature(fv, maxQuad[dim]));
448: }
449: }
450: }
451: for (f = 0; f < Nf; ++f) {
452: PetscObject obj;
453: PetscClassId id;
454: PetscQuadrature q = NULL;
455: PetscInt Nq = 0, Nb, Nc;
457: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
458: if (prob->jetDegree[f] > 1) hasH = PETSC_TRUE;
459: if (!obj) {
460: /* Empty mesh */
461: Nb = Nc = 0;
462: prob->T[f] = prob->Tf[f] = NULL;
463: } else {
464: PetscCall(PetscObjectGetClassId(obj, &id));
465: if (id == PETSCFE_CLASSID) {
466: PetscFE fe = (PetscFE)obj;
468: PetscCall(PetscFEGetQuadrature(fe, &q));
469: {
470: PetscQuadrature fq;
471: PetscInt dim, order;
473: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
474: PetscCall(PetscQuadratureGetOrder(q, &order));
475: if (maxOrder[dim] < 0) maxOrder[dim] = order;
476: 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]);
477: PetscCall(PetscFEGetFaceQuadrature(fe, &fq));
478: if (fq) {
479: PetscCall(PetscQuadratureGetData(fq, &dim, NULL, NULL, NULL, NULL));
480: PetscCall(PetscQuadratureGetOrder(fq, &order));
481: if (maxOrder[dim] < 0) maxOrder[dim] = order;
482: 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]);
483: }
484: }
485: PetscCall(PetscFEGetDimension(fe, &Nb));
486: PetscCall(PetscFEGetNumComponents(fe, &Nc));
487: PetscCall(PetscFEGetCellTabulation(fe, prob->jetDegree[f], &prob->T[f]));
488: PetscCall(PetscFEGetFaceTabulation(fe, prob->jetDegree[f], &prob->Tf[f]));
489: } else if (id == PETSCFV_CLASSID) {
490: PetscFV fv = (PetscFV)obj;
492: PetscCall(PetscFVGetQuadrature(fv, &q));
493: PetscCall(PetscFVGetNumComponents(fv, &Nc));
494: Nb = Nc;
495: PetscCall(PetscFVGetCellTabulation(fv, &prob->T[f]));
496: /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */
497: } else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f);
498: }
499: prob->Nc[f] = Nc;
500: prob->Nb[f] = Nb;
501: prob->off[f + 1] = Nc + prob->off[f];
502: prob->offDer[f + 1] = Nc * dim + prob->offDer[f];
503: prob->offCohesive[0][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[0][f];
504: prob->offDerCohesive[0][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[0][f];
505: prob->offCohesive[1][f] = (prob->cohesive[f] ? 0 : Nc) + prob->offCohesive[0][f];
506: prob->offDerCohesive[1][f] = (prob->cohesive[f] ? 0 : Nc) * dimEmbed + prob->offDerCohesive[0][f];
507: prob->offCohesive[2][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[2][f];
508: prob->offDerCohesive[2][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[2][f];
509: if (q) PetscCall(PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL));
510: NqMax = PetscMax(NqMax, Nq);
511: NbMax = PetscMax(NbMax, Nb);
512: NcMax = PetscMax(NcMax, Nc);
513: prob->totDim += Nb;
514: prob->totComp += Nc;
515: /* There are two faces for all fields on a cohesive cell, except for cohesive fields */
516: if (prob->isCohesive && !prob->cohesive[f]) prob->totDim += Nb;
517: }
518: prob->offCohesive[1][Nf] = prob->offCohesive[0][Nf];
519: prob->offDerCohesive[1][Nf] = prob->offDerCohesive[0][Nf];
520: /* Allocate works space */
521: NsMax = 2; /* A non-cohesive discretizations can be used on a cohesive cell, so we need this extra workspace for all DS */
522: 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));
523: PetscCall(PetscMalloc5(dimEmbed, &prob->x, NbMax * NcMax, &prob->basisReal, NbMax * NcMax * dimEmbed, &prob->basisDerReal, NbMax * NcMax, &prob->testReal, NbMax * NcMax * dimEmbed, &prob->testDerReal));
524: 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,
525: &prob->g2, NsMax * NsMax * NqMax * NcMax * NcMax * dimEmbed * dimEmbed, &prob->g3));
526: PetscTryTypeMethod(prob, setup);
527: prob->setup = PETSC_TRUE;
528: PetscFunctionReturn(PETSC_SUCCESS);
529: }
531: static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob)
532: {
533: PetscFunctionBegin;
534: PetscCall(PetscFree2(prob->Nc, prob->Nb));
535: PetscCall(PetscFree2(prob->off, prob->offDer));
536: PetscCall(PetscFree6(prob->offCohesive[0], prob->offCohesive[1], prob->offCohesive[2], prob->offDerCohesive[0], prob->offDerCohesive[1], prob->offDerCohesive[2]));
537: PetscCall(PetscFree2(prob->T, prob->Tf));
538: PetscCall(PetscFree3(prob->u, prob->u_t, prob->u_x));
539: PetscCall(PetscFree5(prob->x, prob->basisReal, prob->basisDerReal, prob->testReal, prob->testDerReal));
540: PetscCall(PetscFree6(prob->f0, prob->f1, prob->g0, prob->g1, prob->g2, prob->g3));
541: PetscFunctionReturn(PETSC_SUCCESS);
542: }
544: static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew)
545: {
546: PetscObject *tmpd;
547: PetscBool *tmpi;
548: PetscInt *tmpk;
549: PetscBool *tmpc;
550: PetscPointFn **tmpup;
551: PetscSimplePointFn **tmpexactSol, **tmpexactSol_t, **tmplowerBound, **tmpupperBound;
552: void **tmpexactCtx, **tmpexactCtx_t, **tmplowerCtx, **tmpupperCtx;
553: void **tmpctx;
554: PetscInt Nf = prob->Nf, f;
556: PetscFunctionBegin;
557: if (Nf >= NfNew) PetscFunctionReturn(PETSC_SUCCESS);
558: prob->setup = PETSC_FALSE;
559: PetscCall(PetscDSDestroyStructs_Static(prob));
560: PetscCall(PetscMalloc4(NfNew, &tmpd, NfNew, &tmpi, NfNew, &tmpc, NfNew, &tmpk));
561: for (f = 0; f < Nf; ++f) {
562: tmpd[f] = prob->disc[f];
563: tmpi[f] = prob->implicit[f];
564: tmpc[f] = prob->cohesive[f];
565: tmpk[f] = prob->jetDegree[f];
566: }
567: for (f = Nf; f < NfNew; ++f) {
568: tmpd[f] = NULL;
569: tmpi[f] = PETSC_TRUE, tmpc[f] = PETSC_FALSE;
570: tmpk[f] = 1;
571: }
572: PetscCall(PetscFree4(prob->disc, prob->implicit, prob->cohesive, prob->jetDegree));
573: PetscCall(PetscWeakFormSetNumFields(prob->wf, NfNew));
574: prob->Nf = NfNew;
575: prob->disc = tmpd;
576: prob->implicit = tmpi;
577: prob->cohesive = tmpc;
578: prob->jetDegree = tmpk;
579: PetscCall(PetscCalloc2(NfNew, &tmpup, NfNew, &tmpctx));
580: for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f];
581: for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f];
582: for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL;
583: for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL;
584: PetscCall(PetscFree2(prob->update, prob->ctx));
585: prob->update = tmpup;
586: prob->ctx = tmpctx;
587: PetscCall(PetscCalloc4(NfNew, &tmpexactSol, NfNew, &tmpexactCtx, NfNew, &tmpexactSol_t, NfNew, &tmpexactCtx_t));
588: PetscCall(PetscCalloc4(NfNew, &tmplowerBound, NfNew, &tmplowerCtx, NfNew, &tmpupperBound, NfNew, &tmpupperCtx));
589: for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f];
590: for (f = 0; f < Nf; ++f) tmpexactCtx[f] = prob->exactCtx[f];
591: for (f = 0; f < Nf; ++f) tmpexactSol_t[f] = prob->exactSol_t[f];
592: for (f = 0; f < Nf; ++f) tmpexactCtx_t[f] = prob->exactCtx_t[f];
593: for (f = 0; f < Nf; ++f) tmplowerBound[f] = prob->lowerBound[f];
594: for (f = 0; f < Nf; ++f) tmplowerCtx[f] = prob->lowerCtx[f];
595: for (f = 0; f < Nf; ++f) tmpupperBound[f] = prob->upperBound[f];
596: for (f = 0; f < Nf; ++f) tmpupperCtx[f] = prob->upperCtx[f];
597: for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL;
598: for (f = Nf; f < NfNew; ++f) tmpexactCtx[f] = NULL;
599: for (f = Nf; f < NfNew; ++f) tmpexactSol_t[f] = NULL;
600: for (f = Nf; f < NfNew; ++f) tmpexactCtx_t[f] = NULL;
601: for (f = Nf; f < NfNew; ++f) tmplowerBound[f] = NULL;
602: for (f = Nf; f < NfNew; ++f) tmplowerCtx[f] = NULL;
603: for (f = Nf; f < NfNew; ++f) tmpupperBound[f] = NULL;
604: for (f = Nf; f < NfNew; ++f) tmpupperCtx[f] = NULL;
605: PetscCall(PetscFree4(prob->exactSol, prob->exactCtx, prob->exactSol_t, prob->exactCtx_t));
606: PetscCall(PetscFree4(prob->lowerBound, prob->lowerCtx, prob->upperBound, prob->upperCtx));
607: prob->exactSol = tmpexactSol;
608: prob->exactCtx = tmpexactCtx;
609: prob->exactSol_t = tmpexactSol_t;
610: prob->exactCtx_t = tmpexactCtx_t;
611: prob->lowerBound = tmplowerBound;
612: prob->lowerCtx = tmplowerCtx;
613: prob->upperBound = tmpupperBound;
614: prob->upperCtx = tmpupperCtx;
615: PetscFunctionReturn(PETSC_SUCCESS);
616: }
618: /*@
619: PetscDSDestroy - Destroys a `PetscDS` object
621: Collective
623: Input Parameter:
624: . ds - the `PetscDS` object to destroy
626: Level: developer
628: .seealso: `PetscDSView()`
629: @*/
630: PetscErrorCode PetscDSDestroy(PetscDS *ds)
631: {
632: PetscFunctionBegin;
633: if (!*ds) PetscFunctionReturn(PETSC_SUCCESS);
636: if (--((PetscObject)*ds)->refct > 0) {
637: *ds = NULL;
638: PetscFunctionReturn(PETSC_SUCCESS);
639: }
640: ((PetscObject)*ds)->refct = 0;
641: if ((*ds)->subprobs) {
642: PetscInt dim;
644: PetscCall(PetscDSGetSpatialDimension(*ds, &dim));
645: for (PetscInt d = 0; d < dim; ++d) PetscCall(PetscDSDestroy(&(*ds)->subprobs[d]));
646: }
647: PetscCall(PetscFree((*ds)->subprobs));
648: PetscCall(PetscDSDestroyStructs_Static(*ds));
649: for (PetscInt f = 0; f < (*ds)->Nf; ++f) PetscCall(PetscObjectDereference((*ds)->disc[f]));
650: PetscCall(PetscFree4((*ds)->disc, (*ds)->implicit, (*ds)->cohesive, (*ds)->jetDegree));
651: PetscCall(PetscWeakFormDestroy(&(*ds)->wf));
652: PetscCall(PetscFree2((*ds)->update, (*ds)->ctx));
653: PetscCall(PetscFree4((*ds)->exactSol, (*ds)->exactCtx, (*ds)->exactSol_t, (*ds)->exactCtx_t));
654: PetscCall(PetscFree4((*ds)->lowerBound, (*ds)->lowerCtx, (*ds)->upperBound, (*ds)->upperCtx));
655: PetscTryTypeMethod(*ds, destroy);
656: PetscCall(PetscDSDestroyBoundary(*ds));
657: PetscCall(PetscFree((*ds)->constants));
658: for (PetscInt c = 0; c < DM_NUM_POLYTOPES; ++c) {
659: const PetscInt Na = DMPolytopeTypeGetNumArrangements((DMPolytopeType)c);
660: if ((*ds)->quadPerm[c])
661: for (PetscInt o = 0; o < Na; ++o) PetscCall(ISDestroy(&(*ds)->quadPerm[c][o]));
662: PetscCall(PetscFree((*ds)->quadPerm[c]));
663: (*ds)->quadPerm[c] = NULL;
664: }
665: PetscCall(PetscHeaderDestroy(ds));
666: PetscFunctionReturn(PETSC_SUCCESS);
667: }
669: /*@
670: PetscDSCreate - Creates an empty `PetscDS` object. The type can then be set with `PetscDSSetType()`.
672: Collective
674: Input Parameter:
675: . comm - The communicator for the `PetscDS` object
677: Output Parameter:
678: . ds - The `PetscDS` object
680: Level: beginner
682: .seealso: `PetscDS`, `PetscDSSetType()`, `PETSCDSBASIC`, `PetscDSType`, `PetscDSDestroy()`
683: @*/
684: PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *ds)
685: {
686: PetscDS p;
688: PetscFunctionBegin;
689: PetscAssertPointer(ds, 2);
690: PetscCall(PetscDSInitializePackage());
692: PetscCall(PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView));
693: p->Nf = 0;
694: p->setup = PETSC_FALSE;
695: p->numConstants = 0;
696: p->numFuncConstants = 3; // Row and col fields, cell size
697: p->dimEmbed = -1;
698: p->useJacPre = PETSC_TRUE;
699: p->forceQuad = PETSC_TRUE;
700: PetscCall(PetscMalloc1(p->numConstants + p->numFuncConstants, &p->constants));
701: PetscCall(PetscWeakFormCreate(comm, &p->wf));
702: PetscCall(PetscArrayzero(p->quadPerm, DM_NUM_POLYTOPES));
703: *ds = p;
704: PetscFunctionReturn(PETSC_SUCCESS);
705: }
707: /*@
708: PetscDSGetNumFields - Returns the number of fields in the `PetscDS`
710: Not Collective
712: Input Parameter:
713: . prob - The `PetscDS` object
715: Output Parameter:
716: . Nf - The number of fields
718: Level: beginner
720: .seealso: `PetscDS`, `PetscDSGetSpatialDimension()`, `PetscDSCreate()`
721: @*/
722: PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf)
723: {
724: PetscFunctionBegin;
726: PetscAssertPointer(Nf, 2);
727: *Nf = prob->Nf;
728: PetscFunctionReturn(PETSC_SUCCESS);
729: }
731: /*@
732: PetscDSGetSpatialDimension - Returns the spatial dimension of the `PetscDS`, meaning the topological dimension of the discretizations
734: Not Collective
736: Input Parameter:
737: . prob - The `PetscDS` object
739: Output Parameter:
740: . dim - The spatial dimension
742: Level: beginner
744: .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
745: @*/
746: PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim)
747: {
748: PetscFunctionBegin;
750: PetscAssertPointer(dim, 2);
751: *dim = 0;
752: if (prob->Nf) {
753: PetscObject obj;
754: PetscClassId id;
756: PetscCall(PetscDSGetDiscretization(prob, 0, &obj));
757: if (obj) {
758: PetscCall(PetscObjectGetClassId(obj, &id));
759: if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetSpatialDimension((PetscFE)obj, dim));
760: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetSpatialDimension((PetscFV)obj, dim));
761: else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
762: }
763: }
764: PetscFunctionReturn(PETSC_SUCCESS);
765: }
767: /*@
768: PetscDSGetCoordinateDimension - Returns the coordinate dimension of the `PetscDS`, meaning the dimension of the space into which the discretiaztions are embedded
770: Not Collective
772: Input Parameter:
773: . prob - The `PetscDS` object
775: Output Parameter:
776: . dimEmbed - The coordinate dimension
778: Level: beginner
780: .seealso: `PetscDS`, `PetscDSSetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
781: @*/
782: PetscErrorCode PetscDSGetCoordinateDimension(PetscDS prob, PetscInt *dimEmbed)
783: {
784: PetscFunctionBegin;
786: PetscAssertPointer(dimEmbed, 2);
787: PetscCheck(prob->dimEmbed >= 0, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONGSTATE, "No coordinate dimension set for this DS");
788: *dimEmbed = prob->dimEmbed;
789: PetscFunctionReturn(PETSC_SUCCESS);
790: }
792: /*@
793: PetscDSSetCoordinateDimension - Set the coordinate dimension of the `PetscDS`, meaning the dimension of the space into which the discretiaztions are embedded
795: Logically Collective
797: Input Parameters:
798: + prob - The `PetscDS` object
799: - dimEmbed - The coordinate dimension
801: Level: beginner
803: .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
804: @*/
805: PetscErrorCode PetscDSSetCoordinateDimension(PetscDS prob, PetscInt dimEmbed)
806: {
807: PetscFunctionBegin;
809: PetscCheck(dimEmbed >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Coordinate dimension must be non-negative, not %" PetscInt_FMT, dimEmbed);
810: prob->dimEmbed = dimEmbed;
811: PetscFunctionReturn(PETSC_SUCCESS);
812: }
814: /*@
815: PetscDSGetForceQuad - Returns the flag to force matching quadratures among the field discretizations
817: Not collective
819: Input Parameter:
820: . ds - The `PetscDS` object
822: Output Parameter:
823: . forceQuad - The flag
825: Level: intermediate
827: .seealso: `PetscDS`, `PetscDSSetForceQuad()`, `PetscDSGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
828: @*/
829: PetscErrorCode PetscDSGetForceQuad(PetscDS ds, PetscBool *forceQuad)
830: {
831: PetscFunctionBegin;
833: PetscAssertPointer(forceQuad, 2);
834: *forceQuad = ds->forceQuad;
835: PetscFunctionReturn(PETSC_SUCCESS);
836: }
838: /*@
839: PetscDSSetForceQuad - Set the flag to force matching quadratures among the field discretizations
841: Logically collective on ds
843: Input Parameters:
844: + ds - The `PetscDS` object
845: - forceQuad - The flag
847: Level: intermediate
849: .seealso: `PetscDS`, `PetscDSGetForceQuad()`, `PetscDSGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
850: @*/
851: PetscErrorCode PetscDSSetForceQuad(PetscDS ds, PetscBool forceQuad)
852: {
853: PetscFunctionBegin;
855: ds->forceQuad = forceQuad;
856: PetscFunctionReturn(PETSC_SUCCESS);
857: }
859: /*@
860: PetscDSIsCohesive - Returns the flag indicating that this `PetscDS` is for a cohesive cell
862: Not Collective
864: Input Parameter:
865: . ds - The `PetscDS` object
867: Output Parameter:
868: . isCohesive - The flag
870: Level: developer
872: .seealso: `PetscDS`, `PetscDSGetNumCohesive()`, `PetscDSGetCohesive()`, `PetscDSSetCohesive()`, `PetscDSCreate()`
873: @*/
874: PetscErrorCode PetscDSIsCohesive(PetscDS ds, PetscBool *isCohesive)
875: {
876: PetscFunctionBegin;
878: PetscAssertPointer(isCohesive, 2);
879: *isCohesive = ds->isCohesive;
880: PetscFunctionReturn(PETSC_SUCCESS);
881: }
883: /*@
884: PetscDSGetNumCohesive - Returns the number of cohesive fields, meaning those defined on the interior of a cohesive cell
886: Not Collective
888: Input Parameter:
889: . ds - The `PetscDS` object
891: Output Parameter:
892: . numCohesive - The number of cohesive fields
894: Level: developer
896: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSCreate()`
897: @*/
898: PetscErrorCode PetscDSGetNumCohesive(PetscDS ds, PetscInt *numCohesive)
899: {
900: PetscFunctionBegin;
902: PetscAssertPointer(numCohesive, 2);
903: *numCohesive = 0;
904: for (PetscInt f = 0; f < ds->Nf; ++f) *numCohesive += ds->cohesive[f] ? 1 : 0;
905: PetscFunctionReturn(PETSC_SUCCESS);
906: }
908: /*@
909: PetscDSGetCohesive - Returns the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell
911: Not Collective
913: Input Parameters:
914: + ds - The `PetscDS` object
915: - f - The field index
917: Output Parameter:
918: . isCohesive - The flag
920: Level: developer
922: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
923: @*/
924: PetscErrorCode PetscDSGetCohesive(PetscDS ds, PetscInt f, PetscBool *isCohesive)
925: {
926: PetscFunctionBegin;
928: PetscAssertPointer(isCohesive, 3);
929: 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);
930: *isCohesive = ds->cohesive[f];
931: PetscFunctionReturn(PETSC_SUCCESS);
932: }
934: /*@
935: PetscDSSetCohesive - Set the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell
937: Not Collective
939: Input Parameters:
940: + ds - The `PetscDS` object
941: . f - The field index
942: - isCohesive - The flag for a cohesive field
944: Level: developer
946: .seealso: `PetscDS`, `PetscDSGetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
947: @*/
948: PetscErrorCode PetscDSSetCohesive(PetscDS ds, PetscInt f, PetscBool isCohesive)
949: {
950: PetscFunctionBegin;
952: 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);
953: ds->cohesive[f] = isCohesive;
954: ds->isCohesive = PETSC_FALSE;
955: for (PetscInt i = 0; i < ds->Nf; ++i) ds->isCohesive = ds->isCohesive || ds->cohesive[f] ? PETSC_TRUE : PETSC_FALSE;
956: PetscFunctionReturn(PETSC_SUCCESS);
957: }
959: /*@
960: PetscDSGetTotalDimension - Returns the total size of the approximation space for this system
962: Not Collective
964: Input Parameter:
965: . prob - The `PetscDS` object
967: Output Parameter:
968: . dim - The total problem dimension
970: Level: beginner
972: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
973: @*/
974: PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim)
975: {
976: PetscFunctionBegin;
978: PetscCall(PetscDSSetUp(prob));
979: PetscAssertPointer(dim, 2);
980: *dim = prob->totDim;
981: PetscFunctionReturn(PETSC_SUCCESS);
982: }
984: /*@
985: PetscDSGetTotalComponents - Returns the total number of components in this system
987: Not Collective
989: Input Parameter:
990: . prob - The `PetscDS` object
992: Output Parameter:
993: . Nc - The total number of components
995: Level: beginner
997: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
998: @*/
999: PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc)
1000: {
1001: PetscFunctionBegin;
1003: PetscCall(PetscDSSetUp(prob));
1004: PetscAssertPointer(Nc, 2);
1005: *Nc = prob->totComp;
1006: PetscFunctionReturn(PETSC_SUCCESS);
1007: }
1009: /*@
1010: PetscDSGetDiscretization - Returns the discretization object for the given field
1012: Not Collective
1014: Input Parameters:
1015: + prob - The `PetscDS` object
1016: - f - The field number
1018: Output Parameter:
1019: . disc - The discretization object, this can be a `PetscFE` or a `PetscFV`
1021: Level: beginner
1023: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1024: @*/
1025: PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc)
1026: {
1027: PetscFunctionBeginHot;
1029: PetscAssertPointer(disc, 3);
1030: 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);
1031: *disc = prob->disc[f];
1032: PetscFunctionReturn(PETSC_SUCCESS);
1033: }
1035: /*@
1036: PetscDSSetDiscretization - Sets the discretization object for the given field
1038: Not Collective
1040: Input Parameters:
1041: + prob - The `PetscDS` object
1042: . f - The field number
1043: - disc - The discretization object, this can be a `PetscFE` or a `PetscFV`
1045: Level: beginner
1047: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSGetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1048: @*/
1049: PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc)
1050: {
1051: PetscFunctionBegin;
1053: if (disc) PetscAssertPointer(disc, 3);
1054: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1055: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
1056: PetscCall(PetscObjectDereference(prob->disc[f]));
1057: prob->disc[f] = disc;
1058: PetscCall(PetscObjectReference(disc));
1059: if (disc) {
1060: PetscClassId id;
1062: PetscCall(PetscObjectGetClassId(disc, &id));
1063: if (id == PETSCFE_CLASSID) {
1064: PetscCall(PetscDSSetImplicit(prob, f, PETSC_TRUE));
1065: } else if (id == PETSCFV_CLASSID) {
1066: PetscCall(PetscDSSetImplicit(prob, f, PETSC_FALSE));
1067: }
1068: PetscCall(PetscDSSetJetDegree(prob, f, 1));
1069: }
1070: PetscFunctionReturn(PETSC_SUCCESS);
1071: }
1073: /*@
1074: PetscDSGetWeakForm - Returns the weak form object from within the `PetscDS`
1076: Not Collective
1078: Input Parameter:
1079: . ds - The `PetscDS` object
1081: Output Parameter:
1082: . wf - The weak form object
1084: Level: beginner
1086: .seealso: `PetscWeakForm`, `PetscDSSetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1087: @*/
1088: PetscErrorCode PetscDSGetWeakForm(PetscDS ds, PetscWeakForm *wf)
1089: {
1090: PetscFunctionBegin;
1092: PetscAssertPointer(wf, 2);
1093: *wf = ds->wf;
1094: PetscFunctionReturn(PETSC_SUCCESS);
1095: }
1097: /*@
1098: PetscDSSetWeakForm - Sets the weak form object to be used by the `PetscDS`
1100: Not Collective
1102: Input Parameters:
1103: + ds - The `PetscDS` object
1104: - wf - The weak form object
1106: Level: beginner
1108: .seealso: `PetscWeakForm`, `PetscDSGetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1109: @*/
1110: PetscErrorCode PetscDSSetWeakForm(PetscDS ds, PetscWeakForm wf)
1111: {
1112: PetscFunctionBegin;
1115: PetscCall(PetscObjectDereference((PetscObject)ds->wf));
1116: ds->wf = wf;
1117: PetscCall(PetscObjectReference((PetscObject)wf));
1118: PetscCall(PetscWeakFormSetNumFields(wf, ds->Nf));
1119: PetscFunctionReturn(PETSC_SUCCESS);
1120: }
1122: /*@
1123: PetscDSAddDiscretization - Adds a discretization object
1125: Not Collective
1127: Input Parameters:
1128: + prob - The `PetscDS` object
1129: - disc - The discretization object, this can be a `PetscFE` or `PetscFV`
1131: Level: beginner
1133: .seealso: `PetscWeakForm`, `PetscFE`, `PetscFV`, `PetscDSGetDiscretization()`, `PetscDSSetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1134: @*/
1135: PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc)
1136: {
1137: PetscFunctionBegin;
1138: PetscCall(PetscDSSetDiscretization(prob, prob->Nf, disc));
1139: PetscFunctionReturn(PETSC_SUCCESS);
1140: }
1142: /*@
1143: PetscDSGetQuadrature - Returns the quadrature, which must agree for all fields in the `PetscDS`
1145: Not Collective
1147: Input Parameter:
1148: . prob - The `PetscDS` object
1150: Output Parameter:
1151: . q - The quadrature object
1153: Level: intermediate
1155: .seealso: `PetscDS`, `PetscQuadrature`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1156: @*/
1157: PetscErrorCode PetscDSGetQuadrature(PetscDS prob, PetscQuadrature *q)
1158: {
1159: PetscObject obj;
1160: PetscClassId id;
1162: PetscFunctionBegin;
1163: *q = NULL;
1164: if (!prob->Nf) PetscFunctionReturn(PETSC_SUCCESS);
1165: PetscCall(PetscDSGetDiscretization(prob, 0, &obj));
1166: PetscCall(PetscObjectGetClassId(obj, &id));
1167: if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetQuadrature((PetscFE)obj, q));
1168: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetQuadrature((PetscFV)obj, q));
1169: else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
1170: PetscFunctionReturn(PETSC_SUCCESS);
1171: }
1173: /*@
1174: PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for `TSARKIMEX`
1176: Not Collective
1178: Input Parameters:
1179: + prob - The `PetscDS` object
1180: - f - The field number
1182: Output Parameter:
1183: . implicit - The flag indicating what kind of solve to use for this field
1185: Level: developer
1187: .seealso: `TSARKIMEX`, `PetscDS`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1188: @*/
1189: PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit)
1190: {
1191: PetscFunctionBegin;
1193: PetscAssertPointer(implicit, 3);
1194: 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);
1195: *implicit = prob->implicit[f];
1196: PetscFunctionReturn(PETSC_SUCCESS);
1197: }
1199: /*@
1200: PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for `TSARKIMEX`
1202: Not Collective
1204: Input Parameters:
1205: + prob - The `PetscDS` object
1206: . f - The field number
1207: - implicit - The flag indicating what kind of solve to use for this field
1209: Level: developer
1211: .seealso: `TSARKIMEX`, `PetscDSGetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1212: @*/
1213: PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit)
1214: {
1215: PetscFunctionBegin;
1217: 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);
1218: prob->implicit[f] = implicit;
1219: PetscFunctionReturn(PETSC_SUCCESS);
1220: }
1222: /*@
1223: PetscDSGetJetDegree - Returns the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate.
1225: Not Collective
1227: Input Parameters:
1228: + ds - The `PetscDS` object
1229: - f - The field number
1231: Output Parameter:
1232: . k - The highest derivative we need to tabulate
1234: Level: developer
1236: .seealso: `PetscDS`, `PetscDSSetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1237: @*/
1238: PetscErrorCode PetscDSGetJetDegree(PetscDS ds, PetscInt f, PetscInt *k)
1239: {
1240: PetscFunctionBegin;
1242: PetscAssertPointer(k, 3);
1243: 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);
1244: *k = ds->jetDegree[f];
1245: PetscFunctionReturn(PETSC_SUCCESS);
1246: }
1248: /*@
1249: PetscDSSetJetDegree - Set the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate.
1251: Not Collective
1253: Input Parameters:
1254: + ds - The `PetscDS` object
1255: . f - The field number
1256: - k - The highest derivative we need to tabulate
1258: Level: developer
1260: .seealso: `PetscDS`, `PetscDSGetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1261: @*/
1262: PetscErrorCode PetscDSSetJetDegree(PetscDS ds, PetscInt f, PetscInt k)
1263: {
1264: PetscFunctionBegin;
1266: 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);
1267: ds->jetDegree[f] = k;
1268: PetscFunctionReturn(PETSC_SUCCESS);
1269: }
1271: /*@
1272: PetscDSGetObjective - Get the pointwise objective function for a given test field that was provided with `PetscDSSetObjective()`
1274: Not Collective
1276: Input Parameters:
1277: + ds - The `PetscDS`
1278: - f - The test field number
1280: Output Parameter:
1281: . obj - integrand for the test function term, see `PetscPointFn`
1283: Level: intermediate
1285: Note:
1286: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi\,\mathrm{obj}(u, u_t, \nabla u, x, t)$
1288: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSSetObjective()`, `PetscDSGetResidual()`
1289: @*/
1290: PetscErrorCode PetscDSGetObjective(PetscDS ds, PetscInt f, PetscPointFn **obj)
1291: {
1292: PetscPointFn **tmp;
1293: PetscInt n;
1295: PetscFunctionBegin;
1297: PetscAssertPointer(obj, 3);
1298: 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);
1299: PetscCall(PetscWeakFormGetObjective(ds->wf, NULL, 0, f, 0, &n, &tmp));
1300: *obj = tmp ? tmp[0] : NULL;
1301: PetscFunctionReturn(PETSC_SUCCESS);
1302: }
1304: /*@
1305: PetscDSSetObjective - Set the pointwise objective function for a given test field
1307: Not Collective
1309: Input Parameters:
1310: + ds - The `PetscDS`
1311: . f - The test field number
1312: - obj - integrand for the test function term, see `PetscPointFn`
1314: Level: intermediate
1316: Note:
1317: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi\,\mathrm{obj}(u, u_t, \nabla u, x, t)$
1319: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSGetObjective()`, `PetscDSSetResidual()`
1320: @*/
1321: PetscErrorCode PetscDSSetObjective(PetscDS ds, PetscInt f, PetscPointFn *obj)
1322: {
1323: PetscFunctionBegin;
1326: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1327: PetscCall(PetscWeakFormSetIndexObjective(ds->wf, NULL, 0, f, 0, 0, obj));
1328: PetscFunctionReturn(PETSC_SUCCESS);
1329: }
1331: /*@
1332: PetscDSGetResidual - Get the pointwise residual function for a given test field
1334: Not Collective
1336: Input Parameters:
1337: + ds - The `PetscDS`
1338: - f - The test field number
1340: Output Parameters:
1341: + f0 - integrand for the test function term, see `PetscPointFn`
1342: - f1 - integrand for the test function gradient term, see `PetscPointFn`
1344: Level: intermediate
1346: Note:
1347: 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)$
1349: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSSetResidual()`
1350: @*/
1351: PetscErrorCode PetscDSGetResidual(PetscDS ds, PetscInt f, PetscPointFn **f0, PetscPointFn **f1)
1352: {
1353: PetscPointFn **tmp0, **tmp1;
1354: PetscInt n0, n1;
1356: PetscFunctionBegin;
1358: 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);
1359: PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1));
1360: *f0 = tmp0 ? tmp0[0] : NULL;
1361: *f1 = tmp1 ? tmp1[0] : NULL;
1362: PetscFunctionReturn(PETSC_SUCCESS);
1363: }
1365: /*@
1366: PetscDSSetResidual - Set the pointwise residual function for a given test field
1368: Not Collective
1370: Input Parameters:
1371: + ds - The `PetscDS`
1372: . f - The test field number
1373: . f0 - integrand for the test function term, see `PetscPointFn`
1374: - f1 - integrand for the test function gradient term, see `PetscPointFn`
1376: Level: intermediate
1378: Note:
1379: 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)$
1381: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSGetResidual()`
1382: @*/
1383: PetscErrorCode PetscDSSetResidual(PetscDS ds, PetscInt f, PetscPointFn *f0, PetscPointFn *f1)
1384: {
1385: PetscFunctionBegin;
1389: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1390: PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1));
1391: PetscFunctionReturn(PETSC_SUCCESS);
1392: }
1394: /*@
1395: PetscDSGetRHSResidual - Get the pointwise RHS residual function for explicit timestepping for a given test field
1397: Not Collective
1399: Input Parameters:
1400: + ds - The `PetscDS`
1401: - f - The test field number
1403: Output Parameters:
1404: + f0 - integrand for the test function term, see `PetscPointFn`
1405: - f1 - integrand for the test function gradient term, see `PetscPointFn`
1407: Level: intermediate
1409: Note:
1410: 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)$
1412: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSSetRHSResidual()`
1413: @*/
1414: PetscErrorCode PetscDSGetRHSResidual(PetscDS ds, PetscInt f, PetscPointFn **f0, PetscPointFn **f1)
1415: {
1416: PetscPointFn **tmp0, **tmp1;
1417: PetscInt n0, n1;
1419: PetscFunctionBegin;
1421: 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);
1422: PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 100, &n0, &tmp0, &n1, &tmp1));
1423: *f0 = tmp0 ? tmp0[0] : NULL;
1424: *f1 = tmp1 ? tmp1[0] : NULL;
1425: PetscFunctionReturn(PETSC_SUCCESS);
1426: }
1428: /*@
1429: PetscDSSetRHSResidual - Set the pointwise residual function for explicit timestepping for a given test field
1431: Not Collective
1433: Input Parameters:
1434: + ds - The `PetscDS`
1435: . f - The test field number
1436: . f0 - integrand for the test function term, see `PetscPointFn`
1437: - f1 - integrand for the test function gradient term, see `PetscPointFn`
1439: Level: intermediate
1441: Note:
1442: 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)$
1444: .seealso: `PetscDS`, `PetscDSGetResidual()`
1445: @*/
1446: PetscErrorCode PetscDSSetRHSResidual(PetscDS ds, PetscInt f, PetscPointFn *f0, PetscPointFn *f1)
1447: {
1448: PetscFunctionBegin;
1452: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1453: PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 100, 0, f0, 0, f1));
1454: PetscFunctionReturn(PETSC_SUCCESS);
1455: }
1457: /*@
1458: PetscDSHasJacobian - Checks that the Jacobian functions have been set
1460: Not Collective
1462: Input Parameter:
1463: . ds - The `PetscDS`
1465: Output Parameter:
1466: . hasJac - flag that indicates the pointwise function for the Jacobian has been set
1468: Level: intermediate
1470: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1471: @*/
1472: PetscErrorCode PetscDSHasJacobian(PetscDS ds, PetscBool *hasJac)
1473: {
1474: PetscFunctionBegin;
1476: PetscCall(PetscWeakFormHasJacobian(ds->wf, hasJac));
1477: PetscFunctionReturn(PETSC_SUCCESS);
1478: }
1480: /*@
1481: PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field
1483: Not Collective
1485: Input Parameters:
1486: + ds - The `PetscDS`
1487: . f - The test field number
1488: - g - The field number
1490: Output Parameters:
1491: + g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1492: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1493: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1494: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1496: Level: intermediate
1498: Note:
1499: We are using a first order FEM model for the weak form\:
1501: $$
1502: \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
1503: + \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
1504: $$
1506: .seealso: `PetscDS`, `PetscDSSetJacobian()`, `PetscPointJacFn`
1507: @*/
1508: PetscErrorCode PetscDSGetJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn **g0, PetscPointJacFn **g1, PetscPointJacFn **g2, PetscPointJacFn **g3)
1509: {
1510: PetscPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
1511: PetscInt n0, n1, n2, n3;
1513: PetscFunctionBegin;
1515: 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);
1516: 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);
1517: PetscCall(PetscWeakFormGetJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1518: *g0 = tmp0 ? tmp0[0] : NULL;
1519: *g1 = tmp1 ? tmp1[0] : NULL;
1520: *g2 = tmp2 ? tmp2[0] : NULL;
1521: *g3 = tmp3 ? tmp3[0] : NULL;
1522: PetscFunctionReturn(PETSC_SUCCESS);
1523: }
1525: /*@
1526: PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields
1528: Not Collective
1530: Input Parameters:
1531: + ds - The `PetscDS`
1532: . f - The test field number
1533: . g - The field number
1534: . g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1535: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1536: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1537: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1539: Level: intermediate
1541: Note:
1542: We are using a first order FEM model for the weak form\:
1544: $$
1545: \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
1546: + \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
1547: $$
1549: .seealso: `PetscDS`, `PetscDSGetJacobian()`, `PetscPointJacFn`
1550: @*/
1551: PetscErrorCode PetscDSSetJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn *g0, PetscPointJacFn *g1, PetscPointJacFn *g2, PetscPointJacFn *g3)
1552: {
1553: PetscFunctionBegin;
1559: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1560: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1561: PetscCall(PetscWeakFormSetIndexJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1562: PetscFunctionReturn(PETSC_SUCCESS);
1563: }
1565: /*@
1566: PetscDSUseJacobianPreconditioner - Set whether to construct a Jacobian preconditioner
1568: Not Collective
1570: Input Parameters:
1571: + prob - The `PetscDS`
1572: - useJacPre - flag that enables construction of a Jacobian preconditioner
1574: Level: intermediate
1576: Developer Note:
1577: Should be called `PetscDSSetUseJacobianPreconditioner()`
1579: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1580: @*/
1581: PetscErrorCode PetscDSUseJacobianPreconditioner(PetscDS prob, PetscBool useJacPre)
1582: {
1583: PetscFunctionBegin;
1585: prob->useJacPre = useJacPre;
1586: PetscFunctionReturn(PETSC_SUCCESS);
1587: }
1589: /*@
1590: PetscDSHasJacobianPreconditioner - Checks if a Jacobian matrix for constructing a preconditioner has been set
1592: Not Collective
1594: Input Parameter:
1595: . ds - The `PetscDS`
1597: Output Parameter:
1598: . hasJacPre - the flag
1600: Level: intermediate
1602: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1603: @*/
1604: PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS ds, PetscBool *hasJacPre)
1605: {
1606: PetscFunctionBegin;
1608: *hasJacPre = PETSC_FALSE;
1609: if (!ds->useJacPre) PetscFunctionReturn(PETSC_SUCCESS);
1610: PetscCall(PetscWeakFormHasJacobianPreconditioner(ds->wf, hasJacPre));
1611: PetscFunctionReturn(PETSC_SUCCESS);
1612: }
1614: /*@
1615: PetscDSGetJacobianPreconditioner - Get the pointwise Jacobian function for given test and basis field that constructs the matrix used
1616: to compute the preconditioner. If this is missing, the system matrix is used to build the preconditioner.
1618: Not Collective
1620: Input Parameters:
1621: + ds - The `PetscDS`
1622: . f - The test field number
1623: - g - The field number
1625: Output Parameters:
1626: + g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1627: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1628: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1629: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1631: Level: intermediate
1633: Note:
1634: We are using a first order FEM model for the weak form\:
1636: $$
1637: \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
1638: + \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
1639: $$
1641: Developer Note:
1642: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
1644: .seealso: `PetscDS`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`, `PetscPointJacFn`
1645: @*/
1646: PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn **g0, PetscPointJacFn **g1, PetscPointJacFn **g2, PetscPointJacFn **g3)
1647: {
1648: PetscPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
1649: PetscInt n0, n1, n2, n3;
1651: PetscFunctionBegin;
1653: 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);
1654: 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);
1655: PetscCall(PetscWeakFormGetJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1656: *g0 = tmp0 ? tmp0[0] : NULL;
1657: *g1 = tmp1 ? tmp1[0] : NULL;
1658: *g2 = tmp2 ? tmp2[0] : NULL;
1659: *g3 = tmp3 ? tmp3[0] : NULL;
1660: PetscFunctionReturn(PETSC_SUCCESS);
1661: }
1663: /*@
1664: PetscDSSetJacobianPreconditioner - Set the pointwise Jacobian function for given test and basis fields that constructs the matrix used
1665: to compute the preconditioner. If this is missing, the system matrix is used to build the preconditioner.
1667: Not Collective
1669: Input Parameters:
1670: + ds - The `PetscDS`
1671: . f - The test field number
1672: . g - The field number
1673: . g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1674: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1675: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1676: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1678: Level: intermediate
1680: Note:
1681: We are using a first order FEM model for the weak form\:
1683: $$
1684: \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
1685: + \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
1686: $$
1688: Developer Note:
1689: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
1691: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobian()`, `PetscPointJacFn`
1692: @*/
1693: PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn *g0, PetscPointJacFn *g1, PetscPointJacFn *g2, PetscPointJacFn *g3)
1694: {
1695: PetscFunctionBegin;
1701: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1702: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1703: PetscCall(PetscWeakFormSetIndexJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1704: PetscFunctionReturn(PETSC_SUCCESS);
1705: }
1707: /*@
1708: PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, $dF/du_t$, has been set
1710: Not Collective
1712: Input Parameter:
1713: . ds - The `PetscDS`
1715: Output Parameter:
1716: . hasDynJac - flag that pointwise function for dynamic Jacobian has been set
1718: Level: intermediate
1720: .seealso: `PetscDS`, `PetscDSGetDynamicJacobian()`, `PetscDSSetDynamicJacobian()`, `PetscDSGetJacobian()`
1721: @*/
1722: PetscErrorCode PetscDSHasDynamicJacobian(PetscDS ds, PetscBool *hasDynJac)
1723: {
1724: PetscFunctionBegin;
1726: PetscCall(PetscWeakFormHasDynamicJacobian(ds->wf, hasDynJac));
1727: PetscFunctionReturn(PETSC_SUCCESS);
1728: }
1730: /*@
1731: PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, $dF/du_t$, function for given test and basis field
1733: Not Collective
1735: Input Parameters:
1736: + ds - The `PetscDS`
1737: . f - The test field number
1738: - g - The field number
1740: Output Parameters:
1741: + g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1742: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1743: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1744: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1746: Level: intermediate
1748: Note:
1749: We are using a first order FEM model for the weak form\:
1751: $$
1752: \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
1753: + \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
1754: $$
1756: .seealso: `PetscDS`, `PetscDSSetJacobian()`, `PetscDSSetDynamicJacobian()`, `PetscPointJacFn`
1757: @*/
1758: PetscErrorCode PetscDSGetDynamicJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn **g0, PetscPointJacFn **g1, PetscPointJacFn **g2, PetscPointJacFn **g3)
1759: {
1760: PetscPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
1761: PetscInt n0, n1, n2, n3;
1763: PetscFunctionBegin;
1765: 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);
1766: 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);
1767: PetscCall(PetscWeakFormGetDynamicJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1768: *g0 = tmp0 ? tmp0[0] : NULL;
1769: *g1 = tmp1 ? tmp1[0] : NULL;
1770: *g2 = tmp2 ? tmp2[0] : NULL;
1771: *g3 = tmp3 ? tmp3[0] : NULL;
1772: PetscFunctionReturn(PETSC_SUCCESS);
1773: }
1775: /*@
1776: PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, $dF/du_t$, function for given test and basis fields
1778: Not Collective
1780: Input Parameters:
1781: + ds - The `PetscDS`
1782: . f - The test field number
1783: . g - The field number
1784: . g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1785: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1786: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1787: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1789: Level: intermediate
1791: Note:
1792: We are using a first order FEM model for the weak form\:
1794: $$
1795: \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
1796: + \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
1797: $$
1799: .seealso: `PetscDS`, `PetscDSGetDynamicJacobian()`, `PetscDSGetJacobian()`, `PetscPointJacFn`
1800: @*/
1801: PetscErrorCode PetscDSSetDynamicJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn *g0, PetscPointJacFn *g1, PetscPointJacFn *g2, PetscPointJacFn *g3)
1802: {
1803: PetscFunctionBegin;
1809: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1810: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1811: PetscCall(PetscWeakFormSetIndexDynamicJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1812: PetscFunctionReturn(PETSC_SUCCESS);
1813: }
1815: /*@
1816: PetscDSGetRiemannSolver - Returns the Riemann solver for the given field
1818: Not Collective
1820: Input Parameters:
1821: + ds - The `PetscDS` object
1822: - f - The field number
1824: Output Parameter:
1825: . r - Riemann solver, see `PetscRiemannFn`
1827: Level: intermediate
1829: .seealso: `PetscDS`, `PetscRiemannFn`, `PetscDSSetRiemannSolver()`
1830: @*/
1831: PetscErrorCode PetscDSGetRiemannSolver(PetscDS ds, PetscInt f, PetscRiemannFn **r)
1832: {
1833: PetscRiemannFn **tmp;
1834: PetscInt n;
1836: PetscFunctionBegin;
1838: PetscAssertPointer(r, 3);
1839: 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);
1840: PetscCall(PetscWeakFormGetRiemannSolver(ds->wf, NULL, 0, f, 0, &n, &tmp));
1841: *r = tmp ? tmp[0] : NULL;
1842: PetscFunctionReturn(PETSC_SUCCESS);
1843: }
1845: /*@
1846: PetscDSSetRiemannSolver - Sets the Riemann solver for the given field
1848: Not Collective
1850: Input Parameters:
1851: + ds - The `PetscDS` object
1852: . f - The field number
1853: - r - Riemann solver, see `PetscRiemannFn`
1855: Level: intermediate
1857: .seealso: `PetscDS`, `PetscRiemannFn`, `PetscDSGetRiemannSolver()`
1858: @*/
1859: PetscErrorCode PetscDSSetRiemannSolver(PetscDS ds, PetscInt f, PetscRiemannFn *r)
1860: {
1861: PetscFunctionBegin;
1864: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1865: PetscCall(PetscWeakFormSetIndexRiemannSolver(ds->wf, NULL, 0, f, 0, 0, r));
1866: PetscFunctionReturn(PETSC_SUCCESS);
1867: }
1869: /*@
1870: PetscDSGetUpdate - Get the pointwise update function for a given field
1872: Not Collective
1874: Input Parameters:
1875: + ds - The `PetscDS`
1876: - f - The field number
1878: Output Parameter:
1879: . update - update function, see `PetscPointFn`
1881: Level: intermediate
1883: .seealso: `PetscDS`, `PetscPointFn`, `PetscDSSetUpdate()`, `PetscDSSetResidual()`
1884: @*/
1885: PetscErrorCode PetscDSGetUpdate(PetscDS ds, PetscInt f, PetscPointFn **update)
1886: {
1887: PetscFunctionBegin;
1889: 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);
1890: if (update) {
1891: PetscAssertPointer(update, 3);
1892: *update = ds->update[f];
1893: }
1894: PetscFunctionReturn(PETSC_SUCCESS);
1895: }
1897: /*@
1898: PetscDSSetUpdate - Set the pointwise update function for a given field
1900: Not Collective
1902: Input Parameters:
1903: + ds - The `PetscDS`
1904: . f - The field number
1905: - update - update function, see `PetscPointFn`
1907: Level: intermediate
1909: .seealso: `PetscDS`, `PetscPointFn`, `PetscDSGetResidual()`
1910: @*/
1911: PetscErrorCode PetscDSSetUpdate(PetscDS ds, PetscInt f, PetscPointFn *update)
1912: {
1913: PetscFunctionBegin;
1916: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1917: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
1918: ds->update[f] = update;
1919: PetscFunctionReturn(PETSC_SUCCESS);
1920: }
1922: /*@
1923: PetscDSGetContext - Returns the context that was passed by `PetscDSSetContext()`
1925: Not Collective
1927: Input Parameters:
1928: + ds - The `PetscDS`
1929: . f - The field number
1930: - ctx - the context
1932: Level: intermediate
1934: Fortran Notes:
1935: This only works when the context is a Fortran derived type or a `PetscObject`. Define `ctx` with
1936: .vb
1937: type(tUsertype), pointer :: ctx
1938: .ve
1940: .seealso: `PetscDS`, `PetscPointFn`, `PetscDSSetContext()`
1941: @*/
1942: PetscErrorCode PetscDSGetContext(PetscDS ds, PetscInt f, PetscCtxRt ctx)
1943: {
1944: PetscFunctionBegin;
1946: 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);
1947: PetscAssertPointer(ctx, 3);
1948: *(void **)ctx = ds->ctx[f];
1949: PetscFunctionReturn(PETSC_SUCCESS);
1950: }
1952: /*@
1953: PetscDSSetContext - Sets the context that is passed back to some of the pointwise function callbacks used by this `PetscDS`
1955: Not Collective
1957: Input Parameters:
1958: + ds - The `PetscDS`
1959: . f - The field number
1960: - ctx - the context
1962: Level: intermediate
1964: .seealso: `PetscDS`, `PetscPointFn`, `PetscDSGetContext()`
1965: @*/
1966: PetscErrorCode PetscDSSetContext(PetscDS ds, PetscInt f, PetscCtx ctx)
1967: {
1968: PetscFunctionBegin;
1970: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1971: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
1972: ds->ctx[f] = ctx;
1973: PetscFunctionReturn(PETSC_SUCCESS);
1974: }
1976: /*@
1977: PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field
1979: Not Collective
1981: Input Parameters:
1982: + ds - The PetscDS
1983: - f - The test field number
1985: Output Parameters:
1986: + f0 - boundary integrand for the test function term, see `PetscBdPointFn`
1987: - f1 - boundary integrand for the test function gradient term, see `PetscBdPointFn`
1989: Level: intermediate
1991: Note:
1992: We are using a first order FEM model for the weak form\:
1994: $$
1995: \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
1996: $$
1998: .seealso: `PetscDS`, `PetscBdPointFn`, `PetscDSSetBdResidual()`
1999: @*/
2000: PetscErrorCode PetscDSGetBdResidual(PetscDS ds, PetscInt f, PetscBdPointFn **f0, PetscBdPointFn **f1)
2001: {
2002: PetscBdPointFn **tmp0, **tmp1;
2003: PetscInt n0, n1;
2005: PetscFunctionBegin;
2007: 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);
2008: PetscCall(PetscWeakFormGetBdResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1));
2009: *f0 = tmp0 ? tmp0[0] : NULL;
2010: *f1 = tmp1 ? tmp1[0] : NULL;
2011: PetscFunctionReturn(PETSC_SUCCESS);
2012: }
2014: /*@
2015: PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field
2017: Not Collective
2019: Input Parameters:
2020: + ds - The `PetscDS`
2021: . f - The test field number
2022: . f0 - boundary integrand for the test function term, see `PetscBdPointFn`
2023: - f1 - boundary integrand for the test function gradient term, see `PetscBdPointFn`
2025: Level: intermediate
2027: Note:
2028: We are using a first order FEM model for the weak form\:
2030: $$
2031: \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
2032: $$
2034: .seealso: `PetscDS`, `PetscBdPointFn`, `PetscDSGetBdResidual()`
2035: @*/
2036: PetscErrorCode PetscDSSetBdResidual(PetscDS ds, PetscInt f, PetscBdPointFn *f0, PetscBdPointFn *f1)
2037: {
2038: PetscFunctionBegin;
2040: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2041: PetscCall(PetscWeakFormSetIndexBdResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1));
2042: PetscFunctionReturn(PETSC_SUCCESS);
2043: }
2045: /*@
2046: PetscDSHasBdJacobian - Indicates that boundary Jacobian functions have been set
2048: Not Collective
2050: Input Parameter:
2051: . ds - The `PetscDS`
2053: Output Parameter:
2054: . hasBdJac - flag that pointwise function for the boundary Jacobian has been set
2056: Level: intermediate
2058: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2059: @*/
2060: PetscErrorCode PetscDSHasBdJacobian(PetscDS ds, PetscBool *hasBdJac)
2061: {
2062: PetscFunctionBegin;
2064: PetscAssertPointer(hasBdJac, 2);
2065: PetscCall(PetscWeakFormHasBdJacobian(ds->wf, hasBdJac));
2066: PetscFunctionReturn(PETSC_SUCCESS);
2067: }
2069: /*@
2070: PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field
2072: Not Collective
2074: Input Parameters:
2075: + ds - The `PetscDS`
2076: . f - The test field number
2077: - g - The field number
2079: Output Parameters:
2080: + g0 - integrand for the test and basis function term, see `PetscBdPointJacFn`
2081: . g1 - integrand for the test function and basis function gradient term, see `PetscBdPointJacFn`
2082: . g2 - integrand for the test function gradient and basis function term, see `PetscBdPointJacFn`
2083: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscBdPointJacFn`
2085: Level: intermediate
2087: Note:
2088: We are using a first order FEM model for the weak form\:
2090: $$
2091: \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
2092: + \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
2093: $$
2095: .seealso: `PetscDS`, `PetscBdPointJacFn`, `PetscDSSetBdJacobian()`
2096: @*/
2097: PetscErrorCode PetscDSGetBdJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscBdPointJacFn **g0, PetscBdPointJacFn **g1, PetscBdPointJacFn **g2, PetscBdPointJacFn **g3)
2098: {
2099: PetscBdPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
2100: PetscInt n0, n1, n2, n3;
2102: PetscFunctionBegin;
2104: 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);
2105: 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);
2106: PetscCall(PetscWeakFormGetBdJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2107: *g0 = tmp0 ? tmp0[0] : NULL;
2108: *g1 = tmp1 ? tmp1[0] : NULL;
2109: *g2 = tmp2 ? tmp2[0] : NULL;
2110: *g3 = tmp3 ? tmp3[0] : NULL;
2111: PetscFunctionReturn(PETSC_SUCCESS);
2112: }
2114: /*@
2115: PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field
2117: Not Collective
2119: Input Parameters:
2120: + ds - The PetscDS
2121: . f - The test field number
2122: . g - The field number
2123: . g0 - integrand for the test and basis function term, see `PetscBdPointJacFn`
2124: . g1 - integrand for the test function and basis function gradient term, see `PetscBdPointJacFn`
2125: . g2 - integrand for the test function gradient and basis function term, see `PetscBdPointJacFn`
2126: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscBdPointJacFn`
2128: Level: intermediate
2130: Note:
2131: We are using a first order FEM model for the weak form\:
2133: $$
2134: \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
2135: + \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
2136: $$
2138: .seealso: `PetscDS`, `PetscBdPointJacFn`, `PetscDSGetBdJacobian()`
2139: @*/
2140: PetscErrorCode PetscDSSetBdJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscBdPointJacFn *g0, PetscBdPointJacFn *g1, PetscBdPointJacFn *g2, PetscBdPointJacFn *g3)
2141: {
2142: PetscFunctionBegin;
2148: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2149: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2150: PetscCall(PetscWeakFormSetIndexBdJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2151: PetscFunctionReturn(PETSC_SUCCESS);
2152: }
2154: /*@
2155: PetscDSHasBdJacobianPreconditioner - Signals that boundary Jacobian preconditioner functions have been set with `PetscDSSetBdJacobianPreconditioner()`
2157: Not Collective
2159: Input Parameter:
2160: . ds - The `PetscDS`
2162: Output Parameter:
2163: . hasBdJacPre - flag that pointwise function for the boundary Jacobian matrix to construct the preconditioner has been set
2165: Level: intermediate
2167: Developer Note:
2168: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
2170: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2171: @*/
2172: PetscErrorCode PetscDSHasBdJacobianPreconditioner(PetscDS ds, PetscBool *hasBdJacPre)
2173: {
2174: PetscFunctionBegin;
2176: PetscAssertPointer(hasBdJacPre, 2);
2177: PetscCall(PetscWeakFormHasBdJacobianPreconditioner(ds->wf, hasBdJacPre));
2178: PetscFunctionReturn(PETSC_SUCCESS);
2179: }
2181: /*@
2182: PetscDSGetBdJacobianPreconditioner - Get the pointwise boundary Jacobian function for given test and basis field that constructs the
2183: matrix used to construct the preconditioner
2185: Not Collective; No Fortran Support
2187: Input Parameters:
2188: + ds - The `PetscDS`
2189: . f - The test field number
2190: - g - The field number
2192: Output Parameters:
2193: + g0 - integrand for the test and basis function term, see `PetscBdPointJacFn`
2194: . g1 - integrand for the test function and basis function gradient term, see `PetscBdPointJacFn`
2195: . g2 - integrand for the test function gradient and basis function term, see `PetscBdPointJacFn`
2196: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscBdPointJacFn`
2198: Level: intermediate
2200: Note:
2201: We are using a first order FEM model for the weak form\:
2203: $$
2204: \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
2205: + \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
2206: $$
2208: Developer Note:
2209: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
2211: .seealso: `PetscDS`, `PetscBdPointJacFn`, `PetscDSSetBdJacobianPreconditioner()`
2212: @*/
2213: PetscErrorCode PetscDSGetBdJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, PetscBdPointJacFn **g0, PetscBdPointJacFn **g1, PetscBdPointJacFn **g2, PetscBdPointJacFn **g3)
2214: {
2215: PetscBdPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
2216: PetscInt n0, n1, n2, n3;
2218: PetscFunctionBegin;
2220: 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);
2221: 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);
2222: PetscCall(PetscWeakFormGetBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2223: *g0 = tmp0 ? tmp0[0] : NULL;
2224: *g1 = tmp1 ? tmp1[0] : NULL;
2225: *g2 = tmp2 ? tmp2[0] : NULL;
2226: *g3 = tmp3 ? tmp3[0] : NULL;
2227: PetscFunctionReturn(PETSC_SUCCESS);
2228: }
2230: /*@
2231: PetscDSSetBdJacobianPreconditioner - Set the pointwise boundary Jacobian preconditioner function for given test and basis field that constructs the
2232: matrix used to construct the preconditioner
2234: Not Collective; No Fortran Support
2236: Input Parameters:
2237: + ds - The `PetscDS`
2238: . f - The test field number
2239: . g - The field number
2240: . g0 - integrand for the test and basis function term, see `PetscBdPointJacFn`
2241: . g1 - integrand for the test function and basis function gradient term, see `PetscBdPointJacFn`
2242: . g2 - integrand for the test function gradient and basis function term, see `PetscBdPointJacFn`
2243: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscBdPointJacFn`
2245: Level: intermediate
2247: Note:
2248: We are using a first order FEM model for the weak form\:
2250: $$
2251: \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
2252: + \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
2253: $$
2255: Developer Note:
2256: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
2258: .seealso: `PetscDS`, `PetscBdPointJacFn`, `PetscDSGetBdJacobianPreconditioner()`
2259: @*/
2260: PetscErrorCode PetscDSSetBdJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, PetscBdPointJacFn *g0, PetscBdPointJacFn *g1, PetscBdPointJacFn *g2, PetscBdPointJacFn *g3)
2261: {
2262: PetscFunctionBegin;
2268: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2269: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2270: PetscCall(PetscWeakFormSetIndexBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2271: PetscFunctionReturn(PETSC_SUCCESS);
2272: }
2274: /*@
2275: PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field
2277: Not Collective
2279: Input Parameters:
2280: + prob - The `PetscDS`
2281: - f - The test field number
2283: Output Parameters:
2284: + sol - exact solution function for the test field, see `PetscPointExactSolutionFn`
2285: - ctx - exact solution context
2287: Level: intermediate
2289: .seealso: `PetscDS`, `PetscPointExactSolutionFn`, `PetscDSSetExactSolution()`, `PetscDSGetExactSolutionTimeDerivative()`
2290: @*/
2291: PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscPointExactSolutionFn **sol, void **ctx)
2292: {
2293: PetscFunctionBegin;
2295: 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);
2296: if (sol) {
2297: PetscAssertPointer(sol, 3);
2298: *sol = prob->exactSol[f];
2299: }
2300: if (ctx) {
2301: PetscAssertPointer(ctx, 4);
2302: *ctx = prob->exactCtx[f];
2303: }
2304: PetscFunctionReturn(PETSC_SUCCESS);
2305: }
2307: /*@
2308: PetscDSSetExactSolution - Set the pointwise exact solution function for a given test field
2310: Not Collective
2312: Input Parameters:
2313: + prob - The `PetscDS`
2314: . f - The test field number
2315: . sol - solution function for the test fields, see `PetscPointExactSolutionFn`
2316: - ctx - solution context or `NULL`
2318: Level: intermediate
2320: .seealso: `PetscDS`, `PetscPointExactSolutionFn`, `PetscDSGetExactSolution()`
2321: @*/
2322: PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscPointExactSolutionFn *sol, PetscCtx ctx)
2323: {
2324: PetscFunctionBegin;
2326: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2327: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
2328: if (sol) {
2330: prob->exactSol[f] = sol;
2331: }
2332: if (ctx) {
2334: prob->exactCtx[f] = ctx;
2335: }
2336: PetscFunctionReturn(PETSC_SUCCESS);
2337: }
2339: /*@
2340: PetscDSGetExactSolutionTimeDerivative - Get the pointwise time derivative of the exact solution function for a given test field
2342: Not Collective
2344: Input Parameters:
2345: + prob - The `PetscDS`
2346: - f - The test field number
2348: Output Parameters:
2349: + sol - time derivative of the exact solution for the test field, see `PetscPointExactSolutionFn`
2350: - ctx - the exact solution context
2352: Level: intermediate
2354: .seealso: `PetscDS`, `PetscPointExactSolutionFn`, `PetscDSSetExactSolutionTimeDerivative()`, `PetscDSGetExactSolution()`
2355: @*/
2356: PetscErrorCode PetscDSGetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscPointExactSolutionFn **sol, void **ctx)
2357: {
2358: PetscFunctionBegin;
2360: 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);
2361: if (sol) {
2362: PetscAssertPointer(sol, 3);
2363: *sol = prob->exactSol_t[f];
2364: }
2365: if (ctx) {
2366: PetscAssertPointer(ctx, 4);
2367: *ctx = prob->exactCtx_t[f];
2368: }
2369: PetscFunctionReturn(PETSC_SUCCESS);
2370: }
2372: /*@
2373: PetscDSSetExactSolutionTimeDerivative - Set the pointwise time derivative of the exact solution function for a given test field
2375: Not Collective
2377: Input Parameters:
2378: + prob - The `PetscDS`
2379: . f - The test field number
2380: . sol - time derivative of the solution function for the test fields, see `PetscPointExactSolutionFn`
2381: - ctx - the solution context or `NULL`
2383: Level: intermediate
2385: .seealso: `PetscDS`, `PetscPointExactSolutionFn`, `PetscDSGetExactSolutionTimeDerivative()`, `PetscDSSetExactSolution()`
2386: @*/
2387: PetscErrorCode PetscDSSetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscPointExactSolutionFn *sol, PetscCtx ctx)
2388: {
2389: PetscFunctionBegin;
2391: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2392: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
2393: if (sol) {
2395: prob->exactSol_t[f] = sol;
2396: }
2397: if (ctx) {
2399: prob->exactCtx_t[f] = ctx;
2400: }
2401: PetscFunctionReturn(PETSC_SUCCESS);
2402: }
2404: /*@
2405: PetscDSGetLowerBound - Get the pointwise lower bound function for a given field
2407: Not Collective
2409: Input Parameters:
2410: + ds - The PetscDS
2411: - f - The field number
2413: Output Parameters:
2414: + lb - lower bound function for the field, see `PetscPointBoundFn`
2415: - ctx - lower bound context that was set with `PetscDSSetLowerBound()`
2417: Level: intermediate
2419: .seealso: `PetscDS`, `PetscPointBoundFn`, `PetscDSSetLowerBound()`, `PetscDSGetUpperBound()`, `PetscDSGetExactSolution()`
2420: @*/
2421: PetscErrorCode PetscDSGetLowerBound(PetscDS ds, PetscInt f, PetscPointBoundFn **lb, void **ctx)
2422: {
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: if (lb) {
2427: PetscAssertPointer(lb, 3);
2428: *lb = ds->lowerBound[f];
2429: }
2430: if (ctx) {
2431: PetscAssertPointer(ctx, 4);
2432: *ctx = ds->lowerCtx[f];
2433: }
2434: PetscFunctionReturn(PETSC_SUCCESS);
2435: }
2437: /*@
2438: PetscDSSetLowerBound - Set the pointwise lower bound function for a given field
2440: Not Collective
2442: Input Parameters:
2443: + ds - The `PetscDS`
2444: . f - The field number
2445: . lb - lower bound function for the test fields, see `PetscPointBoundFn`
2446: - ctx - lower bound context or `NULL` which will be passed to `lb`
2448: Level: intermediate
2450: .seealso: `PetscDS`, `PetscPointBoundFn`, `PetscDSGetLowerBound()`, `PetscDSGetUpperBound()`, `PetscDSGetExactSolution()`
2451: @*/
2452: PetscErrorCode PetscDSSetLowerBound(PetscDS ds, PetscInt f, PetscPointBoundFn *lb, PetscCtx ctx)
2453: {
2454: PetscFunctionBegin;
2456: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2457: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
2458: if (lb) {
2460: ds->lowerBound[f] = lb;
2461: }
2462: if (ctx) {
2464: ds->lowerCtx[f] = ctx;
2465: }
2466: PetscFunctionReturn(PETSC_SUCCESS);
2467: }
2469: /*@
2470: PetscDSGetUpperBound - Get the pointwise upper bound function for a given field
2472: Not Collective
2474: Input Parameters:
2475: + ds - The `PetscDS`
2476: - f - The field number
2478: Output Parameters:
2479: + ub - upper bound function for the field, see `PetscPointBoundFn`
2480: - ctx - upper bound context that was set with `PetscDSSetUpperBound()`
2482: Level: intermediate
2484: .seealso: `PetscDS`, `PetscPointBoundFn`, `PetscDSSetUpperBound()`, `PetscDSGetLowerBound()`, `PetscDSGetExactSolution()`
2485: @*/
2486: PetscErrorCode PetscDSGetUpperBound(PetscDS ds, PetscInt f, PetscPointBoundFn **ub, void **ctx)
2487: {
2488: PetscFunctionBegin;
2490: 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);
2491: if (ub) {
2492: PetscAssertPointer(ub, 3);
2493: *ub = ds->upperBound[f];
2494: }
2495: if (ctx) {
2496: PetscAssertPointer(ctx, 4);
2497: *ctx = ds->upperCtx[f];
2498: }
2499: PetscFunctionReturn(PETSC_SUCCESS);
2500: }
2502: /*@
2503: PetscDSSetUpperBound - Set the pointwise upper bound function for a given field
2505: Not Collective
2507: Input Parameters:
2508: + ds - The `PetscDS`
2509: . f - The field number
2510: . ub - upper bound function for the test fields, see `PetscPointBoundFn`
2511: - ctx - context or `NULL` that will be passed to `ub`
2513: Level: intermediate
2515: .seealso: `PetscDS`, `PetscPointBoundFn`, `PetscDSGetUpperBound()`, `PetscDSGetLowerBound()`, `PetscDSGetExactSolution()`
2516: @*/
2517: PetscErrorCode PetscDSSetUpperBound(PetscDS ds, PetscInt f, PetscPointBoundFn *ub, PetscCtx ctx)
2518: {
2519: PetscFunctionBegin;
2521: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2522: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
2523: if (ub) {
2525: ds->upperBound[f] = ub;
2526: }
2527: if (ctx) {
2529: ds->upperCtx[f] = ctx;
2530: }
2531: PetscFunctionReturn(PETSC_SUCCESS);
2532: }
2534: /*@
2535: PetscDSGetConstants - Returns the array of constants passed to point functions from a `PetscDS` object
2537: Not Collective
2539: Input Parameter:
2540: . ds - The `PetscDS` object
2542: Output Parameters:
2543: + numConstants - The number of constants, or pass in `NULL` if not required
2544: - constants - The array of constants, `NULL` if there are none
2546: Level: intermediate
2548: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSCreate()`
2549: @*/
2550: PetscErrorCode PetscDSGetConstants(PetscDS ds, PeOp PetscInt *numConstants, PeOp const PetscScalar *constants[])
2551: {
2552: PetscFunctionBegin;
2554: if (numConstants) {
2555: PetscAssertPointer(numConstants, 2);
2556: *numConstants = ds->numConstants;
2557: }
2558: if (constants) {
2559: PetscAssertPointer(constants, 3);
2560: *constants = ds->constants;
2561: }
2562: PetscFunctionReturn(PETSC_SUCCESS);
2563: }
2565: /*@
2566: PetscDSSetConstants - Set the array of constants passed to point functions from a `PetscDS`
2568: Not Collective
2570: Input Parameters:
2571: + ds - The `PetscDS` object
2572: . numConstants - The number of constants
2573: - constants - The array of constants, `NULL` if there are none
2575: Level: intermediate
2577: .seealso: `PetscDS`, `PetscDSGetConstants()`, `PetscDSCreate()`
2578: @*/
2579: PetscErrorCode PetscDSSetConstants(PetscDS ds, PetscInt numConstants, PetscScalar constants[])
2580: {
2581: PetscFunctionBegin;
2583: if (numConstants != ds->numConstants) {
2584: PetscCall(PetscFree(ds->constants));
2585: ds->numConstants = numConstants;
2586: PetscCall(PetscMalloc1(ds->numConstants + ds->numFuncConstants, &ds->constants));
2587: }
2588: if (ds->numConstants) {
2589: PetscAssertPointer(constants, 3);
2590: PetscCall(PetscArraycpy(ds->constants, constants, ds->numConstants));
2591: }
2592: PetscFunctionReturn(PETSC_SUCCESS);
2593: }
2595: /*@
2596: PetscDSSetIntegrationParameters - Set the parameters for a particular integration
2598: Not Collective
2600: Input Parameters:
2601: + ds - The `PetscDS` object
2602: . fieldI - The test field for a given point function, or `PETSC_DETERMINE`
2603: - fieldJ - The basis field for a given point function, or `PETSC_DETERMINE`
2605: Level: intermediate
2607: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSGetConstants()`, `PetscDSCreate()`
2608: @*/
2609: PetscErrorCode PetscDSSetIntegrationParameters(PetscDS ds, PetscInt fieldI, PetscInt fieldJ)
2610: {
2611: PetscFunctionBegin;
2613: ds->constants[ds->numConstants] = fieldI;
2614: ds->constants[ds->numConstants + 1] = fieldJ;
2615: PetscFunctionReturn(PETSC_SUCCESS);
2616: }
2618: /*@
2619: PetscDSSetCellParameters - Set the parameters for a particular cell
2621: Not Collective
2623: Input Parameters:
2624: + ds - The `PetscDS` object
2625: - volume - The cell volume
2627: Level: intermediate
2629: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSGetConstants()`, `PetscDSCreate()`
2630: @*/
2631: PetscErrorCode PetscDSSetCellParameters(PetscDS ds, PetscReal volume)
2632: {
2633: PetscFunctionBegin;
2635: ds->constants[ds->numConstants + 2] = volume;
2636: PetscFunctionReturn(PETSC_SUCCESS);
2637: }
2639: /*@
2640: PetscDSGetFieldIndex - Returns the index of the given field
2642: Not Collective
2644: Input Parameters:
2645: + prob - The `PetscDS` object
2646: - disc - The discretization object
2648: Output Parameter:
2649: . f - The field number
2651: Level: beginner
2653: .seealso: `PetscDS`, `PetscGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2654: @*/
2655: PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f)
2656: {
2657: PetscInt g;
2659: PetscFunctionBegin;
2661: PetscAssertPointer(f, 3);
2662: *f = -1;
2663: for (g = 0; g < prob->Nf; ++g) {
2664: if (disc == prob->disc[g]) break;
2665: }
2666: PetscCheck(g != prob->Nf, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS.");
2667: *f = g;
2668: PetscFunctionReturn(PETSC_SUCCESS);
2669: }
2671: /*@
2672: PetscDSGetFieldSize - Returns the size of the given field in the full space basis
2674: Not Collective
2676: Input Parameters:
2677: + prob - The `PetscDS` object
2678: - f - The field number
2680: Output Parameter:
2681: . size - The size
2683: Level: beginner
2685: .seealso: `PetscDS`, `PetscDSGetFieldOffset()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2686: @*/
2687: PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size)
2688: {
2689: PetscFunctionBegin;
2691: PetscAssertPointer(size, 3);
2692: 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);
2693: PetscCall(PetscDSSetUp(prob));
2694: *size = prob->Nb[f];
2695: PetscFunctionReturn(PETSC_SUCCESS);
2696: }
2698: /*@
2699: PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis
2701: Not Collective
2703: Input Parameters:
2704: + prob - The `PetscDS` object
2705: - f - The field number
2707: Output Parameter:
2708: . off - The offset
2710: Level: beginner
2712: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2713: @*/
2714: PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off)
2715: {
2716: PetscInt size;
2718: PetscFunctionBegin;
2720: PetscAssertPointer(off, 3);
2721: 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);
2722: *off = 0;
2723: for (PetscInt g = 0; g < f; ++g) {
2724: PetscCall(PetscDSGetFieldSize(prob, g, &size));
2725: *off += size;
2726: }
2727: PetscFunctionReturn(PETSC_SUCCESS);
2728: }
2730: /*@
2731: PetscDSGetFieldOffsetCohesive - Returns the offset of the given field in the full space basis on a cohesive cell
2733: Not Collective
2735: Input Parameters:
2736: + ds - The `PetscDS` object
2737: - f - The field number
2739: Output Parameter:
2740: . off - The offset
2742: Level: beginner
2744: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2745: @*/
2746: PetscErrorCode PetscDSGetFieldOffsetCohesive(PetscDS ds, PetscInt f, PetscInt *off)
2747: {
2748: PetscInt size;
2750: PetscFunctionBegin;
2752: PetscAssertPointer(off, 3);
2753: 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);
2754: *off = 0;
2755: for (PetscInt g = 0; g < f; ++g) {
2756: PetscBool cohesive;
2758: PetscCall(PetscDSGetCohesive(ds, g, &cohesive));
2759: PetscCall(PetscDSGetFieldSize(ds, g, &size));
2760: *off += cohesive ? size : size * 2;
2761: }
2762: PetscFunctionReturn(PETSC_SUCCESS);
2763: }
2765: /*@
2766: PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point
2768: Not Collective
2770: Input Parameter:
2771: . prob - The `PetscDS` object
2773: Output Parameter:
2774: . dimensions - The number of dimensions
2776: Level: beginner
2778: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2779: @*/
2780: PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[])
2781: {
2782: PetscFunctionBegin;
2784: PetscCall(PetscDSSetUp(prob));
2785: PetscAssertPointer(dimensions, 2);
2786: *dimensions = prob->Nb;
2787: PetscFunctionReturn(PETSC_SUCCESS);
2788: }
2790: /*@
2791: PetscDSGetComponents - Returns the number of components for each field on an evaluation point
2793: Not Collective
2795: Input Parameter:
2796: . prob - The `PetscDS` object
2798: Output Parameter:
2799: . components - The number of components
2801: Level: beginner
2803: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2804: @*/
2805: PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[])
2806: {
2807: PetscFunctionBegin;
2809: PetscCall(PetscDSSetUp(prob));
2810: PetscAssertPointer(components, 2);
2811: *components = prob->Nc;
2812: PetscFunctionReturn(PETSC_SUCCESS);
2813: }
2815: /*@
2816: PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point
2818: Not Collective
2820: Input Parameters:
2821: + prob - The `PetscDS` object
2822: - f - The field number
2824: Output Parameter:
2825: . off - The offset
2827: Level: beginner
2829: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2830: @*/
2831: PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off)
2832: {
2833: PetscFunctionBegin;
2835: PetscAssertPointer(off, 3);
2836: 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);
2837: PetscCall(PetscDSSetUp(prob));
2838: *off = prob->off[f];
2839: PetscFunctionReturn(PETSC_SUCCESS);
2840: }
2842: /*@
2843: PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point
2845: Not Collective
2847: Input Parameter:
2848: . prob - The `PetscDS` object
2850: Output Parameter:
2851: . offsets - The offsets
2853: Level: beginner
2855: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2856: @*/
2857: PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[])
2858: {
2859: PetscFunctionBegin;
2861: PetscAssertPointer(offsets, 2);
2862: PetscCall(PetscDSSetUp(prob));
2863: *offsets = prob->off;
2864: PetscFunctionReturn(PETSC_SUCCESS);
2865: }
2867: /*@
2868: PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point
2870: Not Collective
2872: Input Parameter:
2873: . prob - The `PetscDS` object
2875: Output Parameter:
2876: . offsets - The offsets
2878: Level: beginner
2880: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2881: @*/
2882: PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[])
2883: {
2884: PetscFunctionBegin;
2886: PetscAssertPointer(offsets, 2);
2887: PetscCall(PetscDSSetUp(prob));
2888: *offsets = prob->offDer;
2889: PetscFunctionReturn(PETSC_SUCCESS);
2890: }
2892: /*@
2893: PetscDSGetComponentOffsetsCohesive - Returns the offset of each field on an evaluation point
2895: Not Collective
2897: Input Parameters:
2898: + ds - The `PetscDS` object
2899: - s - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive
2901: Output Parameter:
2902: . offsets - The offsets
2904: Level: beginner
2906: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2907: @*/
2908: PetscErrorCode PetscDSGetComponentOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
2909: {
2910: PetscFunctionBegin;
2912: PetscAssertPointer(offsets, 3);
2913: PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS");
2914: PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s);
2915: PetscCall(PetscDSSetUp(ds));
2916: *offsets = ds->offCohesive[s];
2917: PetscFunctionReturn(PETSC_SUCCESS);
2918: }
2920: /*@
2921: PetscDSGetComponentDerivativeOffsetsCohesive - Returns the offset of each field derivative on an evaluation point
2923: Not Collective
2925: Input Parameters:
2926: + ds - The `PetscDS` object
2927: - s - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive
2929: Output Parameter:
2930: . offsets - The offsets
2932: Level: beginner
2934: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2935: @*/
2936: PetscErrorCode PetscDSGetComponentDerivativeOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
2937: {
2938: PetscFunctionBegin;
2940: PetscAssertPointer(offsets, 3);
2941: PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS");
2942: PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s);
2943: PetscCall(PetscDSSetUp(ds));
2944: *offsets = ds->offDerCohesive[s];
2945: PetscFunctionReturn(PETSC_SUCCESS);
2946: }
2948: /*@
2949: PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization
2951: Not Collective
2953: Input Parameter:
2954: . prob - The `PetscDS` object
2956: Output Parameter:
2957: . T - The basis function and derivatives tabulation at quadrature points for each field, see `PetscTabulation` for its details
2959: Level: intermediate
2961: Note:
2962: The tabulation is only valid so long as the `PetscDS` has not be destroyed. There is no `PetscDSRestoreTabulation()` in C.
2964: Fortran Note:
2965: Use the declaration
2966: .vb
2967: PetscTabulation, pointer :: tab(:)
2968: .ve
2969: and access the values using, for example,
2970: .vb
2971: tab(i)%ptr%K
2972: tab(i)%ptr%T(j)%ptr
2973: .ve
2974: where $ i = 1, 2, ..., Nf $ and $ j = 1, 2, ..., tab(i)%ptr%K+1 $.
2976: Use `PetscDSRestoreTabulation()` to restore the array
2978: Developer Note:
2979: The Fortran language syntax does not directly support arrays of pointers, the '%ptr' notation allows mimicking their use in Fortran.
2981: .seealso: `PetscDS`, `PetscTabulation`, `PetscDSCreate()`
2982: @*/
2983: PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscTabulation *T[]) PeNS
2984: {
2985: PetscFunctionBegin;
2987: PetscAssertPointer(T, 2);
2988: PetscCall(PetscDSSetUp(prob));
2989: *T = prob->T;
2990: PetscFunctionReturn(PETSC_SUCCESS);
2991: }
2993: /*@
2994: PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces
2996: Not Collective
2998: Input Parameter:
2999: . prob - The `PetscDS` object
3001: Output Parameter:
3002: . Tf - The basis function and derivative tabulation on each local face at quadrature points for each field
3004: Level: intermediate
3006: Note:
3007: The tabulation is only valid so long as the `PetscDS` has not be destroyed. There is no `PetscDSRestoreFaceTabulation()` in C.
3009: .seealso: `PetscTabulation`, `PetscDS`, `PetscDSGetTabulation()`, `PetscDSCreate()`
3010: @*/
3011: PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscTabulation *Tf[])
3012: {
3013: PetscFunctionBegin;
3015: PetscAssertPointer(Tf, 2);
3016: PetscCall(PetscDSSetUp(prob));
3017: *Tf = prob->Tf;
3018: PetscFunctionReturn(PETSC_SUCCESS);
3019: }
3021: /*@
3022: PetscDSGetEvaluationArrays - Get scratch arrays used to evaluate fields, time derivatives, and field gradients at quadrature points.
3024: Not Collective
3026: Input Parameter:
3027: . prob - the `PetscDS`
3029: Output Parameters:
3030: + u - array for the field values, or `NULL` if not needed
3031: . u_t - array for the field time derivatives, or `NULL` if not needed
3032: - u_x - array for the field gradients, or `NULL` if not needed
3034: Level: developer
3036: Note:
3037: The returned arrays are owned by the `PetscDS` and must not be freed by the caller.
3039: .seealso: `PetscDS`, `PetscDSGetWeakFormArrays()`, `PetscDSGetWorkspace()`
3040: @*/
3041: PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar *u[], PetscScalar *u_t[], PetscScalar *u_x[])
3042: {
3043: PetscFunctionBegin;
3045: PetscCall(PetscDSSetUp(prob));
3046: if (u) {
3047: PetscAssertPointer(u, 2);
3048: *u = prob->u;
3049: }
3050: if (u_t) {
3051: PetscAssertPointer(u_t, 3);
3052: *u_t = prob->u_t;
3053: }
3054: if (u_x) {
3055: PetscAssertPointer(u_x, 4);
3056: *u_x = prob->u_x;
3057: }
3058: PetscFunctionReturn(PETSC_SUCCESS);
3059: }
3061: PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar *f0[], PetscScalar *f1[], PetscScalar *g0[], PetscScalar *g1[], PetscScalar *g2[], PetscScalar *g3[])
3062: {
3063: PetscFunctionBegin;
3065: PetscCall(PetscDSSetUp(prob));
3066: if (f0) {
3067: PetscAssertPointer(f0, 2);
3068: *f0 = prob->f0;
3069: }
3070: if (f1) {
3071: PetscAssertPointer(f1, 3);
3072: *f1 = prob->f1;
3073: }
3074: if (g0) {
3075: PetscAssertPointer(g0, 4);
3076: *g0 = prob->g0;
3077: }
3078: if (g1) {
3079: PetscAssertPointer(g1, 5);
3080: *g1 = prob->g1;
3081: }
3082: if (g2) {
3083: PetscAssertPointer(g2, 6);
3084: *g2 = prob->g2;
3085: }
3086: if (g3) {
3087: PetscAssertPointer(g3, 7);
3088: *g3 = prob->g3;
3089: }
3090: PetscFunctionReturn(PETSC_SUCCESS);
3091: }
3093: /*@
3094: PetscDSGetWorkspace - Get scratch storage used during discretization computations.
3096: Not Collective
3098: Input Parameter:
3099: . prob - the `PetscDS`
3101: Output Parameters:
3102: + x - array for real-valued quadrature point coordinates, or `NULL` if not needed
3103: . basisReal - array for the real-valued basis function values, or `NULL` if not needed
3104: . basisDerReal - array for the real-valued basis function derivatives, or `NULL` if not needed
3105: . testReal - array for the real-valued test function values, or `NULL` if not needed
3106: - testDerReal - array for the real-valued test function derivatives, or `NULL` if not needed
3108: Level: developer
3110: Note:
3111: The returned arrays are owned by the `PetscDS` and must not be freed by the caller.
3113: .seealso: `PetscDS`, `PetscDSGetEvaluationArrays()`, `PetscDSGetWeakFormArrays()`
3114: @*/
3115: PetscErrorCode PetscDSGetWorkspace(PetscDS prob, PetscReal **x, PetscScalar **basisReal, PetscScalar **basisDerReal, PetscScalar **testReal, PetscScalar **testDerReal)
3116: {
3117: PetscFunctionBegin;
3119: PetscCall(PetscDSSetUp(prob));
3120: if (x) {
3121: PetscAssertPointer(x, 2);
3122: *x = prob->x;
3123: }
3124: if (basisReal) {
3125: PetscAssertPointer(basisReal, 3);
3126: *basisReal = prob->basisReal;
3127: }
3128: if (basisDerReal) {
3129: PetscAssertPointer(basisDerReal, 4);
3130: *basisDerReal = prob->basisDerReal;
3131: }
3132: if (testReal) {
3133: PetscAssertPointer(testReal, 5);
3134: *testReal = prob->testReal;
3135: }
3136: if (testDerReal) {
3137: PetscAssertPointer(testDerReal, 6);
3138: *testDerReal = prob->testDerReal;
3139: }
3140: PetscFunctionReturn(PETSC_SUCCESS);
3141: }
3143: /*@
3144: PetscDSAddBoundary - Add a boundary condition to the model.
3146: Collective
3148: Input Parameters:
3149: + ds - The `PetscDS` object
3150: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3151: . name - The name for the boundary condition
3152: . label - The label defining constrained points
3153: . Nv - The number of `DMLabel` values for constrained points
3154: . values - An array of label values for constrained points
3155: . field - The field to constrain
3156: . Nc - The number of constrained field components (0 will constrain all fields)
3157: . comps - An array of constrained component numbers
3158: . bcFunc - A pointwise function giving boundary values
3159: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or `NULL`
3160: - ctx - An optional application context for `bcFunc`
3162: Output Parameter:
3163: . bd - The boundary number
3165: Options Database Keys:
3166: + -bc_NAME values - comma separated list of values for the boundary condition NAME
3167: - -bc_NAME_comp comps - comma separated list of components for the boundary condition NAME
3169: Level: developer
3171: Note:
3172: Both `bcFunc` and `bcFunc_t` will depend on the boundary condition type. If the type if `DM_BC_ESSENTIAL`, then the calling sequence is\:
3173: .vb
3174: void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])
3175: .ve
3177: If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value, then the calling sequence is\:
3178: .vb
3179: void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3180: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3181: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3182: PetscReal time, const PetscReal x[], PetscScalar bcval[])
3183: .ve
3184: + dim - the coordinate dimension
3185: . Nf - the number of fields
3186: . uOff - the offset into u[] and u_t[] for each field
3187: . uOff_x - the offset into u_x[] for each field
3188: . u - each field evaluated at the current point
3189: . u_t - the time derivative of each field evaluated at the current point
3190: . u_x - the gradient of each field evaluated at the current point
3191: . aOff - the offset into a[] and a_t[] for each auxiliary field
3192: . aOff_x - the offset into a_x[] for each auxiliary field
3193: . a - each auxiliary field evaluated at the current point
3194: . a_t - the time derivative of each auxiliary field evaluated at the current point
3195: . a_x - the gradient of auxiliary each field evaluated at the current point
3196: . t - current time
3197: . x - coordinates of the current point
3198: . numConstants - number of constant parameters
3199: . constants - constant parameters
3200: - bcval - output values at the current point
3202: Notes:
3203: The pointwise functions are used to provide boundary values for essential boundary
3204: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3205: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3206: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3208: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundaryByName()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3209: @*/
3210: PetscErrorCode PetscDSAddBoundary(PetscDS ds, DMBoundaryConditionType type, const char name[], DMLabel label, PetscInt Nv, const PetscInt values[], PetscInt field, PetscInt Nc, const PetscInt comps[], PetscVoidFn *bcFunc, PetscVoidFn *bcFunc_t, PetscCtx ctx, PetscInt *bd)
3211: {
3212: DSBoundary head = ds->boundary, b;
3213: PetscInt n = 0;
3214: const char *lname;
3216: PetscFunctionBegin;
3219: PetscAssertPointer(name, 3);
3224: 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);
3225: if (Nc > 0) {
3226: PetscInt *fcomps;
3228: PetscCall(PetscDSGetComponents(ds, &fcomps));
3229: 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);
3230: for (PetscInt c = 0; c < Nc; ++c) {
3231: 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);
3232: }
3233: }
3234: PetscCall(PetscNew(&b));
3235: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3236: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf));
3237: PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf));
3238: PetscCall(PetscMalloc1(Nv, &b->values));
3239: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3240: PetscCall(PetscMalloc1(Nc, &b->comps));
3241: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3242: PetscCall(PetscObjectGetName((PetscObject)label, &lname));
3243: PetscCall(PetscStrallocpy(lname, (char **)&b->lname));
3244: b->type = type;
3245: b->label = label;
3246: b->Nv = Nv;
3247: b->field = field;
3248: b->Nc = Nc;
3249: b->func = bcFunc;
3250: b->func_t = bcFunc_t;
3251: b->ctx = ctx;
3252: b->next = NULL;
3253: /* Append to linked list so that we can preserve the order */
3254: if (!head) ds->boundary = b;
3255: while (head) {
3256: if (!head->next) {
3257: head->next = b;
3258: head = b;
3259: }
3260: head = head->next;
3261: ++n;
3262: }
3263: if (bd) {
3264: PetscAssertPointer(bd, 13);
3265: *bd = n;
3266: }
3267: PetscFunctionReturn(PETSC_SUCCESS);
3268: }
3270: // PetscClangLinter pragma ignore: -fdoc-section-header-unknown
3271: /*@
3272: PetscDSAddBoundaryByName - Add a boundary condition to the model.
3274: Collective
3276: Input Parameters:
3277: + ds - The `PetscDS` object
3278: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3279: . name - The boundary condition name
3280: . lname - The name of the label defining constrained points
3281: . Nv - The number of `DMLabel` values for constrained points
3282: . values - An array of label values for constrained points
3283: . field - The field to constrain
3284: . Nc - The number of constrained field components (0 will constrain all fields)
3285: . comps - An array of constrained component numbers
3286: . bcFunc - A pointwise function giving boundary values
3287: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or `NULL`
3288: - ctx - An optional application context for `bcFunc`
3290: Output Parameter:
3291: . bd - The boundary number
3293: Options Database Keys:
3294: + -bc_NAME values - comma separated list of values for the boundary condition NAME
3295: - -bc_NAME_comp comps - comma separated list of components for the boundary condition NAME
3297: Calling Sequence of `bcFunc` and `bcFunc_t`:
3298: If the type is `DM_BC_ESSENTIAL`
3299: .vb
3300: void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])
3301: .ve
3302: If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value,
3303: .vb
3304: void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3305: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3306: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3307: PetscReal time, const PetscReal x[], PetscScalar bcval[])
3308: .ve
3309: + dim - the coordinate dimension
3310: . Nf - the number of fields
3311: . uOff - the offset into `u`[] and `u_t`[] for each field
3312: . uOff_x - the offset into `u_x`[] for each field
3313: . u - each field evaluated at the current point
3314: . u_t - the time derivative of each field evaluated at the current point
3315: . u_x - the gradient of each field evaluated at the current point
3316: . aOff - the offset into `a`[] and `a_t`[] for each auxiliary field
3317: . aOff_x - the offset into `a_x`[] for each auxiliary field
3318: . a - each auxiliary field evaluated at the current point
3319: . a_t - the time derivative of each auxiliary field evaluated at the current point
3320: . a_x - the gradient of auxiliary each field evaluated at the current point
3321: . t - current time
3322: . x - coordinates of the current point
3323: . numConstants - number of constant parameters
3324: . constants - constant parameters
3325: - bcval - output values at the current point
3327: Level: developer
3329: Notes:
3330: The pointwise functions are used to provide boundary values for essential boundary
3331: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3332: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3333: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3335: This function should only be used with `DMFOREST` currently, since labels cannot be defined before the underlying `DMPLEX` is built.
3337: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3338: @*/
3339: PetscErrorCode PetscDSAddBoundaryByName(PetscDS ds, DMBoundaryConditionType type, const char name[], const char lname[], PetscInt Nv, const PetscInt values[], PetscInt field, PetscInt Nc, const PetscInt comps[], PetscVoidFn *bcFunc, PetscVoidFn *bcFunc_t, PetscCtx ctx, PetscInt *bd)
3340: {
3341: DSBoundary head = ds->boundary, b;
3342: PetscInt n = 0;
3344: PetscFunctionBegin;
3347: PetscAssertPointer(name, 3);
3348: PetscAssertPointer(lname, 4);
3352: PetscCall(PetscNew(&b));
3353: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3354: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf));
3355: PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf));
3356: PetscCall(PetscMalloc1(Nv, &b->values));
3357: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3358: PetscCall(PetscMalloc1(Nc, &b->comps));
3359: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3360: PetscCall(PetscStrallocpy(lname, (char **)&b->lname));
3361: b->type = type;
3362: b->label = NULL;
3363: b->Nv = Nv;
3364: b->field = field;
3365: b->Nc = Nc;
3366: b->func = bcFunc;
3367: b->func_t = bcFunc_t;
3368: b->ctx = ctx;
3369: b->next = NULL;
3370: /* Append to linked list so that we can preserve the order */
3371: if (!head) ds->boundary = b;
3372: while (head) {
3373: if (!head->next) {
3374: head->next = b;
3375: head = b;
3376: }
3377: head = head->next;
3378: ++n;
3379: }
3380: if (bd) {
3381: PetscAssertPointer(bd, 13);
3382: *bd = n;
3383: }
3384: PetscFunctionReturn(PETSC_SUCCESS);
3385: }
3387: /*@
3388: PetscDSUpdateBoundary - Change a boundary condition for the model.
3390: Input Parameters:
3391: + ds - The `PetscDS` object
3392: . bd - The boundary condition number
3393: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3394: . name - The boundary condition name
3395: . label - The label defining constrained points
3396: . Nv - The number of `DMLabel` ids for constrained points
3397: . values - An array of ids for constrained points
3398: . field - The field to constrain
3399: . Nc - The number of constrained field components
3400: . comps - An array of constrained component numbers
3401: . bcFunc - A pointwise function giving boundary values
3402: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or `NULL`
3403: - ctx - An optional application context for `bcFunc`
3405: Level: developer
3407: Notes:
3408: The pointwise functions are used to provide boundary values for essential boundary
3409: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3410: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3411: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3413: The boundary condition number is the order in which it was registered. The user can get the number of boundary conditions from `PetscDSGetNumBoundary()`.
3414: See `PetscDSAddBoundary()` for a description of the calling sequences for the callbacks.
3416: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSGetNumBoundary()`, `DMLabel`
3417: @*/
3418: 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[], PetscVoidFn *bcFunc, PetscVoidFn *bcFunc_t, PetscCtx ctx)
3419: {
3420: DSBoundary b = ds->boundary;
3421: PetscInt n = 0;
3423: PetscFunctionBegin;
3425: while (b) {
3426: if (n == bd) break;
3427: b = b->next;
3428: ++n;
3429: }
3430: PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n);
3431: if (name) {
3432: PetscCall(PetscFree(b->name));
3433: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3434: }
3435: b->type = type;
3436: if (label) {
3437: const char *name;
3439: b->label = label;
3440: PetscCall(PetscFree(b->lname));
3441: PetscCall(PetscObjectGetName((PetscObject)label, &name));
3442: PetscCall(PetscStrallocpy(name, (char **)&b->lname));
3443: }
3444: if (Nv >= 0) {
3445: b->Nv = Nv;
3446: PetscCall(PetscFree(b->values));
3447: PetscCall(PetscMalloc1(Nv, &b->values));
3448: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3449: }
3450: if (field >= 0) b->field = field;
3451: if (Nc >= 0) {
3452: b->Nc = Nc;
3453: PetscCall(PetscFree(b->comps));
3454: PetscCall(PetscMalloc1(Nc, &b->comps));
3455: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3456: }
3457: if (bcFunc) b->func = bcFunc;
3458: if (bcFunc_t) b->func_t = bcFunc_t;
3459: if (ctx) b->ctx = ctx;
3460: PetscFunctionReturn(PETSC_SUCCESS);
3461: }
3463: /*@
3464: PetscDSGetNumBoundary - Get the number of registered boundary conditions
3466: Input Parameter:
3467: . ds - The `PetscDS` object
3469: Output Parameter:
3470: . numBd - The number of boundary conditions
3472: Level: intermediate
3474: .seealso: `PetscDS`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`
3475: @*/
3476: PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd)
3477: {
3478: DSBoundary b = ds->boundary;
3480: PetscFunctionBegin;
3482: PetscAssertPointer(numBd, 2);
3483: *numBd = 0;
3484: while (b) {
3485: ++(*numBd);
3486: b = b->next;
3487: }
3488: PetscFunctionReturn(PETSC_SUCCESS);
3489: }
3491: /*@
3492: PetscDSGetBoundary - Gets a boundary condition from the model
3494: Input Parameters:
3495: + ds - The `PetscDS` object
3496: - bd - The boundary condition number
3498: Output Parameters:
3499: + wf - The `PetscWeakForm` holding the pointwise functions
3500: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3501: . name - The boundary condition name
3502: . label - The label defining constrained points
3503: . Nv - The number of `DMLabel` ids for constrained points
3504: . values - An array of ids for constrained points
3505: . field - The field to constrain
3506: . Nc - The number of constrained field components
3507: . comps - An array of constrained component numbers
3508: . func - A pointwise function giving boundary values
3509: . func_t - A pointwise function giving the time derivative of the boundary values
3510: - ctx - An optional application context for `func`
3512: Level: developer
3514: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `DMLabel`
3515: @*/
3516: 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[], PetscVoidFn **func, PetscVoidFn **func_t, PetscCtxRt ctx)
3517: {
3518: DSBoundary b = ds->boundary;
3519: PetscInt n = 0;
3521: PetscFunctionBegin;
3523: while (b) {
3524: if (n == bd) break;
3525: b = b->next;
3526: ++n;
3527: }
3528: PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n);
3529: if (wf) {
3530: PetscAssertPointer(wf, 3);
3531: *wf = b->wf;
3532: }
3533: if (type) {
3534: PetscAssertPointer(type, 4);
3535: *type = b->type;
3536: }
3537: if (name) {
3538: PetscAssertPointer(name, 5);
3539: *name = b->name;
3540: }
3541: if (label) {
3542: PetscAssertPointer(label, 6);
3543: *label = b->label;
3544: }
3545: if (Nv) {
3546: PetscAssertPointer(Nv, 7);
3547: *Nv = b->Nv;
3548: }
3549: if (values) {
3550: PetscAssertPointer(values, 8);
3551: *values = b->values;
3552: }
3553: if (field) {
3554: PetscAssertPointer(field, 9);
3555: *field = b->field;
3556: }
3557: if (Nc) {
3558: PetscAssertPointer(Nc, 10);
3559: *Nc = b->Nc;
3560: }
3561: if (comps) {
3562: PetscAssertPointer(comps, 11);
3563: *comps = b->comps;
3564: }
3565: if (func) {
3566: PetscAssertPointer(func, 12);
3567: *func = b->func;
3568: }
3569: if (func_t) {
3570: PetscAssertPointer(func_t, 13);
3571: *func_t = b->func_t;
3572: }
3573: if (ctx) {
3574: PetscAssertPointer(ctx, 14);
3575: *(void **)ctx = b->ctx;
3576: }
3577: PetscFunctionReturn(PETSC_SUCCESS);
3578: }
3580: /*@
3581: PetscDSUpdateBoundaryLabels - Update `DMLabel` in each boundary condition using the label name and the input `DM`
3583: Not Collective
3585: Input Parameters:
3586: + ds - The source `PetscDS` object
3587: - dm - The `DM` holding labels
3589: Level: intermediate
3591: .seealso: `PetscDS`, `DMBoundary`, `DM`, `PetscDSCopyBoundary()`, `PetscDSCreate()`, `DMGetLabel()`
3592: @*/
3593: PetscErrorCode PetscDSUpdateBoundaryLabels(PetscDS ds, DM dm)
3594: {
3595: DSBoundary b;
3597: PetscFunctionBegin;
3600: for (b = ds->boundary; b; b = b->next) {
3601: if (b->lname) PetscCall(DMGetLabel(dm, b->lname, &b->label));
3602: }
3603: PetscFunctionReturn(PETSC_SUCCESS);
3604: }
3606: static PetscErrorCode DSBoundaryDuplicate_Internal(DSBoundary b, DSBoundary *bNew)
3607: {
3608: PetscFunctionBegin;
3609: PetscCall(PetscNew(bNew));
3610: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &(*bNew)->wf));
3611: PetscCall(PetscWeakFormCopy(b->wf, (*bNew)->wf));
3612: PetscCall(PetscStrallocpy(b->name, (char **)&((*bNew)->name)));
3613: PetscCall(PetscStrallocpy(b->lname, (char **)&((*bNew)->lname)));
3614: (*bNew)->type = b->type;
3615: (*bNew)->label = b->label;
3616: (*bNew)->Nv = b->Nv;
3617: PetscCall(PetscMalloc1(b->Nv, &(*bNew)->values));
3618: PetscCall(PetscArraycpy((*bNew)->values, b->values, b->Nv));
3619: (*bNew)->field = b->field;
3620: (*bNew)->Nc = b->Nc;
3621: PetscCall(PetscMalloc1(b->Nc, &(*bNew)->comps));
3622: PetscCall(PetscArraycpy((*bNew)->comps, b->comps, b->Nc));
3623: (*bNew)->func = b->func;
3624: (*bNew)->func_t = b->func_t;
3625: (*bNew)->ctx = b->ctx;
3626: PetscFunctionReturn(PETSC_SUCCESS);
3627: }
3629: /*@
3630: PetscDSCopyBoundary - Copy all boundary condition objects to the new `PetscDS`
3632: Not Collective
3634: Input Parameters:
3635: + ds - The source `PetscDS` object
3636: . numFields - The number of selected fields, or `PETSC_DEFAULT` for all fields
3637: - fields - The selected fields, or `NULL` for all fields
3639: Output Parameter:
3640: . newds - The target `PetscDS`, now with a copy of the boundary conditions
3642: Level: intermediate
3644: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3645: @*/
3646: PetscErrorCode PetscDSCopyBoundary(PetscDS ds, PetscInt numFields, const PetscInt fields[], PetscDS newds)
3647: {
3648: DSBoundary b, *lastnext;
3650: PetscFunctionBegin;
3653: if (ds == newds) PetscFunctionReturn(PETSC_SUCCESS);
3654: PetscCall(PetscDSDestroyBoundary(newds));
3655: lastnext = &newds->boundary;
3656: for (b = ds->boundary; b; b = b->next) {
3657: DSBoundary bNew;
3658: PetscInt fieldNew = -1;
3660: if (numFields > 0 && fields) {
3661: PetscInt f;
3663: for (f = 0; f < numFields; ++f)
3664: if (b->field == fields[f]) break;
3665: if (f == numFields) continue;
3666: fieldNew = f;
3667: }
3668: PetscCall(DSBoundaryDuplicate_Internal(b, &bNew));
3669: bNew->field = fieldNew < 0 ? b->field : fieldNew;
3670: *lastnext = bNew;
3671: lastnext = &bNew->next;
3672: }
3673: PetscFunctionReturn(PETSC_SUCCESS);
3674: }
3676: /*@
3677: PetscDSDestroyBoundary - Remove all `DMBoundary` objects from the `PetscDS`
3679: Not Collective
3681: Input Parameter:
3682: . ds - The `PetscDS` object
3684: Level: intermediate
3686: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`
3687: @*/
3688: PetscErrorCode PetscDSDestroyBoundary(PetscDS ds)
3689: {
3690: DSBoundary next = ds->boundary;
3692: PetscFunctionBegin;
3693: while (next) {
3694: DSBoundary b = next;
3696: next = b->next;
3697: PetscCall(PetscWeakFormDestroy(&b->wf));
3698: PetscCall(PetscFree(b->name));
3699: PetscCall(PetscFree(b->lname));
3700: PetscCall(PetscFree(b->values));
3701: PetscCall(PetscFree(b->comps));
3702: PetscCall(PetscFree(b));
3703: }
3704: PetscFunctionReturn(PETSC_SUCCESS);
3705: }
3707: /*@
3708: PetscDSSelectDiscretizations - Copy discretizations to the new `PetscDS` with different field layout
3710: Not Collective
3712: Input Parameters:
3713: + prob - The `PetscDS` object
3714: . numFields - Number of new fields
3715: . fields - Old field number for each new field
3716: . minDegree - Minimum degree for a discretization, or `PETSC_DETERMINE` for no limit
3717: - maxDegree - Maximum degree for a discretization, or `PETSC_DETERMINE` for no limit
3719: Output Parameter:
3720: . newprob - The `PetscDS` copy
3722: Level: intermediate
3724: .seealso: `PetscDS`, `PetscDSSelectEquations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3725: @*/
3726: PetscErrorCode PetscDSSelectDiscretizations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscInt minDegree, PetscInt maxDegree, PetscDS newprob)
3727: {
3728: PetscInt Nf, Nfn, fn;
3730: PetscFunctionBegin;
3732: if (fields) PetscAssertPointer(fields, 3);
3734: PetscCall(PetscDSGetNumFields(prob, &Nf));
3735: PetscCall(PetscDSGetNumFields(newprob, &Nfn));
3736: numFields = numFields < 0 ? Nf : numFields;
3737: for (fn = 0; fn < numFields; ++fn) {
3738: const PetscInt f = fields ? fields[fn] : fn;
3739: PetscObject disc;
3740: PetscClassId id;
3742: if (f >= Nf) continue;
3743: PetscCall(PetscDSGetDiscretization(prob, f, &disc));
3744: PetscCallContinue(PetscObjectGetClassId(disc, &id));
3745: if (id == PETSCFE_CLASSID) {
3746: PetscFE fe;
3748: PetscCall(PetscFELimitDegree((PetscFE)disc, minDegree, maxDegree, &fe));
3749: PetscCall(PetscDSSetDiscretization(newprob, fn, (PetscObject)fe));
3750: PetscCall(PetscFEDestroy(&fe));
3751: } else {
3752: PetscCall(PetscDSSetDiscretization(newprob, fn, disc));
3753: }
3754: }
3755: PetscFunctionReturn(PETSC_SUCCESS);
3756: }
3758: /*@
3759: PetscDSSelectEquations - Copy pointwise function pointers to the new `PetscDS` with different field layout
3761: Not Collective
3763: Input Parameters:
3764: + prob - The `PetscDS` object
3765: . numFields - Number of new fields
3766: - fields - Old field number for each new field
3768: Output Parameter:
3769: . newprob - The `PetscDS` copy
3771: Level: intermediate
3773: .seealso: `PetscDS`, `PetscDSSelectDiscretizations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3774: @*/
3775: PetscErrorCode PetscDSSelectEquations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
3776: {
3777: PetscInt Nf, Nfn, fn, gn;
3779: PetscFunctionBegin;
3781: if (fields) PetscAssertPointer(fields, 3);
3783: PetscCall(PetscDSGetNumFields(prob, &Nf));
3784: PetscCall(PetscDSGetNumFields(newprob, &Nfn));
3785: PetscCheck(numFields <= Nfn, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_SIZ, "Number of fields %" PetscInt_FMT " to transfer must not be greater than the total number of fields %" PetscInt_FMT, numFields, Nfn);
3786: for (fn = 0; fn < numFields; ++fn) {
3787: const PetscInt f = fields ? fields[fn] : fn;
3788: PetscPointFn *obj;
3789: PetscPointFn *f0, *f1;
3790: PetscBdPointFn *f0Bd, *f1Bd;
3791: PetscRiemannFn *r;
3793: if (f >= Nf) continue;
3794: PetscCall(PetscDSGetObjective(prob, f, &obj));
3795: PetscCall(PetscDSGetResidual(prob, f, &f0, &f1));
3796: PetscCall(PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd));
3797: PetscCall(PetscDSGetRiemannSolver(prob, f, &r));
3798: PetscCall(PetscDSSetObjective(newprob, fn, obj));
3799: PetscCall(PetscDSSetResidual(newprob, fn, f0, f1));
3800: PetscCall(PetscDSSetBdResidual(newprob, fn, f0Bd, f1Bd));
3801: PetscCall(PetscDSSetRiemannSolver(newprob, fn, r));
3802: for (gn = 0; gn < numFields; ++gn) {
3803: const PetscInt g = fields ? fields[gn] : gn;
3804: PetscPointJacFn *g0, *g1, *g2, *g3;
3805: PetscPointJacFn *g0p, *g1p, *g2p, *g3p;
3806: PetscBdPointJacFn *g0Bd, *g1Bd, *g2Bd, *g3Bd;
3808: if (g >= Nf) continue;
3809: PetscCall(PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3));
3810: PetscCall(PetscDSGetJacobianPreconditioner(prob, f, g, &g0p, &g1p, &g2p, &g3p));
3811: PetscCall(PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd));
3812: PetscCall(PetscDSSetJacobian(newprob, fn, gn, g0, g1, g2, g3));
3813: PetscCall(PetscDSSetJacobianPreconditioner(newprob, fn, gn, g0p, g1p, g2p, g3p));
3814: PetscCall(PetscDSSetBdJacobian(newprob, fn, gn, g0Bd, g1Bd, g2Bd, g3Bd));
3815: }
3816: }
3817: PetscFunctionReturn(PETSC_SUCCESS);
3818: }
3820: /*@
3821: PetscDSCopyEquations - Copy all pointwise function pointers to another `PetscDS`
3823: Not Collective
3825: Input Parameter:
3826: . prob - The `PetscDS` object
3828: Output Parameter:
3829: . newprob - The `PetscDS` copy
3831: Level: intermediate
3833: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3834: @*/
3835: PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob)
3836: {
3837: PetscWeakForm wf, newwf;
3838: PetscInt Nf, Ng;
3840: PetscFunctionBegin;
3843: PetscCall(PetscDSGetNumFields(prob, &Nf));
3844: PetscCall(PetscDSGetNumFields(newprob, &Ng));
3845: PetscCheck(Nf == Ng, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %" PetscInt_FMT " != %" PetscInt_FMT, Nf, Ng);
3846: PetscCall(PetscDSGetWeakForm(prob, &wf));
3847: PetscCall(PetscDSGetWeakForm(newprob, &newwf));
3848: PetscCall(PetscWeakFormCopy(wf, newwf));
3849: PetscFunctionReturn(PETSC_SUCCESS);
3850: }
3852: /*@
3853: PetscDSCopyConstants - Copy all constants set with `PetscDSSetConstants()` to another `PetscDS`
3855: Not Collective
3857: Input Parameter:
3858: . prob - The `PetscDS` object
3860: Output Parameter:
3861: . newprob - The `PetscDS` copy
3863: Level: intermediate
3865: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3866: @*/
3867: PetscErrorCode PetscDSCopyConstants(PetscDS prob, PetscDS newprob)
3868: {
3869: PetscInt Nc;
3870: const PetscScalar *constants;
3872: PetscFunctionBegin;
3875: PetscCall(PetscDSGetConstants(prob, &Nc, &constants));
3876: PetscCall(PetscDSSetConstants(newprob, Nc, (PetscScalar *)constants));
3877: PetscFunctionReturn(PETSC_SUCCESS);
3878: }
3880: /*@
3881: PetscDSCopyExactSolutions - Copy all exact solutions set with `PetscDSSetExactSolution()` and `PetscDSSetExactSolutionTimeDerivative()` to another `PetscDS`
3883: Not Collective
3885: Input Parameter:
3886: . ds - The `PetscDS` object
3888: Output Parameter:
3889: . newds - The `PetscDS` copy
3891: Level: intermediate
3893: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSCopyBounds()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3894: @*/
3895: PetscErrorCode PetscDSCopyExactSolutions(PetscDS ds, PetscDS newds)
3896: {
3897: PetscSimplePointFn *sol;
3898: void *ctx;
3899: PetscInt Nf, f;
3901: PetscFunctionBegin;
3904: PetscCall(PetscDSGetNumFields(ds, &Nf));
3905: for (f = 0; f < Nf; ++f) {
3906: PetscCall(PetscDSGetExactSolution(ds, f, &sol, &ctx));
3907: PetscCall(PetscDSSetExactSolution(newds, f, sol, ctx));
3908: PetscCall(PetscDSGetExactSolutionTimeDerivative(ds, f, &sol, &ctx));
3909: PetscCall(PetscDSSetExactSolutionTimeDerivative(newds, f, sol, ctx));
3910: }
3911: PetscFunctionReturn(PETSC_SUCCESS);
3912: }
3914: /*@
3915: PetscDSCopyBounds - Copy lower and upper solution bounds set with `PetscDSSetLowerBound()` and `PetscDSSetLowerBound()` to another `PetscDS`
3917: Not Collective
3919: Input Parameter:
3920: . ds - The `PetscDS` object
3922: Output Parameter:
3923: . newds - The `PetscDS` copy
3925: Level: intermediate
3927: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSCopyExactSolutions()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3928: @*/
3929: PetscErrorCode PetscDSCopyBounds(PetscDS ds, PetscDS newds)
3930: {
3931: PetscSimplePointFn *bound;
3932: void *ctx;
3933: PetscInt Nf, f;
3935: PetscFunctionBegin;
3938: PetscCall(PetscDSGetNumFields(ds, &Nf));
3939: for (f = 0; f < Nf; ++f) {
3940: PetscCall(PetscDSGetLowerBound(ds, f, &bound, &ctx));
3941: PetscCall(PetscDSSetLowerBound(newds, f, bound, ctx));
3942: PetscCall(PetscDSGetUpperBound(ds, f, &bound, &ctx));
3943: PetscCall(PetscDSSetUpperBound(newds, f, bound, ctx));
3944: }
3945: PetscFunctionReturn(PETSC_SUCCESS);
3946: }
3948: /*@
3949: PetscDSCopy - Copy the contents of a `PetscDS` into another `PetscDS` on a new `DM`.
3951: Collective
3953: Input Parameters:
3954: + ds - the source `PetscDS`
3955: . minDegree - the minimum polynomial degree to consider when selecting discretizations, or `PETSC_DETERMINE`
3956: . maxDegree - the maximum polynomial degree to consider when selecting discretizations, or `PETSC_DETERMINE`
3957: - dmNew - the target `DM` used to resolve boundary condition labels for the copied boundaries
3959: Output Parameter:
3960: . dsNew - the destination `PetscDS`
3962: Level: developer
3964: Note:
3965: This copies constants, exact solutions, bounds, discretizations, equations, field contexts,
3966: cohesive flags, jet degrees, and boundary conditions.
3968: .seealso: `PetscDS`, `PetscDSCopyEquations()`, `PetscDSCopyConstants()`, `PetscDSCopyExactSolutions()`, `PetscDSCopyBounds()`, `PetscDSCopyBoundary()`
3969: @*/
3970: PetscErrorCode PetscDSCopy(PetscDS ds, PetscInt minDegree, PetscInt maxDegree, DM dmNew, PetscDS dsNew)
3971: {
3972: DSBoundary b;
3973: PetscInt cdim, Nf, f, d;
3974: PetscBool isCohesive;
3975: void *ctx;
3977: PetscFunctionBegin;
3978: PetscCall(PetscDSCopyConstants(ds, dsNew));
3979: PetscCall(PetscDSCopyExactSolutions(ds, dsNew));
3980: PetscCall(PetscDSCopyBounds(ds, dsNew));
3981: PetscCall(PetscDSSelectDiscretizations(ds, PETSC_DETERMINE, NULL, minDegree, maxDegree, dsNew));
3982: PetscCall(PetscDSCopyEquations(ds, dsNew));
3983: PetscCall(PetscDSGetNumFields(ds, &Nf));
3984: for (f = 0; f < Nf; ++f) {
3985: PetscCall(PetscDSGetContext(ds, f, &ctx));
3986: PetscCall(PetscDSSetContext(dsNew, f, ctx));
3987: PetscCall(PetscDSGetCohesive(ds, f, &isCohesive));
3988: PetscCall(PetscDSSetCohesive(dsNew, f, isCohesive));
3989: PetscCall(PetscDSGetJetDegree(ds, f, &d));
3990: PetscCall(PetscDSSetJetDegree(dsNew, f, d));
3991: }
3992: if (Nf) {
3993: PetscCall(PetscDSGetCoordinateDimension(ds, &cdim));
3994: PetscCall(PetscDSSetCoordinateDimension(dsNew, cdim));
3995: }
3996: PetscCall(PetscDSCopyBoundary(ds, PETSC_DETERMINE, NULL, dsNew));
3997: for (b = dsNew->boundary; b; b = b->next) {
3998: PetscCall(DMGetLabel(dmNew, b->lname, &b->label));
3999: /* Do not check if label exists here, since p4est calls this for the reference tree which does not have the labels */
4000: //PetscCheck(b->label,PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Label %s missing in new DM", name);
4001: }
4002: PetscFunctionReturn(PETSC_SUCCESS);
4003: }
4005: /*@
4006: PetscDSGetHeightSubspace - Get the `PetscDS` for the trace subspace at a given height in the mesh.
4008: Not Collective
4010: Input Parameters:
4011: + prob - the `PetscDS`
4012: - height - the height (0 for the ambient cell, 1 for faces, etc.)
4014: Output Parameter:
4015: . subprob - the `PetscDS` for the trace subspace; `prob` itself is returned when `height` is 0
4017: Level: developer
4019: Note:
4020: Only `PetscFE` discretizations are currently supported.
4022: .seealso: `PetscDS`, `PetscFE`, `PetscFEGetHeightSubspace()`, `PetscDSGetSpatialDimension()`
4023: @*/
4024: PetscErrorCode PetscDSGetHeightSubspace(PetscDS prob, PetscInt height, PetscDS *subprob)
4025: {
4026: PetscInt dim, Nf, f;
4028: PetscFunctionBegin;
4030: PetscAssertPointer(subprob, 3);
4031: if (height == 0) {
4032: *subprob = prob;
4033: PetscFunctionReturn(PETSC_SUCCESS);
4034: }
4035: PetscCall(PetscDSGetNumFields(prob, &Nf));
4036: PetscCall(PetscDSGetSpatialDimension(prob, &dim));
4037: PetscCheck(height <= dim, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_OUTOFRANGE, "DS can only handle height in [0, %" PetscInt_FMT "], not %" PetscInt_FMT, dim, height);
4038: if (!prob->subprobs) PetscCall(PetscCalloc1(dim, &prob->subprobs));
4039: if (!prob->subprobs[height - 1]) {
4040: PetscInt cdim;
4042: PetscCall(PetscDSCreate(PetscObjectComm((PetscObject)prob), &prob->subprobs[height - 1]));
4043: PetscCall(PetscDSGetCoordinateDimension(prob, &cdim));
4044: PetscCall(PetscDSSetCoordinateDimension(prob->subprobs[height - 1], cdim));
4045: for (f = 0; f < Nf; ++f) {
4046: PetscFE subfe;
4047: PetscObject obj;
4048: PetscClassId id;
4050: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
4051: PetscCall(PetscObjectGetClassId(obj, &id));
4052: PetscCheck(id == PETSCFE_CLASSID, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unsupported discretization type for field %" PetscInt_FMT, f);
4053: PetscCall(PetscFEGetHeightSubspace((PetscFE)obj, height, &subfe));
4054: PetscCall(PetscDSSetDiscretization(prob->subprobs[height - 1], f, (PetscObject)subfe));
4055: }
4056: }
4057: *subprob = prob->subprobs[height - 1];
4058: PetscFunctionReturn(PETSC_SUCCESS);
4059: }
4061: /*@
4062: PetscDSPermuteQuadPoint - Permute a quadrature point index according to a cell orientation.
4064: Not Collective
4066: Input Parameters:
4067: + ds - the `PetscDS`
4068: . ornt - the cell orientation, in `[-Na, Na)` where `Na` is half the number of arrangements for the cell type
4069: . field - the field number whose quadrature is used
4070: - q - the input quadrature point index in `[0, Nq)`
4072: Output Parameter:
4073: . qperm - the permuted quadrature point index
4075: Level: developer
4077: .seealso: `PetscDS`, `PetscQuadrature`, `PetscQuadratureComputePermutations()`, `DMPolytopeTypeGetNumArrangements()`
4078: @*/
4079: PetscErrorCode PetscDSPermuteQuadPoint(PetscDS ds, PetscInt ornt, PetscInt field, PetscInt q, PetscInt *qperm)
4080: {
4081: IS permIS;
4082: PetscQuadrature quad;
4083: DMPolytopeType ct;
4084: const PetscInt *perm;
4085: PetscInt Na, Nq;
4087: PetscFunctionBeginHot;
4088: PetscCall(PetscFEGetQuadrature((PetscFE)ds->disc[field], &quad));
4089: PetscCall(PetscQuadratureGetData(quad, NULL, NULL, &Nq, NULL, NULL));
4090: PetscCall(PetscQuadratureGetCellType(quad, &ct));
4091: PetscCheck(q >= 0 && q < Nq, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Quadrature point %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", q, Nq);
4092: Na = DMPolytopeTypeGetNumArrangements(ct) / 2;
4093: 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);
4094: if (!ds->quadPerm[(PetscInt)ct]) PetscCall(PetscQuadratureComputePermutations(quad, NULL, &ds->quadPerm[(PetscInt)ct]));
4095: permIS = ds->quadPerm[(PetscInt)ct][ornt + Na];
4096: PetscCall(ISGetIndices(permIS, &perm));
4097: *qperm = perm[q];
4098: PetscCall(ISRestoreIndices(permIS, &perm));
4099: PetscFunctionReturn(PETSC_SUCCESS);
4100: }
4102: PetscErrorCode PetscDSGetDiscType_Internal(PetscDS ds, PetscInt f, PetscDiscType *disctype)
4103: {
4104: PetscObject obj;
4105: PetscClassId id;
4106: PetscInt Nf;
4108: PetscFunctionBegin;
4110: PetscAssertPointer(disctype, 3);
4111: *disctype = PETSC_DISC_NONE;
4112: PetscCall(PetscDSGetNumFields(ds, &Nf));
4113: PetscCheck(f < Nf, PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_SIZ, "Field %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, Nf);
4114: PetscCall(PetscDSGetDiscretization(ds, f, &obj));
4115: if (obj) {
4116: PetscCall(PetscObjectGetClassId(obj, &id));
4117: if (id == PETSCFE_CLASSID) *disctype = PETSC_DISC_FE;
4118: else *disctype = PETSC_DISC_FV;
4119: }
4120: PetscFunctionReturn(PETSC_SUCCESS);
4121: }
4123: static PetscErrorCode PetscDSDestroy_Basic(PetscDS ds)
4124: {
4125: PetscFunctionBegin;
4126: PetscCall(PetscFree(ds->data));
4127: PetscFunctionReturn(PETSC_SUCCESS);
4128: }
4130: static PetscErrorCode PetscDSInitialize_Basic(PetscDS ds)
4131: {
4132: PetscFunctionBegin;
4133: ds->ops->setfromoptions = NULL;
4134: ds->ops->setup = NULL;
4135: ds->ops->view = NULL;
4136: ds->ops->destroy = PetscDSDestroy_Basic;
4137: PetscFunctionReturn(PETSC_SUCCESS);
4138: }
4140: /*MC
4141: PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions
4143: Level: intermediate
4145: .seealso: `PetscDSType`, `PetscDSCreate()`, `PetscDSSetType()`
4146: M*/
4148: PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS ds)
4149: {
4150: PetscDS_Basic *b;
4152: PetscFunctionBegin;
4154: PetscCall(PetscNew(&b));
4155: ds->data = b;
4157: PetscCall(PetscDSInitialize_Basic(ds));
4158: PetscFunctionReturn(PETSC_SUCCESS);
4159: }