Actual source code: dtds.c
1: #include <petsc/private/petscdsimpl.h>
3: PetscClassId PETSCDS_CLASSID = 0;
5: PetscFunctionList PetscDSList = NULL;
6: PetscBool PetscDSRegisterAllCalled = PETSC_FALSE;
8: /* A PetscDS (Discrete System) encodes a set of equations posed in a discrete space, which represents a set of
9: nonlinear continuum equations. The equations can have multiple fields, each field having a different
10: discretization. In addition, different pieces of the domain can have different field combinations and equations.
12: The DS provides the user a description of the approximation space on any given cell. It also gives pointwise
13: functions representing the equations.
15: Each field is associated with a label, marking the cells on which it is supported. Note that a field can be
16: supported on the closure of a cell not in the label due to overlap of the boundary of neighboring cells. The DM
17: then creates a DS for each set of cells with identical approximation spaces. When assembling, the user asks for
18: the space associated with a given cell. DMPlex uses the labels associated with each DS in the default integration loop.
19: */
21: /*@C
22: PetscDSRegister - Adds a new `PetscDS` implementation
24: Not Collective; No Fortran Support
26: Input Parameters:
27: + sname - The name of a new user-defined creation routine
28: - function - The creation routine itself
30: Example Usage:
31: .vb
32: PetscDSRegister("my_ds", MyPetscDSCreate);
33: .ve
35: Then, your PetscDS type can be chosen with the procedural interface via
36: .vb
37: PetscDSCreate(MPI_Comm, PetscDS *);
38: PetscDSSetType(PetscDS, "my_ds");
39: .ve
40: or at runtime via the option
41: .vb
42: -petscds_type my_ds
43: .ve
45: Level: advanced
47: Note:
48: `PetscDSRegister()` may be called multiple times to add several user-defined `PetscDSs`
50: .seealso: `PetscDSType`, `PetscDS`, `PetscDSRegisterAll()`, `PetscDSRegisterDestroy()`
51: @*/
52: PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS))
53: {
54: PetscFunctionBegin;
55: PetscCall(PetscFunctionListAdd(&PetscDSList, sname, function));
56: PetscFunctionReturn(PETSC_SUCCESS);
57: }
59: /*@
60: PetscDSSetType - Builds a particular `PetscDS`
62: Collective; No Fortran Support
64: Input Parameters:
65: + prob - The `PetscDS` object
66: - name - The `PetscDSType`
68: Options Database Key:
69: . -petscds_type <type> - Sets the PetscDS type; use -help for a list of available types
71: Level: intermediate
73: .seealso: `PetscDSType`, `PetscDS`, `PetscDSGetType()`, `PetscDSCreate()`
74: @*/
75: PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name)
76: {
77: PetscErrorCode (*r)(PetscDS);
78: PetscBool match;
80: PetscFunctionBegin;
82: PetscCall(PetscObjectTypeCompare((PetscObject)prob, name, &match));
83: if (match) PetscFunctionReturn(PETSC_SUCCESS);
85: PetscCall(PetscDSRegisterAll());
86: PetscCall(PetscFunctionListFind(PetscDSList, name, &r));
87: PetscCheck(r, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscDS type: %s", name);
89: PetscTryTypeMethod(prob, destroy);
90: prob->ops->destroy = NULL;
92: PetscCall((*r)(prob));
93: PetscCall(PetscObjectChangeTypeName((PetscObject)prob, name));
94: PetscFunctionReturn(PETSC_SUCCESS);
95: }
97: /*@
98: PetscDSGetType - Gets the `PetscDSType` name (as a string) from the `PetscDS`
100: Not Collective; No Fortran Support
102: Input Parameter:
103: . prob - The `PetscDS`
105: Output Parameter:
106: . name - The `PetscDSType` name
108: Level: intermediate
110: .seealso: `PetscDSType`, `PetscDS`, `PetscDSSetType()`, `PetscDSCreate()`
111: @*/
112: PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name)
113: {
114: PetscFunctionBegin;
116: PetscAssertPointer(name, 2);
117: PetscCall(PetscDSRegisterAll());
118: *name = ((PetscObject)prob)->type_name;
119: PetscFunctionReturn(PETSC_SUCCESS);
120: }
122: static PetscErrorCode PetscDSView_Ascii(PetscDS ds, PetscViewer viewer)
123: {
124: PetscViewerFormat format;
125: const PetscScalar *constants;
126: PetscInt Nf, numConstants, f;
128: PetscFunctionBegin;
129: PetscCall(PetscDSGetNumFields(ds, &Nf));
130: PetscCall(PetscViewerGetFormat(viewer, &format));
131: PetscCall(PetscViewerASCIIPrintf(viewer, "Discrete System with %" PetscInt_FMT " fields\n", Nf));
132: PetscCall(PetscViewerASCIIPushTab(viewer));
133: PetscCall(PetscViewerASCIIPrintf(viewer, " cell total dim %" PetscInt_FMT " total comp %" PetscInt_FMT "\n", ds->totDim, ds->totComp));
134: if (ds->isCohesive) PetscCall(PetscViewerASCIIPrintf(viewer, " cohesive cell\n"));
135: for (f = 0; f < Nf; ++f) {
136: DSBoundary b;
137: PetscObject obj;
138: PetscClassId id;
139: PetscQuadrature q;
140: const char *name;
141: PetscInt Nc, Nq, Nqc;
143: PetscCall(PetscDSGetDiscretization(ds, f, &obj));
144: PetscCall(PetscObjectGetClassId(obj, &id));
145: PetscCall(PetscObjectGetName(obj, &name));
146: PetscCall(PetscViewerASCIIPrintf(viewer, "Field %s", name ? name : "<unknown>"));
147: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
148: if (id == PETSCFE_CLASSID) {
149: PetscCall(PetscFEGetNumComponents((PetscFE)obj, &Nc));
150: PetscCall(PetscFEGetQuadrature((PetscFE)obj, &q));
151: PetscCall(PetscViewerASCIIPrintf(viewer, " FEM"));
152: } else if (id == PETSCFV_CLASSID) {
153: PetscCall(PetscFVGetNumComponents((PetscFV)obj, &Nc));
154: PetscCall(PetscFVGetQuadrature((PetscFV)obj, &q));
155: PetscCall(PetscViewerASCIIPrintf(viewer, " FVM"));
156: } else SETERRQ(PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f);
157: if (Nc > 1) PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " components", Nc));
158: else PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " component ", Nc));
159: if (ds->implicit[f]) PetscCall(PetscViewerASCIIPrintf(viewer, " (implicit)"));
160: else PetscCall(PetscViewerASCIIPrintf(viewer, " (explicit)"));
161: if (q) {
162: PetscCall(PetscQuadratureGetData(q, NULL, &Nqc, &Nq, NULL, NULL));
163: PetscCall(PetscViewerASCIIPrintf(viewer, " (Nq %" PetscInt_FMT " Nqc %" PetscInt_FMT ")", Nq, Nqc));
164: }
165: PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT "-jet", ds->jetDegree[f]));
166: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
167: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
168: PetscCall(PetscViewerASCIIPushTab(viewer));
169: if (id == PETSCFE_CLASSID) PetscCall(PetscFEView((PetscFE)obj, viewer));
170: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVView((PetscFV)obj, viewer));
171: PetscCall(PetscViewerASCIIPopTab(viewer));
173: for (b = ds->boundary; b; b = b->next) {
174: char *name;
175: PetscInt c, i;
177: if (b->field != f) continue;
178: PetscCall(PetscViewerASCIIPushTab(viewer));
179: PetscCall(PetscViewerASCIIPrintf(viewer, "Boundary %s (%s) %s\n", b->name, b->lname, DMBoundaryConditionTypes[b->type]));
180: if (!b->Nc) {
181: PetscCall(PetscViewerASCIIPrintf(viewer, " all components\n"));
182: } else {
183: PetscCall(PetscViewerASCIIPrintf(viewer, " components: "));
184: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
185: for (c = 0; c < b->Nc; ++c) {
186: if (c > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", "));
187: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->comps[c]));
188: }
189: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
190: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
191: }
192: PetscCall(PetscViewerASCIIPrintf(viewer, " values: "));
193: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
194: for (i = 0; i < b->Nv; ++i) {
195: if (i > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", "));
196: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->values[i]));
197: }
198: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
199: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
200: #if defined(__clang__)
201: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat-pedantic")
202: #elif defined(__GNUC__) || defined(__GNUG__)
203: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat")
204: #endif
205: if (b->func) {
206: PetscCall(PetscDLAddr(b->func, &name));
207: if (name) PetscCall(PetscViewerASCIIPrintf(viewer, " func: %s\n", name));
208: else PetscCall(PetscViewerASCIIPrintf(viewer, " func: %p\n", b->func));
209: PetscCall(PetscFree(name));
210: }
211: if (b->func_t) {
212: PetscCall(PetscDLAddr(b->func_t, &name));
213: if (name) PetscCall(PetscViewerASCIIPrintf(viewer, " func_t: %s\n", name));
214: else PetscCall(PetscViewerASCIIPrintf(viewer, " func_t: %p\n", b->func_t));
215: PetscCall(PetscFree(name));
216: }
217: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_END()
218: PetscCall(PetscWeakFormView(b->wf, viewer));
219: PetscCall(PetscViewerASCIIPopTab(viewer));
220: }
221: }
222: PetscCall(PetscDSGetConstants(ds, &numConstants, &constants));
223: if (numConstants) {
224: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " constants\n", numConstants));
225: PetscCall(PetscViewerASCIIPushTab(viewer));
226: for (f = 0; f < numConstants; ++f) PetscCall(PetscViewerASCIIPrintf(viewer, "%g\n", (double)PetscRealPart(constants[f])));
227: PetscCall(PetscViewerASCIIPopTab(viewer));
228: }
229: PetscCall(PetscWeakFormView(ds->wf, viewer));
230: PetscCall(PetscViewerASCIIPopTab(viewer));
231: PetscFunctionReturn(PETSC_SUCCESS);
232: }
234: /*@
235: PetscDSViewFromOptions - View a `PetscDS` based on values in the options database
237: Collective
239: Input Parameters:
240: + A - the `PetscDS` object
241: . obj - Optional object that provides the options prefix used in the search
242: - name - command line option
244: Level: intermediate
246: .seealso: `PetscDSType`, `PetscDS`, `PetscDSView()`, `PetscObjectViewFromOptions()`, `PetscDSCreate()`
247: @*/
248: PetscErrorCode PetscDSViewFromOptions(PetscDS A, PetscObject obj, const char name[])
249: {
250: PetscFunctionBegin;
252: PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
253: PetscFunctionReturn(PETSC_SUCCESS);
254: }
256: /*@
257: PetscDSView - Views a `PetscDS`
259: Collective
261: Input Parameters:
262: + prob - the `PetscDS` object to view
263: - v - the viewer
265: Level: developer
267: .seealso: `PetscDSType`, `PetscDS`, `PetscViewer`, `PetscDSDestroy()`, `PetscDSViewFromOptions()`
268: @*/
269: PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v)
270: {
271: PetscBool iascii;
273: PetscFunctionBegin;
275: if (!v) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)prob), &v));
277: PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &iascii));
278: if (iascii) PetscCall(PetscDSView_Ascii(prob, v));
279: PetscTryTypeMethod(prob, view, v);
280: PetscFunctionReturn(PETSC_SUCCESS);
281: }
283: /*@
284: PetscDSSetFromOptions - sets parameters in a `PetscDS` from the options database
286: Collective
288: Input Parameter:
289: . prob - the `PetscDS` object to set options for
291: Options Database Keys:
292: + -petscds_type <type> - Set the `PetscDS` type
293: . -petscds_view <view opt> - View the `PetscDS`
294: . -petscds_jac_pre - Turn formation of a separate Jacobian preconditioner on or off
295: . -bc_<name> <ids> - Specify a list of label ids for a boundary condition
296: - -bc_<name>_comp <comps> - Specify a list of field components to constrain for a boundary condition
298: Level: intermediate
300: .seealso: `PetscDS`, `PetscDSView()`
301: @*/
302: PetscErrorCode PetscDSSetFromOptions(PetscDS prob)
303: {
304: DSBoundary b;
305: const char *defaultType;
306: char name[256];
307: PetscBool flg;
309: PetscFunctionBegin;
311: if (!((PetscObject)prob)->type_name) {
312: defaultType = PETSCDSBASIC;
313: } else {
314: defaultType = ((PetscObject)prob)->type_name;
315: }
316: PetscCall(PetscDSRegisterAll());
318: PetscObjectOptionsBegin((PetscObject)prob);
319: for (b = prob->boundary; b; b = b->next) {
320: char optname[1024];
321: PetscInt ids[1024], len = 1024;
322: PetscBool flg;
324: PetscCall(PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name));
325: PetscCall(PetscMemzero(ids, sizeof(ids)));
326: PetscCall(PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg));
327: if (flg) {
328: b->Nv = len;
329: PetscCall(PetscFree(b->values));
330: PetscCall(PetscMalloc1(len, &b->values));
331: PetscCall(PetscArraycpy(b->values, ids, len));
332: PetscCall(PetscWeakFormRewriteKeys(b->wf, b->label, len, b->values));
333: }
334: len = 1024;
335: PetscCall(PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name));
336: PetscCall(PetscMemzero(ids, sizeof(ids)));
337: PetscCall(PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg));
338: if (flg) {
339: b->Nc = len;
340: PetscCall(PetscFree(b->comps));
341: PetscCall(PetscMalloc1(len, &b->comps));
342: PetscCall(PetscArraycpy(b->comps, ids, len));
343: }
344: }
345: PetscCall(PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg));
346: if (flg) {
347: PetscCall(PetscDSSetType(prob, name));
348: } else if (!((PetscObject)prob)->type_name) {
349: PetscCall(PetscDSSetType(prob, defaultType));
350: }
351: PetscCall(PetscOptionsBool("-petscds_jac_pre", "Discrete System", "PetscDSUseJacobianPreconditioner", prob->useJacPre, &prob->useJacPre, &flg));
352: PetscCall(PetscOptionsBool("-petscds_force_quad", "Discrete System", "PetscDSSetForceQuad", prob->forceQuad, &prob->forceQuad, &flg));
353: PetscCall(PetscOptionsInt("-petscds_print_integrate", "Discrete System", "", prob->printIntegrate, &prob->printIntegrate, NULL));
354: PetscTryTypeMethod(prob, setfromoptions);
355: /* process any options handlers added with PetscObjectAddOptionsHandler() */
356: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)prob, PetscOptionsObject));
357: PetscOptionsEnd();
358: if (prob->Nf) PetscCall(PetscDSViewFromOptions(prob, NULL, "-petscds_view"));
359: PetscFunctionReturn(PETSC_SUCCESS);
360: }
362: /*@
363: PetscDSSetUp - Construct data structures for the `PetscDS`
365: Collective
367: Input Parameter:
368: . prob - the `PetscDS` object to setup
370: Level: developer
372: .seealso: `PetscDS`, `PetscDSView()`, `PetscDSDestroy()`
373: @*/
374: PetscErrorCode PetscDSSetUp(PetscDS prob)
375: {
376: const PetscInt Nf = prob->Nf;
377: PetscBool hasH = PETSC_FALSE;
378: PetscInt maxOrder[4] = {-2, -2, -2, -2};
379: PetscInt dim, dimEmbed, NbMax = 0, NcMax = 0, NqMax = 0, NsMax = 1, f;
381: PetscFunctionBegin;
383: if (prob->setup) PetscFunctionReturn(PETSC_SUCCESS);
384: /* Calculate sizes */
385: PetscCall(PetscDSGetSpatialDimension(prob, &dim));
386: PetscCall(PetscDSGetCoordinateDimension(prob, &dimEmbed));
387: prob->totDim = prob->totComp = 0;
388: PetscCall(PetscMalloc2(Nf, &prob->Nc, Nf, &prob->Nb));
389: PetscCall(PetscCalloc2(Nf + 1, &prob->off, Nf + 1, &prob->offDer));
390: 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]));
391: PetscCall(PetscMalloc2(Nf, &prob->T, Nf, &prob->Tf));
392: if (prob->forceQuad) {
393: // Note: This assumes we have one kind of cell at each dimension.
394: // We can fix this by having quadrature hold the celltype
395: PetscQuadrature maxQuad[4] = {NULL, NULL, NULL, NULL};
397: for (f = 0; f < Nf; ++f) {
398: PetscObject obj;
399: PetscClassId id;
400: PetscQuadrature q = NULL, fq = NULL;
401: PetscInt dim = -1, order = -1, forder = -1;
403: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
404: if (!obj) continue;
405: PetscCall(PetscObjectGetClassId(obj, &id));
406: if (id == PETSCFE_CLASSID) {
407: PetscFE fe = (PetscFE)obj;
409: PetscCall(PetscFEGetQuadrature(fe, &q));
410: PetscCall(PetscFEGetFaceQuadrature(fe, &fq));
411: } else if (id == PETSCFV_CLASSID) {
412: PetscFV fv = (PetscFV)obj;
414: PetscCall(PetscFVGetQuadrature(fv, &q));
415: }
416: if (q) {
417: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
418: PetscCall(PetscQuadratureGetOrder(q, &order));
419: if (order > maxOrder[dim]) {
420: maxOrder[dim] = order;
421: maxQuad[dim] = q;
422: }
423: }
424: if (fq) {
425: PetscCall(PetscQuadratureGetData(fq, &dim, NULL, NULL, NULL, NULL));
426: PetscCall(PetscQuadratureGetOrder(fq, &forder));
427: if (forder > maxOrder[dim]) {
428: maxOrder[dim] = forder;
429: maxQuad[dim] = fq;
430: }
431: }
432: }
433: for (f = 0; f < Nf; ++f) {
434: PetscObject obj;
435: PetscClassId id;
436: PetscQuadrature q;
437: PetscInt dim;
439: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
440: if (!obj) continue;
441: PetscCall(PetscObjectGetClassId(obj, &id));
442: if (id == PETSCFE_CLASSID) {
443: PetscFE fe = (PetscFE)obj;
445: PetscCall(PetscFEGetQuadrature(fe, &q));
446: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
447: PetscCall(PetscFESetQuadrature(fe, maxQuad[dim]));
448: PetscCall(PetscFESetFaceQuadrature(fe, dim ? maxQuad[dim - 1] : NULL));
449: } else if (id == PETSCFV_CLASSID) {
450: PetscFV fv = (PetscFV)obj;
452: PetscCall(PetscFVGetQuadrature(fv, &q));
453: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
454: PetscCall(PetscFVSetQuadrature(fv, maxQuad[dim]));
455: }
456: }
457: }
458: for (f = 0; f < Nf; ++f) {
459: PetscObject obj;
460: PetscClassId id;
461: PetscQuadrature q = NULL;
462: PetscInt Nq = 0, Nb, Nc;
464: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
465: if (prob->jetDegree[f] > 1) hasH = PETSC_TRUE;
466: if (!obj) {
467: /* Empty mesh */
468: Nb = Nc = 0;
469: prob->T[f] = prob->Tf[f] = NULL;
470: } else {
471: PetscCall(PetscObjectGetClassId(obj, &id));
472: if (id == PETSCFE_CLASSID) {
473: PetscFE fe = (PetscFE)obj;
475: PetscCall(PetscFEGetQuadrature(fe, &q));
476: {
477: PetscQuadrature fq;
478: PetscInt dim, order;
480: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
481: PetscCall(PetscQuadratureGetOrder(q, &order));
482: if (maxOrder[dim] < 0) maxOrder[dim] = order;
483: 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]);
484: PetscCall(PetscFEGetFaceQuadrature(fe, &fq));
485: if (fq) {
486: PetscCall(PetscQuadratureGetData(fq, &dim, NULL, NULL, NULL, NULL));
487: PetscCall(PetscQuadratureGetOrder(fq, &order));
488: if (maxOrder[dim] < 0) maxOrder[dim] = order;
489: 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]);
490: }
491: }
492: PetscCall(PetscFEGetDimension(fe, &Nb));
493: PetscCall(PetscFEGetNumComponents(fe, &Nc));
494: PetscCall(PetscFEGetCellTabulation(fe, prob->jetDegree[f], &prob->T[f]));
495: PetscCall(PetscFEGetFaceTabulation(fe, prob->jetDegree[f], &prob->Tf[f]));
496: } else if (id == PETSCFV_CLASSID) {
497: PetscFV fv = (PetscFV)obj;
499: PetscCall(PetscFVGetQuadrature(fv, &q));
500: PetscCall(PetscFVGetNumComponents(fv, &Nc));
501: Nb = Nc;
502: PetscCall(PetscFVGetCellTabulation(fv, &prob->T[f]));
503: /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */
504: } else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f);
505: }
506: prob->Nc[f] = Nc;
507: prob->Nb[f] = Nb;
508: prob->off[f + 1] = Nc + prob->off[f];
509: prob->offDer[f + 1] = Nc * dim + prob->offDer[f];
510: prob->offCohesive[0][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[0][f];
511: prob->offDerCohesive[0][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[0][f];
512: prob->offCohesive[1][f] = (prob->cohesive[f] ? 0 : Nc) + prob->offCohesive[0][f];
513: prob->offDerCohesive[1][f] = (prob->cohesive[f] ? 0 : Nc) * dimEmbed + prob->offDerCohesive[0][f];
514: prob->offCohesive[2][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[2][f];
515: prob->offDerCohesive[2][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[2][f];
516: if (q) PetscCall(PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL));
517: NqMax = PetscMax(NqMax, Nq);
518: NbMax = PetscMax(NbMax, Nb);
519: NcMax = PetscMax(NcMax, Nc);
520: prob->totDim += Nb;
521: prob->totComp += Nc;
522: /* There are two faces for all fields on a cohesive cell, except for cohesive fields */
523: if (prob->isCohesive && !prob->cohesive[f]) prob->totDim += Nb;
524: }
525: prob->offCohesive[1][Nf] = prob->offCohesive[0][Nf];
526: prob->offDerCohesive[1][Nf] = prob->offDerCohesive[0][Nf];
527: /* Allocate works space */
528: NsMax = 2; /* A non-cohesive discretizations can be used on a cohesive cell, so we need this extra workspace for all DS */
529: 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));
530: PetscCall(PetscMalloc5(dimEmbed, &prob->x, NbMax * NcMax, &prob->basisReal, NbMax * NcMax * dimEmbed, &prob->basisDerReal, NbMax * NcMax, &prob->testReal, NbMax * NcMax * dimEmbed, &prob->testDerReal));
531: 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,
532: &prob->g2, NsMax * NsMax * NqMax * NcMax * NcMax * dimEmbed * dimEmbed, &prob->g3));
533: PetscTryTypeMethod(prob, setup);
534: prob->setup = PETSC_TRUE;
535: PetscFunctionReturn(PETSC_SUCCESS);
536: }
538: static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob)
539: {
540: PetscFunctionBegin;
541: PetscCall(PetscFree2(prob->Nc, prob->Nb));
542: PetscCall(PetscFree2(prob->off, prob->offDer));
543: PetscCall(PetscFree6(prob->offCohesive[0], prob->offCohesive[1], prob->offCohesive[2], prob->offDerCohesive[0], prob->offDerCohesive[1], prob->offDerCohesive[2]));
544: PetscCall(PetscFree2(prob->T, prob->Tf));
545: PetscCall(PetscFree3(prob->u, prob->u_t, prob->u_x));
546: PetscCall(PetscFree5(prob->x, prob->basisReal, prob->basisDerReal, prob->testReal, prob->testDerReal));
547: PetscCall(PetscFree6(prob->f0, prob->f1, prob->g0, prob->g1, prob->g2, prob->g3));
548: PetscFunctionReturn(PETSC_SUCCESS);
549: }
551: static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew)
552: {
553: PetscObject *tmpd;
554: PetscBool *tmpi;
555: PetscInt *tmpk;
556: PetscBool *tmpc;
557: PetscPointFunc *tmpup;
558: PetscSimplePointFn **tmpexactSol, **tmpexactSol_t;
559: void **tmpexactCtx, **tmpexactCtx_t;
560: void **tmpctx;
561: PetscInt Nf = prob->Nf, f;
563: PetscFunctionBegin;
564: if (Nf >= NfNew) PetscFunctionReturn(PETSC_SUCCESS);
565: prob->setup = PETSC_FALSE;
566: PetscCall(PetscDSDestroyStructs_Static(prob));
567: PetscCall(PetscMalloc4(NfNew, &tmpd, NfNew, &tmpi, NfNew, &tmpc, NfNew, &tmpk));
568: for (f = 0; f < Nf; ++f) {
569: tmpd[f] = prob->disc[f];
570: tmpi[f] = prob->implicit[f];
571: tmpc[f] = prob->cohesive[f];
572: tmpk[f] = prob->jetDegree[f];
573: }
574: for (f = Nf; f < NfNew; ++f) {
575: tmpd[f] = NULL;
576: tmpi[f] = PETSC_TRUE, tmpc[f] = PETSC_FALSE;
577: tmpk[f] = 1;
578: }
579: PetscCall(PetscFree4(prob->disc, prob->implicit, prob->cohesive, prob->jetDegree));
580: PetscCall(PetscWeakFormSetNumFields(prob->wf, NfNew));
581: prob->Nf = NfNew;
582: prob->disc = tmpd;
583: prob->implicit = tmpi;
584: prob->cohesive = tmpc;
585: prob->jetDegree = tmpk;
586: PetscCall(PetscCalloc2(NfNew, &tmpup, NfNew, &tmpctx));
587: for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f];
588: for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f];
589: for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL;
590: for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL;
591: PetscCall(PetscFree2(prob->update, prob->ctx));
592: prob->update = tmpup;
593: prob->ctx = tmpctx;
594: PetscCall(PetscCalloc4(NfNew, &tmpexactSol, NfNew, &tmpexactCtx, NfNew, &tmpexactSol_t, NfNew, &tmpexactCtx_t));
595: for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f];
596: for (f = 0; f < Nf; ++f) tmpexactCtx[f] = prob->exactCtx[f];
597: for (f = 0; f < Nf; ++f) tmpexactSol_t[f] = prob->exactSol_t[f];
598: for (f = 0; f < Nf; ++f) tmpexactCtx_t[f] = prob->exactCtx_t[f];
599: for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL;
600: for (f = Nf; f < NfNew; ++f) tmpexactCtx[f] = NULL;
601: for (f = Nf; f < NfNew; ++f) tmpexactSol_t[f] = NULL;
602: for (f = Nf; f < NfNew; ++f) tmpexactCtx_t[f] = NULL;
603: PetscCall(PetscFree4(prob->exactSol, prob->exactCtx, prob->exactSol_t, prob->exactCtx_t));
604: prob->exactSol = tmpexactSol;
605: prob->exactCtx = tmpexactCtx;
606: prob->exactSol_t = tmpexactSol_t;
607: prob->exactCtx_t = tmpexactCtx_t;
608: PetscFunctionReturn(PETSC_SUCCESS);
609: }
611: /*@
612: PetscDSDestroy - Destroys a `PetscDS` object
614: Collective
616: Input Parameter:
617: . ds - the `PetscDS` object to destroy
619: Level: developer
621: .seealso: `PetscDSView()`
622: @*/
623: PetscErrorCode PetscDSDestroy(PetscDS *ds)
624: {
625: PetscInt f;
627: PetscFunctionBegin;
628: if (!*ds) PetscFunctionReturn(PETSC_SUCCESS);
631: if (--((PetscObject)*ds)->refct > 0) {
632: *ds = NULL;
633: PetscFunctionReturn(PETSC_SUCCESS);
634: }
635: ((PetscObject)*ds)->refct = 0;
636: if ((*ds)->subprobs) {
637: PetscInt dim, d;
639: PetscCall(PetscDSGetSpatialDimension(*ds, &dim));
640: for (d = 0; d < dim; ++d) PetscCall(PetscDSDestroy(&(*ds)->subprobs[d]));
641: }
642: PetscCall(PetscFree((*ds)->subprobs));
643: PetscCall(PetscDSDestroyStructs_Static(*ds));
644: for (f = 0; f < (*ds)->Nf; ++f) PetscCall(PetscObjectDereference((*ds)->disc[f]));
645: PetscCall(PetscFree4((*ds)->disc, (*ds)->implicit, (*ds)->cohesive, (*ds)->jetDegree));
646: PetscCall(PetscWeakFormDestroy(&(*ds)->wf));
647: PetscCall(PetscFree2((*ds)->update, (*ds)->ctx));
648: PetscCall(PetscFree4((*ds)->exactSol, (*ds)->exactCtx, (*ds)->exactSol_t, (*ds)->exactCtx_t));
649: PetscTryTypeMethod(*ds, destroy);
650: PetscCall(PetscDSDestroyBoundary(*ds));
651: PetscCall(PetscFree((*ds)->constants));
652: for (PetscInt c = 0; c < DM_NUM_POLYTOPES; ++c) {
653: const PetscInt Na = DMPolytopeTypeGetNumArrangements((DMPolytopeType)c);
654: if ((*ds)->quadPerm[c])
655: for (PetscInt o = 0; o < Na; ++o) PetscCall(ISDestroy(&(*ds)->quadPerm[c][o]));
656: PetscCall(PetscFree((*ds)->quadPerm[c]));
657: (*ds)->quadPerm[c] = NULL;
658: }
659: PetscCall(PetscHeaderDestroy(ds));
660: PetscFunctionReturn(PETSC_SUCCESS);
661: }
663: /*@
664: PetscDSCreate - Creates an empty `PetscDS` object. The type can then be set with `PetscDSSetType()`.
666: Collective
668: Input Parameter:
669: . comm - The communicator for the `PetscDS` object
671: Output Parameter:
672: . ds - The `PetscDS` object
674: Level: beginner
676: .seealso: `PetscDS`, `PetscDSSetType()`, `PETSCDSBASIC`, `PetscDSType`
677: @*/
678: PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *ds)
679: {
680: PetscDS p;
682: PetscFunctionBegin;
683: PetscAssertPointer(ds, 2);
684: PetscCall(PetscDSInitializePackage());
686: PetscCall(PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView));
687: p->Nf = 0;
688: p->setup = PETSC_FALSE;
689: p->numConstants = 0;
690: p->numFuncConstants = 3; // Row and col fields, cell size
691: p->dimEmbed = -1;
692: p->useJacPre = PETSC_TRUE;
693: p->forceQuad = PETSC_TRUE;
694: PetscCall(PetscMalloc1(p->numConstants + p->numFuncConstants, &p->constants));
695: PetscCall(PetscWeakFormCreate(comm, &p->wf));
696: PetscCall(PetscArrayzero(p->quadPerm, DM_NUM_POLYTOPES));
697: *ds = p;
698: PetscFunctionReturn(PETSC_SUCCESS);
699: }
701: /*@
702: PetscDSGetNumFields - Returns the number of fields in the `PetscDS`
704: Not Collective
706: Input Parameter:
707: . prob - The `PetscDS` object
709: Output Parameter:
710: . Nf - The number of fields
712: Level: beginner
714: .seealso: `PetscDS`, `PetscDSGetSpatialDimension()`, `PetscDSCreate()`
715: @*/
716: PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf)
717: {
718: PetscFunctionBegin;
720: PetscAssertPointer(Nf, 2);
721: *Nf = prob->Nf;
722: PetscFunctionReturn(PETSC_SUCCESS);
723: }
725: /*@
726: PetscDSGetSpatialDimension - Returns the spatial dimension of the `PetscDS`, meaning the topological dimension of the discretizations
728: Not Collective
730: Input Parameter:
731: . prob - The `PetscDS` object
733: Output Parameter:
734: . dim - The spatial dimension
736: Level: beginner
738: .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
739: @*/
740: PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim)
741: {
742: PetscFunctionBegin;
744: PetscAssertPointer(dim, 2);
745: *dim = 0;
746: if (prob->Nf) {
747: PetscObject obj;
748: PetscClassId id;
750: PetscCall(PetscDSGetDiscretization(prob, 0, &obj));
751: if (obj) {
752: PetscCall(PetscObjectGetClassId(obj, &id));
753: if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetSpatialDimension((PetscFE)obj, dim));
754: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetSpatialDimension((PetscFV)obj, dim));
755: else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
756: }
757: }
758: PetscFunctionReturn(PETSC_SUCCESS);
759: }
761: /*@
762: PetscDSGetCoordinateDimension - Returns the coordinate dimension of the `PetscDS`, meaning the dimension of the space into which the discretiaztions are embedded
764: Not Collective
766: Input Parameter:
767: . prob - The `PetscDS` object
769: Output Parameter:
770: . dimEmbed - The coordinate dimension
772: Level: beginner
774: .seealso: `PetscDS`, `PetscDSSetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
775: @*/
776: PetscErrorCode PetscDSGetCoordinateDimension(PetscDS prob, PetscInt *dimEmbed)
777: {
778: PetscFunctionBegin;
780: PetscAssertPointer(dimEmbed, 2);
781: PetscCheck(prob->dimEmbed >= 0, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONGSTATE, "No coordinate dimension set for this DS");
782: *dimEmbed = prob->dimEmbed;
783: PetscFunctionReturn(PETSC_SUCCESS);
784: }
786: /*@
787: PetscDSSetCoordinateDimension - Set the coordinate dimension of the `PetscDS`, meaning the dimension of the space into which the discretiaztions are embedded
789: Logically Collective
791: Input Parameters:
792: + prob - The `PetscDS` object
793: - dimEmbed - The coordinate dimension
795: Level: beginner
797: .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
798: @*/
799: PetscErrorCode PetscDSSetCoordinateDimension(PetscDS prob, PetscInt dimEmbed)
800: {
801: PetscFunctionBegin;
803: PetscCheck(dimEmbed >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Coordinate dimension must be non-negative, not %" PetscInt_FMT, dimEmbed);
804: prob->dimEmbed = dimEmbed;
805: PetscFunctionReturn(PETSC_SUCCESS);
806: }
808: /*@
809: PetscDSGetForceQuad - Returns the flag to force matching quadratures among the field discretizations
811: Not collective
813: Input Parameter:
814: . ds - The `PetscDS` object
816: Output Parameter:
817: . forceQuad - The flag
819: Level: intermediate
821: .seealso: `PetscDS`, `PetscDSSetForceQuad()`, `PetscDSGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
822: @*/
823: PetscErrorCode PetscDSGetForceQuad(PetscDS ds, PetscBool *forceQuad)
824: {
825: PetscFunctionBegin;
827: PetscAssertPointer(forceQuad, 2);
828: *forceQuad = ds->forceQuad;
829: PetscFunctionReturn(PETSC_SUCCESS);
830: }
832: /*@
833: PetscDSSetForceQuad - Set the flag to force matching quadratures among the field discretizations
835: Logically collective on ds
837: Input Parameters:
838: + ds - The `PetscDS` object
839: - forceQuad - The flag
841: Level: intermediate
843: .seealso: `PetscDS`, `PetscDSGetForceQuad()`, `PetscDSGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
844: @*/
845: PetscErrorCode PetscDSSetForceQuad(PetscDS ds, PetscBool forceQuad)
846: {
847: PetscFunctionBegin;
849: ds->forceQuad = forceQuad;
850: PetscFunctionReturn(PETSC_SUCCESS);
851: }
853: /*@
854: PetscDSIsCohesive - Returns the flag indicating that this `PetscDS` is for a cohesive cell
856: Not Collective
858: Input Parameter:
859: . ds - The `PetscDS` object
861: Output Parameter:
862: . isCohesive - The flag
864: Level: developer
866: .seealso: `PetscDS`, `PetscDSGetNumCohesive()`, `PetscDSGetCohesive()`, `PetscDSSetCohesive()`, `PetscDSCreate()`
867: @*/
868: PetscErrorCode PetscDSIsCohesive(PetscDS ds, PetscBool *isCohesive)
869: {
870: PetscFunctionBegin;
872: PetscAssertPointer(isCohesive, 2);
873: *isCohesive = ds->isCohesive;
874: PetscFunctionReturn(PETSC_SUCCESS);
875: }
877: /*@
878: PetscDSGetNumCohesive - Returns the number of cohesive fields, meaning those defined on the interior of a cohesive cell
880: Not Collective
882: Input Parameter:
883: . ds - The `PetscDS` object
885: Output Parameter:
886: . numCohesive - The number of cohesive fields
888: Level: developer
890: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSCreate()`
891: @*/
892: PetscErrorCode PetscDSGetNumCohesive(PetscDS ds, PetscInt *numCohesive)
893: {
894: PetscInt f;
896: PetscFunctionBegin;
898: PetscAssertPointer(numCohesive, 2);
899: *numCohesive = 0;
900: for (f = 0; f < ds->Nf; ++f) *numCohesive += ds->cohesive[f] ? 1 : 0;
901: PetscFunctionReturn(PETSC_SUCCESS);
902: }
904: /*@
905: PetscDSGetCohesive - Returns the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell
907: Not Collective
909: Input Parameters:
910: + ds - The `PetscDS` object
911: - f - The field index
913: Output Parameter:
914: . isCohesive - The flag
916: Level: developer
918: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
919: @*/
920: PetscErrorCode PetscDSGetCohesive(PetscDS ds, PetscInt f, PetscBool *isCohesive)
921: {
922: PetscFunctionBegin;
924: PetscAssertPointer(isCohesive, 3);
925: 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);
926: *isCohesive = ds->cohesive[f];
927: PetscFunctionReturn(PETSC_SUCCESS);
928: }
930: /*@
931: PetscDSSetCohesive - Set the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell
933: Not Collective
935: Input Parameters:
936: + ds - The `PetscDS` object
937: . f - The field index
938: - isCohesive - The flag for a cohesive field
940: Level: developer
942: .seealso: `PetscDS`, `PetscDSGetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
943: @*/
944: PetscErrorCode PetscDSSetCohesive(PetscDS ds, PetscInt f, PetscBool isCohesive)
945: {
946: PetscInt i;
948: PetscFunctionBegin;
950: 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);
951: ds->cohesive[f] = isCohesive;
952: ds->isCohesive = PETSC_FALSE;
953: for (i = 0; i < ds->Nf; ++i) ds->isCohesive = ds->isCohesive || ds->cohesive[f] ? PETSC_TRUE : PETSC_FALSE;
954: PetscFunctionReturn(PETSC_SUCCESS);
955: }
957: /*@
958: PetscDSGetTotalDimension - Returns the total size of the approximation space for this system
960: Not Collective
962: Input Parameter:
963: . prob - The `PetscDS` object
965: Output Parameter:
966: . dim - The total problem dimension
968: Level: beginner
970: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
971: @*/
972: PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim)
973: {
974: PetscFunctionBegin;
976: PetscCall(PetscDSSetUp(prob));
977: PetscAssertPointer(dim, 2);
978: *dim = prob->totDim;
979: PetscFunctionReturn(PETSC_SUCCESS);
980: }
982: /*@
983: PetscDSGetTotalComponents - Returns the total number of components in this system
985: Not Collective
987: Input Parameter:
988: . prob - The `PetscDS` object
990: Output Parameter:
991: . Nc - The total number of components
993: Level: beginner
995: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
996: @*/
997: PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc)
998: {
999: PetscFunctionBegin;
1001: PetscCall(PetscDSSetUp(prob));
1002: PetscAssertPointer(Nc, 2);
1003: *Nc = prob->totComp;
1004: PetscFunctionReturn(PETSC_SUCCESS);
1005: }
1007: /*@
1008: PetscDSGetDiscretization - Returns the discretization object for the given field
1010: Not Collective
1012: Input Parameters:
1013: + prob - The `PetscDS` object
1014: - f - The field number
1016: Output Parameter:
1017: . disc - The discretization object
1019: Level: beginner
1021: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1022: @*/
1023: PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc)
1024: {
1025: PetscFunctionBeginHot;
1027: PetscAssertPointer(disc, 3);
1028: 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);
1029: *disc = prob->disc[f];
1030: PetscFunctionReturn(PETSC_SUCCESS);
1031: }
1033: /*@
1034: PetscDSSetDiscretization - Sets the discretization object for the given field
1036: Not Collective
1038: Input Parameters:
1039: + prob - The `PetscDS` object
1040: . f - The field number
1041: - disc - The discretization object
1043: Level: beginner
1045: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSGetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1046: @*/
1047: PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc)
1048: {
1049: PetscFunctionBegin;
1051: if (disc) PetscAssertPointer(disc, 3);
1052: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1053: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
1054: PetscCall(PetscObjectDereference(prob->disc[f]));
1055: prob->disc[f] = disc;
1056: PetscCall(PetscObjectReference(disc));
1057: if (disc) {
1058: PetscClassId id;
1060: PetscCall(PetscObjectGetClassId(disc, &id));
1061: if (id == PETSCFE_CLASSID) {
1062: PetscCall(PetscDSSetImplicit(prob, f, PETSC_TRUE));
1063: } else if (id == PETSCFV_CLASSID) {
1064: PetscCall(PetscDSSetImplicit(prob, f, PETSC_FALSE));
1065: }
1066: PetscCall(PetscDSSetJetDegree(prob, f, 1));
1067: }
1068: PetscFunctionReturn(PETSC_SUCCESS);
1069: }
1071: /*@
1072: PetscDSGetWeakForm - Returns the weak form object
1074: Not Collective
1076: Input Parameter:
1077: . ds - The `PetscDS` object
1079: Output Parameter:
1080: . wf - The weak form object
1082: Level: beginner
1084: .seealso: `PetscWeakForm`, `PetscDSSetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1085: @*/
1086: PetscErrorCode PetscDSGetWeakForm(PetscDS ds, PetscWeakForm *wf)
1087: {
1088: PetscFunctionBegin;
1090: PetscAssertPointer(wf, 2);
1091: *wf = ds->wf;
1092: PetscFunctionReturn(PETSC_SUCCESS);
1093: }
1095: /*@
1096: PetscDSSetWeakForm - Sets the weak form object
1098: Not Collective
1100: Input Parameters:
1101: + ds - The `PetscDS` object
1102: - wf - The weak form object
1104: Level: beginner
1106: .seealso: `PetscWeakForm`, `PetscDSGetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1107: @*/
1108: PetscErrorCode PetscDSSetWeakForm(PetscDS ds, PetscWeakForm wf)
1109: {
1110: PetscFunctionBegin;
1113: PetscCall(PetscObjectDereference((PetscObject)ds->wf));
1114: ds->wf = wf;
1115: PetscCall(PetscObjectReference((PetscObject)wf));
1116: PetscCall(PetscWeakFormSetNumFields(wf, ds->Nf));
1117: PetscFunctionReturn(PETSC_SUCCESS);
1118: }
1120: /*@
1121: PetscDSAddDiscretization - Adds a discretization object
1123: Not Collective
1125: Input Parameters:
1126: + prob - The `PetscDS` object
1127: - disc - The boundary discretization object
1129: Level: beginner
1131: .seealso: `PetscWeakForm`, `PetscDSGetDiscretization()`, `PetscDSSetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1132: @*/
1133: PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc)
1134: {
1135: PetscFunctionBegin;
1136: PetscCall(PetscDSSetDiscretization(prob, prob->Nf, disc));
1137: PetscFunctionReturn(PETSC_SUCCESS);
1138: }
1140: /*@
1141: PetscDSGetQuadrature - Returns the quadrature, which must agree for all fields in the `PetscDS`
1143: Not Collective
1145: Input Parameter:
1146: . prob - The `PetscDS` object
1148: Output Parameter:
1149: . q - The quadrature object
1151: Level: intermediate
1153: .seealso: `PetscDS`, `PetscQuadrature`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1154: @*/
1155: PetscErrorCode PetscDSGetQuadrature(PetscDS prob, PetscQuadrature *q)
1156: {
1157: PetscObject obj;
1158: PetscClassId id;
1160: PetscFunctionBegin;
1161: *q = NULL;
1162: if (!prob->Nf) PetscFunctionReturn(PETSC_SUCCESS);
1163: PetscCall(PetscDSGetDiscretization(prob, 0, &obj));
1164: PetscCall(PetscObjectGetClassId(obj, &id));
1165: if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetQuadrature((PetscFE)obj, q));
1166: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetQuadrature((PetscFV)obj, q));
1167: else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
1168: PetscFunctionReturn(PETSC_SUCCESS);
1169: }
1171: /*@
1172: PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for `TSARKIMEX`
1174: Not Collective
1176: Input Parameters:
1177: + prob - The `PetscDS` object
1178: - f - The field number
1180: Output Parameter:
1181: . implicit - The flag indicating what kind of solve to use for this field
1183: Level: developer
1185: .seealso: `TSARKIMEX`, `PetscDS`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1186: @*/
1187: PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit)
1188: {
1189: PetscFunctionBegin;
1191: PetscAssertPointer(implicit, 3);
1192: 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);
1193: *implicit = prob->implicit[f];
1194: PetscFunctionReturn(PETSC_SUCCESS);
1195: }
1197: /*@
1198: PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for `TSARKIMEX`
1200: Not Collective
1202: Input Parameters:
1203: + prob - The `PetscDS` object
1204: . f - The field number
1205: - implicit - The flag indicating what kind of solve to use for this field
1207: Level: developer
1209: .seealso: `TSARKIMEX`, `PetscDSGetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1210: @*/
1211: PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit)
1212: {
1213: PetscFunctionBegin;
1215: 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);
1216: prob->implicit[f] = implicit;
1217: PetscFunctionReturn(PETSC_SUCCESS);
1218: }
1220: /*@
1221: PetscDSGetJetDegree - Returns the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate.
1223: Not Collective
1225: Input Parameters:
1226: + ds - The `PetscDS` object
1227: - f - The field number
1229: Output Parameter:
1230: . k - The highest derivative we need to tabulate
1232: Level: developer
1234: .seealso: `PetscDS`, `PetscDSSetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1235: @*/
1236: PetscErrorCode PetscDSGetJetDegree(PetscDS ds, PetscInt f, PetscInt *k)
1237: {
1238: PetscFunctionBegin;
1240: PetscAssertPointer(k, 3);
1241: 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);
1242: *k = ds->jetDegree[f];
1243: PetscFunctionReturn(PETSC_SUCCESS);
1244: }
1246: /*@
1247: PetscDSSetJetDegree - Set the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate.
1249: Not Collective
1251: Input Parameters:
1252: + ds - The `PetscDS` object
1253: . f - The field number
1254: - k - The highest derivative we need to tabulate
1256: Level: developer
1258: .seealso: `PetscDS`, `PetscDSGetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1259: @*/
1260: PetscErrorCode PetscDSSetJetDegree(PetscDS ds, PetscInt f, PetscInt k)
1261: {
1262: PetscFunctionBegin;
1264: 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);
1265: ds->jetDegree[f] = k;
1266: PetscFunctionReturn(PETSC_SUCCESS);
1267: }
1269: /*@C
1270: PetscDSGetObjective - Get the pointwise objective function for a given test field
1272: Not Collective
1274: Input Parameters:
1275: + ds - The `PetscDS`
1276: - f - The test field number
1278: Output Parameter:
1279: . obj - integrand for the test function term
1281: Calling sequence of `obj`:
1282: + dim - the spatial dimension
1283: . Nf - the number of fields
1284: . NfAux - the number of auxiliary fields
1285: . uOff - the offset into u[] and u_t[] for each field
1286: . uOff_x - the offset into u_x[] for each field
1287: . u - each field evaluated at the current point
1288: . u_t - the time derivative of each field evaluated at the current point
1289: . u_x - the gradient of each field evaluated at the current point
1290: . aOff - the offset into a[] and a_t[] for each auxiliary field
1291: . aOff_x - the offset into a_x[] for each auxiliary field
1292: . a - each auxiliary field evaluated at the current point
1293: . a_t - the time derivative of each auxiliary field evaluated at the current point
1294: . a_x - the gradient of auxiliary each field evaluated at the current point
1295: . t - current time
1296: . x - coordinates of the current point
1297: . numConstants - number of constant parameters
1298: . constants - constant parameters
1299: - obj - output values at the current point
1301: Level: intermediate
1303: Note:
1304: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi\,\mathrm{obj}(u, u_t, \nabla u, x, t)$
1306: .seealso: `PetscDS`, `PetscDSSetObjective()`, `PetscDSGetResidual()`
1307: @*/
1308: PetscErrorCode PetscDSGetObjective(PetscDS ds, PetscInt f, void (**obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
1309: {
1310: PetscPointFunc *tmp;
1311: PetscInt n;
1313: PetscFunctionBegin;
1315: PetscAssertPointer(obj, 3);
1316: 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);
1317: PetscCall(PetscWeakFormGetObjective(ds->wf, NULL, 0, f, 0, &n, &tmp));
1318: *obj = tmp ? tmp[0] : NULL;
1319: PetscFunctionReturn(PETSC_SUCCESS);
1320: }
1322: /*@C
1323: PetscDSSetObjective - Set the pointwise objective function for a given test field
1325: Not Collective
1327: Input Parameters:
1328: + ds - The `PetscDS`
1329: . f - The test field number
1330: - obj - integrand for the test function term
1332: Calling sequence of `obj`:
1333: + dim - the spatial dimension
1334: . Nf - the number of fields
1335: . NfAux - the number of auxiliary fields
1336: . uOff - the offset into u[] and u_t[] for each field
1337: . uOff_x - the offset into u_x[] for each field
1338: . u - each field evaluated at the current point
1339: . u_t - the time derivative of each field evaluated at the current point
1340: . u_x - the gradient of each field evaluated at the current point
1341: . aOff - the offset into a[] and a_t[] for each auxiliary field
1342: . aOff_x - the offset into a_x[] for each auxiliary field
1343: . a - each auxiliary field evaluated at the current point
1344: . a_t - the time derivative of each auxiliary field evaluated at the current point
1345: . a_x - the gradient of auxiliary each field evaluated at the current point
1346: . t - current time
1347: . x - coordinates of the current point
1348: . numConstants - number of constant parameters
1349: . constants - constant parameters
1350: - obj - output values at the current point
1352: Level: intermediate
1354: Note:
1355: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi\,\mathrm{obj}(u, u_t, \nabla u, x, t)$
1357: .seealso: `PetscDS`, `PetscDSGetObjective()`, `PetscDSSetResidual()`
1358: @*/
1359: PetscErrorCode PetscDSSetObjective(PetscDS ds, PetscInt f, void (*obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
1360: {
1361: PetscFunctionBegin;
1364: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1365: PetscCall(PetscWeakFormSetIndexObjective(ds->wf, NULL, 0, f, 0, 0, obj));
1366: PetscFunctionReturn(PETSC_SUCCESS);
1367: }
1369: /*@C
1370: PetscDSGetResidual - Get the pointwise residual function for a given test field
1372: Not Collective
1374: Input Parameters:
1375: + ds - The `PetscDS`
1376: - f - The test field number
1378: Output Parameters:
1379: + f0 - integrand for the test function term
1380: - f1 - integrand for the test function gradient term
1382: Calling sequence of `f0`:
1383: + dim - the spatial dimension
1384: . Nf - the number of fields
1385: . NfAux - the number of auxiliary fields
1386: . uOff - the offset into u[] and u_t[] for each field
1387: . uOff_x - the offset into u_x[] for each field
1388: . u - each field evaluated at the current point
1389: . u_t - the time derivative of each field evaluated at the current point
1390: . u_x - the gradient of each field evaluated at the current point
1391: . aOff - the offset into a[] and a_t[] for each auxiliary field
1392: . aOff_x - the offset into a_x[] for each auxiliary field
1393: . a - each auxiliary field evaluated at the current point
1394: . a_t - the time derivative of each auxiliary field evaluated at the current point
1395: . a_x - the gradient of auxiliary each field evaluated at the current point
1396: . t - current time
1397: . x - coordinates of the current point
1398: . numConstants - number of constant parameters
1399: . constants - constant parameters
1400: - f0 - output values at the current point
1402: Level: intermediate
1404: Note:
1405: `f1` has an identical form and is omitted for brevity.
1407: 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)$
1409: .seealso: `PetscDS`, `PetscDSSetResidual()`
1410: @*/
1411: PetscErrorCode PetscDSGetResidual(PetscDS ds, PetscInt f, void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (**f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1412: {
1413: PetscPointFunc *tmp0, *tmp1;
1414: PetscInt n0, n1;
1416: PetscFunctionBegin;
1418: 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);
1419: PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1));
1420: *f0 = tmp0 ? tmp0[0] : NULL;
1421: *f1 = tmp1 ? tmp1[0] : NULL;
1422: PetscFunctionReturn(PETSC_SUCCESS);
1423: }
1425: /*@C
1426: PetscDSSetResidual - Set the pointwise residual function for a given test field
1428: Not Collective
1430: Input Parameters:
1431: + ds - The `PetscDS`
1432: . f - The test field number
1433: . f0 - integrand for the test function term
1434: - f1 - integrand for the test function gradient term
1436: Calling sequence of `f0`:
1437: + dim - the spatial dimension
1438: . Nf - the number of fields
1439: . NfAux - the number of auxiliary fields
1440: . uOff - the offset into u[] and u_t[] for each field
1441: . uOff_x - the offset into u_x[] for each field
1442: . u - each field evaluated at the current point
1443: . u_t - the time derivative of each field evaluated at the current point
1444: . u_x - the gradient of each field evaluated at the current point
1445: . aOff - the offset into a[] and a_t[] for each auxiliary field
1446: . aOff_x - the offset into a_x[] for each auxiliary field
1447: . a - each auxiliary field evaluated at the current point
1448: . a_t - the time derivative of each auxiliary field evaluated at the current point
1449: . a_x - the gradient of auxiliary each field evaluated at the current point
1450: . t - current time
1451: . x - coordinates of the current point
1452: . numConstants - number of constant parameters
1453: . constants - constant parameters
1454: - f0 - output values at the current point
1456: Level: intermediate
1458: Note:
1459: `f1` has an identical form and is omitted for brevity.
1461: 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)$
1463: .seealso: `PetscDS`, `PetscDSGetResidual()`
1464: @*/
1465: PetscErrorCode PetscDSSetResidual(PetscDS ds, PetscInt f, void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (*f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1466: {
1467: PetscFunctionBegin;
1471: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1472: PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1));
1473: PetscFunctionReturn(PETSC_SUCCESS);
1474: }
1476: /*@C
1477: PetscDSGetRHSResidual - Get the pointwise RHS residual function for explicit timestepping for a given test field
1479: Not Collective
1481: Input Parameters:
1482: + ds - The `PetscDS`
1483: - f - The test field number
1485: Output Parameters:
1486: + f0 - integrand for the test function term
1487: - f1 - integrand for the test function gradient term
1489: Calling sequence of `f0`:
1490: + dim - the spatial dimension
1491: . Nf - the number of fields
1492: . NfAux - the number of auxiliary fields
1493: . uOff - the offset into u[] and u_t[] for each field
1494: . uOff_x - the offset into u_x[] for each field
1495: . u - each field evaluated at the current point
1496: . u_t - the time derivative of each field evaluated at the current point
1497: . u_x - the gradient of each field evaluated at the current point
1498: . aOff - the offset into a[] and a_t[] for each auxiliary field
1499: . aOff_x - the offset into a_x[] for each auxiliary field
1500: . a - each auxiliary field evaluated at the current point
1501: . a_t - the time derivative of each auxiliary field evaluated at the current point
1502: . a_x - the gradient of auxiliary each field evaluated at the current point
1503: . t - current time
1504: . x - coordinates of the current point
1505: . numConstants - number of constant parameters
1506: . constants - constant parameters
1507: - f0 - output values at the current point
1509: Level: intermediate
1511: Note:
1512: `f1` has an identical form and is omitted for brevity.
1514: 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)$
1516: .seealso: `PetscDS`, `PetscDSSetRHSResidual()`
1517: @*/
1518: PetscErrorCode PetscDSGetRHSResidual(PetscDS ds, PetscInt f, void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (**f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1519: {
1520: PetscPointFunc *tmp0, *tmp1;
1521: PetscInt n0, n1;
1523: PetscFunctionBegin;
1525: 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);
1526: PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 100, &n0, &tmp0, &n1, &tmp1));
1527: *f0 = tmp0 ? tmp0[0] : NULL;
1528: *f1 = tmp1 ? tmp1[0] : NULL;
1529: PetscFunctionReturn(PETSC_SUCCESS);
1530: }
1532: /*@C
1533: PetscDSSetRHSResidual - Set the pointwise residual function for explicit timestepping for a given test field
1535: Not Collective
1537: Input Parameters:
1538: + ds - The `PetscDS`
1539: . f - The test field number
1540: . f0 - integrand for the test function term
1541: - f1 - integrand for the test function gradient term
1543: Calling sequence for the callbacks `f0`:
1544: + dim - the spatial dimension
1545: . Nf - the number of fields
1546: . NfAux - the number of auxiliary fields
1547: . uOff - the offset into u[] and u_t[] for each field
1548: . uOff_x - the offset into u_x[] for each field
1549: . u - each field evaluated at the current point
1550: . u_t - the time derivative of each field evaluated at the current point
1551: . u_x - the gradient of each field evaluated at the current point
1552: . aOff - the offset into a[] and a_t[] for each auxiliary field
1553: . aOff_x - the offset into a_x[] for each auxiliary field
1554: . a - each auxiliary field evaluated at the current point
1555: . a_t - the time derivative of each auxiliary field evaluated at the current point
1556: . a_x - the gradient of auxiliary each field evaluated at the current point
1557: . t - current time
1558: . x - coordinates of the current point
1559: . numConstants - number of constant parameters
1560: . constants - constant parameters
1561: - f0 - output values at the current point
1563: Level: intermediate
1565: Note:
1566: `f1` has an identical form and is omitted for brevity.
1568: 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)$
1570: .seealso: `PetscDS`, `PetscDSGetResidual()`
1571: @*/
1572: PetscErrorCode PetscDSSetRHSResidual(PetscDS ds, PetscInt f, void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (*f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1573: {
1574: PetscFunctionBegin;
1578: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1579: PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 100, 0, f0, 0, f1));
1580: PetscFunctionReturn(PETSC_SUCCESS);
1581: }
1583: /*@
1584: PetscDSHasJacobian - Checks that the Jacobian functions have been set
1586: Not Collective
1588: Input Parameter:
1589: . ds - The `PetscDS`
1591: Output Parameter:
1592: . hasJac - flag that pointwise function for the Jacobian has been set
1594: Level: intermediate
1596: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1597: @*/
1598: PetscErrorCode PetscDSHasJacobian(PetscDS ds, PetscBool *hasJac)
1599: {
1600: PetscFunctionBegin;
1602: PetscCall(PetscWeakFormHasJacobian(ds->wf, hasJac));
1603: PetscFunctionReturn(PETSC_SUCCESS);
1604: }
1606: /*@C
1607: PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field
1609: Not Collective
1611: Input Parameters:
1612: + ds - The `PetscDS`
1613: . f - The test field number
1614: - g - The field number
1616: Output Parameters:
1617: + g0 - integrand for the test and basis function term
1618: . g1 - integrand for the test function and basis function gradient term
1619: . g2 - integrand for the test function gradient and basis function term
1620: - g3 - integrand for the test function gradient and basis function gradient term
1622: Calling sequence of `g0`:
1623: + dim - the spatial dimension
1624: . Nf - the number of fields
1625: . NfAux - the number of auxiliary fields
1626: . uOff - the offset into u[] and u_t[] for each field
1627: . uOff_x - the offset into u_x[] for each field
1628: . u - each field evaluated at the current point
1629: . u_t - the time derivative of each field evaluated at the current point
1630: . u_x - the gradient of each field evaluated at the current point
1631: . aOff - the offset into a[] and a_t[] for each auxiliary field
1632: . aOff_x - the offset into a_x[] for each auxiliary field
1633: . a - each auxiliary field evaluated at the current point
1634: . a_t - the time derivative of each auxiliary field evaluated at the current point
1635: . a_x - the gradient of auxiliary each field evaluated at the current point
1636: . t - current time
1637: . u_tShift - the multiplier a for dF/dU_t
1638: . x - coordinates of the current point
1639: . numConstants - number of constant parameters
1640: . constants - constant parameters
1641: - g0 - output values at the current point
1643: Level: intermediate
1645: Note:
1646: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1648: We are using a first order FEM model for the weak form\:
1650: $$
1651: \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
1652: + \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
1653: $$
1655: .seealso: `PetscDS`, `PetscDSSetJacobian()`
1656: @*/
1657: PetscErrorCode PetscDSGetJacobian(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1658: {
1659: PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1660: PetscInt n0, n1, n2, n3;
1662: PetscFunctionBegin;
1664: 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);
1665: 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);
1666: PetscCall(PetscWeakFormGetJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1667: *g0 = tmp0 ? tmp0[0] : NULL;
1668: *g1 = tmp1 ? tmp1[0] : NULL;
1669: *g2 = tmp2 ? tmp2[0] : NULL;
1670: *g3 = tmp3 ? tmp3[0] : NULL;
1671: PetscFunctionReturn(PETSC_SUCCESS);
1672: }
1674: /*@C
1675: PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields
1677: Not Collective
1679: Input Parameters:
1680: + ds - The `PetscDS`
1681: . f - The test field number
1682: . g - The field number
1683: . g0 - integrand for the test and basis function term
1684: . g1 - integrand for the test function and basis function gradient term
1685: . g2 - integrand for the test function gradient and basis function term
1686: - g3 - integrand for the test function gradient and basis function gradient term
1688: Calling sequence of `g0`:
1689: + dim - the spatial dimension
1690: . Nf - the number of fields
1691: . NfAux - the number of auxiliary fields
1692: . uOff - the offset into u[] and u_t[] for each field
1693: . uOff_x - the offset into u_x[] for each field
1694: . u - each field evaluated at the current point
1695: . u_t - the time derivative of each field evaluated at the current point
1696: . u_x - the gradient of each field evaluated at the current point
1697: . aOff - the offset into a[] and a_t[] for each auxiliary field
1698: . aOff_x - the offset into a_x[] for each auxiliary field
1699: . a - each auxiliary field evaluated at the current point
1700: . a_t - the time derivative of each auxiliary field evaluated at the current point
1701: . a_x - the gradient of auxiliary each field evaluated at the current point
1702: . t - current time
1703: . u_tShift - the multiplier a for dF/dU_t
1704: . x - coordinates of the current point
1705: . numConstants - number of constant parameters
1706: . constants - constant parameters
1707: - g0 - output values at the current point
1709: Level: intermediate
1711: Note:
1712: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1714: We are using a first order FEM model for the weak form\:
1716: $$
1717: \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
1718: + \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
1719: $$
1721: .seealso: `PetscDS`, `PetscDSGetJacobian()`
1722: @*/
1723: PetscErrorCode PetscDSSetJacobian(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1724: {
1725: PetscFunctionBegin;
1731: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1732: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1733: PetscCall(PetscWeakFormSetIndexJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1734: PetscFunctionReturn(PETSC_SUCCESS);
1735: }
1737: /*@
1738: PetscDSUseJacobianPreconditioner - Set whether to construct a Jacobian preconditioner
1740: Not Collective
1742: Input Parameters:
1743: + prob - The `PetscDS`
1744: - useJacPre - flag that enables construction of a Jacobian preconditioner
1746: Level: intermediate
1748: Developer Note:
1749: Should be called `PetscDSSetUseJacobianPreconditioner()`
1751: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1752: @*/
1753: PetscErrorCode PetscDSUseJacobianPreconditioner(PetscDS prob, PetscBool useJacPre)
1754: {
1755: PetscFunctionBegin;
1757: prob->useJacPre = useJacPre;
1758: PetscFunctionReturn(PETSC_SUCCESS);
1759: }
1761: /*@
1762: PetscDSHasJacobianPreconditioner - Checks if a Jacobian preconditioner matrix has been set
1764: Not Collective
1766: Input Parameter:
1767: . ds - The `PetscDS`
1769: Output Parameter:
1770: . hasJacPre - flag that pointwise function for Jacobian preconditioner matrix has been set
1772: Level: intermediate
1774: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1775: @*/
1776: PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS ds, PetscBool *hasJacPre)
1777: {
1778: PetscFunctionBegin;
1780: *hasJacPre = PETSC_FALSE;
1781: if (!ds->useJacPre) PetscFunctionReturn(PETSC_SUCCESS);
1782: PetscCall(PetscWeakFormHasJacobianPreconditioner(ds->wf, hasJacPre));
1783: PetscFunctionReturn(PETSC_SUCCESS);
1784: }
1786: /*@C
1787: PetscDSGetJacobianPreconditioner - Get the pointwise Jacobian preconditioner function for given test and basis field. If this is missing,
1788: the system matrix is used to build the preconditioner.
1790: Not Collective
1792: Input Parameters:
1793: + ds - The `PetscDS`
1794: . f - The test field number
1795: - g - The field number
1797: Output Parameters:
1798: + g0 - integrand for the test and basis function term
1799: . g1 - integrand for the test function and basis function gradient term
1800: . g2 - integrand for the test function gradient and basis function term
1801: - g3 - integrand for the test function gradient and basis function gradient term
1803: Calling sequence of `g0`:
1804: + dim - the spatial dimension
1805: . Nf - the number of fields
1806: . NfAux - the number of auxiliary fields
1807: . uOff - the offset into u[] and u_t[] for each field
1808: . uOff_x - the offset into u_x[] for each field
1809: . u - each field evaluated at the current point
1810: . u_t - the time derivative of each field evaluated at the current point
1811: . u_x - the gradient of each field evaluated at the current point
1812: . aOff - the offset into a[] and a_t[] for each auxiliary field
1813: . aOff_x - the offset into a_x[] for each auxiliary field
1814: . a - each auxiliary field evaluated at the current point
1815: . a_t - the time derivative of each auxiliary field evaluated at the current point
1816: . a_x - the gradient of auxiliary each field evaluated at the current point
1817: . t - current time
1818: . u_tShift - the multiplier a for dF/dU_t
1819: . x - coordinates of the current point
1820: . numConstants - number of constant parameters
1821: . constants - constant parameters
1822: - g0 - output values at the current point
1824: Level: intermediate
1826: Note:
1827: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1828: We are using a first order FEM model for the weak form\:
1830: $$
1831: \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
1832: + \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
1833: $$
1835: .seealso: `PetscDS`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1836: @*/
1837: PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1838: {
1839: PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1840: PetscInt n0, n1, n2, n3;
1842: PetscFunctionBegin;
1844: 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);
1845: 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);
1846: PetscCall(PetscWeakFormGetJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1847: *g0 = tmp0 ? tmp0[0] : NULL;
1848: *g1 = tmp1 ? tmp1[0] : NULL;
1849: *g2 = tmp2 ? tmp2[0] : NULL;
1850: *g3 = tmp3 ? tmp3[0] : NULL;
1851: PetscFunctionReturn(PETSC_SUCCESS);
1852: }
1854: /*@C
1855: PetscDSSetJacobianPreconditioner - Set the pointwise Jacobian preconditioner function for given test and basis fields.
1856: If this is missing, the system matrix is used to build the preconditioner.
1858: Not Collective
1860: Input Parameters:
1861: + ds - The `PetscDS`
1862: . f - The test field number
1863: . g - The field number
1864: . g0 - integrand for the test and basis function term
1865: . g1 - integrand for the test function and basis function gradient term
1866: . g2 - integrand for the test function gradient and basis function term
1867: - g3 - integrand for the test function gradient and basis function gradient term
1869: Calling sequence of `g0`:
1870: + dim - the spatial dimension
1871: . Nf - the number of fields
1872: . NfAux - the number of auxiliary fields
1873: . uOff - the offset into u[] and u_t[] for each field
1874: . uOff_x - the offset into u_x[] for each field
1875: . u - each field evaluated at the current point
1876: . u_t - the time derivative of each field evaluated at the current point
1877: . u_x - the gradient of each field evaluated at the current point
1878: . aOff - the offset into a[] and a_t[] for each auxiliary field
1879: . aOff_x - the offset into a_x[] for each auxiliary field
1880: . a - each auxiliary field evaluated at the current point
1881: . a_t - the time derivative of each auxiliary field evaluated at the current point
1882: . a_x - the gradient of auxiliary each field evaluated at the current point
1883: . t - current time
1884: . u_tShift - the multiplier a for dF/dU_t
1885: . x - coordinates of the current point
1886: . numConstants - number of constant parameters
1887: . constants - constant parameters
1888: - g0 - output values at the current point
1890: Level: intermediate
1892: Note:
1893: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1895: We are using a first order FEM model for the weak form\:
1897: $$
1898: \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
1899: + \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
1900: $$
1902: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobian()`
1903: @*/
1904: PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1905: {
1906: PetscFunctionBegin;
1912: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1913: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1914: PetscCall(PetscWeakFormSetIndexJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1915: PetscFunctionReturn(PETSC_SUCCESS);
1916: }
1918: /*@
1919: PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, dF/du_t, has been set
1921: Not Collective
1923: Input Parameter:
1924: . ds - The `PetscDS`
1926: Output Parameter:
1927: . hasDynJac - flag that pointwise function for dynamic Jacobian has been set
1929: Level: intermediate
1931: .seealso: `PetscDS`, `PetscDSGetDynamicJacobian()`, `PetscDSSetDynamicJacobian()`, `PetscDSGetJacobian()`
1932: @*/
1933: PetscErrorCode PetscDSHasDynamicJacobian(PetscDS ds, PetscBool *hasDynJac)
1934: {
1935: PetscFunctionBegin;
1937: PetscCall(PetscWeakFormHasDynamicJacobian(ds->wf, hasDynJac));
1938: PetscFunctionReturn(PETSC_SUCCESS);
1939: }
1941: /*@C
1942: PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, dF/du_t, function for given test and basis field
1944: Not Collective
1946: Input Parameters:
1947: + ds - The `PetscDS`
1948: . f - The test field number
1949: - g - The field number
1951: Output Parameters:
1952: + g0 - integrand for the test and basis function term
1953: . g1 - integrand for the test function and basis function gradient term
1954: . g2 - integrand for the test function gradient and basis function term
1955: - g3 - integrand for the test function gradient and basis function gradient term
1957: Calling sequence of `g0`:
1958: + dim - the spatial dimension
1959: . Nf - the number of fields
1960: . NfAux - the number of auxiliary fields
1961: . uOff - the offset into u[] and u_t[] for each field
1962: . uOff_x - the offset into u_x[] for each field
1963: . u - each field evaluated at the current point
1964: . u_t - the time derivative of each field evaluated at the current point
1965: . u_x - the gradient of each field evaluated at the current point
1966: . aOff - the offset into a[] and a_t[] for each auxiliary field
1967: . aOff_x - the offset into a_x[] for each auxiliary field
1968: . a - each auxiliary field evaluated at the current point
1969: . a_t - the time derivative of each auxiliary field evaluated at the current point
1970: . a_x - the gradient of auxiliary each field evaluated at the current point
1971: . t - current time
1972: . u_tShift - the multiplier a for dF/dU_t
1973: . x - coordinates of the current point
1974: . numConstants - number of constant parameters
1975: . constants - constant parameters
1976: - g0 - output values at the current point
1978: Level: intermediate
1980: Note:
1981: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1983: We are using a first order FEM model for the weak form\:
1985: $$
1986: \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
1987: + \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
1988: $$
1990: .seealso: `PetscDS`, `PetscDSSetJacobian()`
1991: @*/
1992: PetscErrorCode PetscDSGetDynamicJacobian(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1993: {
1994: PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1995: PetscInt n0, n1, n2, n3;
1997: PetscFunctionBegin;
1999: 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);
2000: 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);
2001: PetscCall(PetscWeakFormGetDynamicJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2002: *g0 = tmp0 ? tmp0[0] : NULL;
2003: *g1 = tmp1 ? tmp1[0] : NULL;
2004: *g2 = tmp2 ? tmp2[0] : NULL;
2005: *g3 = tmp3 ? tmp3[0] : NULL;
2006: PetscFunctionReturn(PETSC_SUCCESS);
2007: }
2009: /*@C
2010: PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, dF/du_t, function for given test and basis fields
2012: Not Collective
2014: Input Parameters:
2015: + ds - The `PetscDS`
2016: . f - The test field number
2017: . g - The field number
2018: . g0 - integrand for the test and basis function term
2019: . g1 - integrand for the test function and basis function gradient term
2020: . g2 - integrand for the test function gradient and basis function term
2021: - g3 - integrand for the test function gradient and basis function gradient term
2023: Calling sequence of `g0`:
2024: + dim - the spatial dimension
2025: . Nf - the number of fields
2026: . NfAux - the number of auxiliary fields
2027: . uOff - the offset into u[] and u_t[] for each field
2028: . uOff_x - the offset into u_x[] for each field
2029: . u - each field evaluated at the current point
2030: . u_t - the time derivative of each field evaluated at the current point
2031: . u_x - the gradient of each field evaluated at the current point
2032: . aOff - the offset into a[] and a_t[] for each auxiliary field
2033: . aOff_x - the offset into a_x[] for each auxiliary field
2034: . a - each auxiliary field evaluated at the current point
2035: . a_t - the time derivative of each auxiliary field evaluated at the current point
2036: . a_x - the gradient of auxiliary each field evaluated at the current point
2037: . t - current time
2038: . u_tShift - the multiplier a for dF/dU_t
2039: . x - coordinates of the current point
2040: . numConstants - number of constant parameters
2041: . constants - constant parameters
2042: - g0 - output values at the current point
2044: Level: intermediate
2046: Note:
2047: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2049: We are using a first order FEM model for the weak form\:
2051: $$
2052: \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
2053: + \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
2054: $$
2056: .seealso: `PetscDS`, `PetscDSGetJacobian()`
2057: @*/
2058: PetscErrorCode PetscDSSetDynamicJacobian(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
2059: {
2060: PetscFunctionBegin;
2066: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2067: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2068: PetscCall(PetscWeakFormSetIndexDynamicJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2069: PetscFunctionReturn(PETSC_SUCCESS);
2070: }
2072: /*@C
2073: PetscDSGetRiemannSolver - Returns the Riemann solver for the given field
2075: Not Collective
2077: Input Parameters:
2078: + ds - The `PetscDS` object
2079: - f - The field number
2081: Output Parameter:
2082: . r - Riemann solver
2084: Calling sequence of `r`:
2085: + dim - The spatial dimension
2086: . Nf - The number of fields
2087: . x - The coordinates at a point on the interface
2088: . n - The normal vector to the interface
2089: . uL - The state vector to the left of the interface
2090: . uR - The state vector to the right of the interface
2091: . flux - output array of flux through the interface
2092: . numConstants - number of constant parameters
2093: . constants - constant parameters
2094: - ctx - optional user context
2096: Level: intermediate
2098: .seealso: `PetscDS`, `PetscDSSetRiemannSolver()`
2099: @*/
2100: PetscErrorCode PetscDSGetRiemannSolver(PetscDS ds, PetscInt f, void (**r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscInt numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx))
2101: {
2102: PetscRiemannFunc *tmp;
2103: PetscInt n;
2105: PetscFunctionBegin;
2107: PetscAssertPointer(r, 3);
2108: 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);
2109: PetscCall(PetscWeakFormGetRiemannSolver(ds->wf, NULL, 0, f, 0, &n, &tmp));
2110: *r = tmp ? tmp[0] : NULL;
2111: PetscFunctionReturn(PETSC_SUCCESS);
2112: }
2114: /*@C
2115: PetscDSSetRiemannSolver - Sets the Riemann solver for the given field
2117: Not Collective
2119: Input Parameters:
2120: + ds - The `PetscDS` object
2121: . f - The field number
2122: - r - Riemann solver
2124: Calling sequence of `r`:
2125: + dim - The spatial dimension
2126: . Nf - The number of fields
2127: . x - The coordinates at a point on the interface
2128: . n - The normal vector to the interface
2129: . uL - The state vector to the left of the interface
2130: . uR - The state vector to the right of the interface
2131: . flux - output array of flux through the interface
2132: . numConstants - number of constant parameters
2133: . constants - constant parameters
2134: - ctx - optional user context
2136: Level: intermediate
2138: .seealso: `PetscDS`, `PetscDSGetRiemannSolver()`
2139: @*/
2140: PetscErrorCode PetscDSSetRiemannSolver(PetscDS ds, PetscInt f, void (*r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscInt numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx))
2141: {
2142: PetscFunctionBegin;
2145: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2146: PetscCall(PetscWeakFormSetIndexRiemannSolver(ds->wf, NULL, 0, f, 0, 0, r));
2147: PetscFunctionReturn(PETSC_SUCCESS);
2148: }
2150: /*@C
2151: PetscDSGetUpdate - Get the pointwise update function for a given field
2153: Not Collective
2155: Input Parameters:
2156: + ds - The `PetscDS`
2157: - f - The field number
2159: Output Parameter:
2160: . update - update function
2162: Calling sequence of `update`:
2163: + dim - the spatial dimension
2164: . Nf - the number of fields
2165: . NfAux - the number of auxiliary fields
2166: . uOff - the offset into u[] and u_t[] for each field
2167: . uOff_x - the offset into u_x[] for each field
2168: . u - each field evaluated at the current point
2169: . u_t - the time derivative of each field evaluated at the current point
2170: . u_x - the gradient of each field evaluated at the current point
2171: . aOff - the offset into a[] and a_t[] for each auxiliary field
2172: . aOff_x - the offset into a_x[] for each auxiliary field
2173: . a - each auxiliary field evaluated at the current point
2174: . a_t - the time derivative of each auxiliary field evaluated at the current point
2175: . a_x - the gradient of auxiliary each field evaluated at the current point
2176: . t - current time
2177: . x - coordinates of the current point
2178: . numConstants - number of constant parameters
2179: . constants - constant parameters
2180: - uNew - new value for field at the current point
2182: Level: intermediate
2184: .seealso: `PetscDS`, `PetscDSSetUpdate()`, `PetscDSSetResidual()`
2185: @*/
2186: PetscErrorCode PetscDSGetUpdate(PetscDS ds, PetscInt f, void (**update)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
2187: {
2188: PetscFunctionBegin;
2190: 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);
2191: if (update) {
2192: PetscAssertPointer(update, 3);
2193: *update = ds->update[f];
2194: }
2195: PetscFunctionReturn(PETSC_SUCCESS);
2196: }
2198: /*@C
2199: PetscDSSetUpdate - Set the pointwise update function for a given field
2201: Not Collective
2203: Input Parameters:
2204: + ds - The `PetscDS`
2205: . f - The field number
2206: - update - update function
2208: Calling sequence of `update`:
2209: + dim - the spatial dimension
2210: . Nf - the number of fields
2211: . NfAux - the number of auxiliary fields
2212: . uOff - the offset into u[] and u_t[] for each field
2213: . uOff_x - the offset into u_x[] for each field
2214: . u - each field evaluated at the current point
2215: . u_t - the time derivative of each field evaluated at the current point
2216: . u_x - the gradient of each field evaluated at the current point
2217: . aOff - the offset into a[] and a_t[] for each auxiliary field
2218: . aOff_x - the offset into a_x[] for each auxiliary field
2219: . a - each auxiliary field evaluated at the current point
2220: . a_t - the time derivative of each auxiliary field evaluated at the current point
2221: . a_x - the gradient of auxiliary each field evaluated at the current point
2222: . t - current time
2223: . x - coordinates of the current point
2224: . numConstants - number of constant parameters
2225: . constants - constant parameters
2226: - uNew - new field values at the current point
2228: Level: intermediate
2230: .seealso: `PetscDS`, `PetscDSGetResidual()`
2231: @*/
2232: PetscErrorCode PetscDSSetUpdate(PetscDS ds, PetscInt f, void (*update)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
2233: {
2234: PetscFunctionBegin;
2237: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2238: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
2239: ds->update[f] = update;
2240: PetscFunctionReturn(PETSC_SUCCESS);
2241: }
2243: PetscErrorCode PetscDSGetContext(PetscDS ds, PetscInt f, void *ctx)
2244: {
2245: PetscFunctionBegin;
2247: 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);
2248: PetscAssertPointer(ctx, 3);
2249: *(void **)ctx = ds->ctx[f];
2250: PetscFunctionReturn(PETSC_SUCCESS);
2251: }
2253: PetscErrorCode PetscDSSetContext(PetscDS ds, PetscInt f, void *ctx)
2254: {
2255: PetscFunctionBegin;
2257: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2258: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
2259: ds->ctx[f] = ctx;
2260: PetscFunctionReturn(PETSC_SUCCESS);
2261: }
2263: /*@C
2264: PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field
2266: Not Collective
2268: Input Parameters:
2269: + ds - The PetscDS
2270: - f - The test field number
2272: Output Parameters:
2273: + f0 - boundary integrand for the test function term
2274: - f1 - boundary integrand for the test function gradient term
2276: Calling sequence of `f0`:
2277: + dim - the spatial dimension
2278: . Nf - the number of fields
2279: . NfAux - the number of auxiliary fields
2280: . uOff - the offset into u[] and u_t[] for each field
2281: . uOff_x - the offset into u_x[] for each field
2282: . u - each field evaluated at the current point
2283: . u_t - the time derivative of each field evaluated at the current point
2284: . u_x - the gradient of each field evaluated at the current point
2285: . aOff - the offset into a[] and a_t[] for each auxiliary field
2286: . aOff_x - the offset into a_x[] for each auxiliary field
2287: . a - each auxiliary field evaluated at the current point
2288: . a_t - the time derivative of each auxiliary field evaluated at the current point
2289: . a_x - the gradient of auxiliary each field evaluated at the current point
2290: . t - current time
2291: . x - coordinates of the current point
2292: . n - unit normal at the current point
2293: . numConstants - number of constant parameters
2294: . constants - constant parameters
2295: - f0 - output values at the current point
2297: Level: intermediate
2299: Note:
2300: The calling sequence of `f1` is identical, and therefore omitted for brevity.
2302: We are using a first order FEM model for the weak form\:
2304: $$
2305: \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
2306: $$
2308: .seealso: `PetscDS`, `PetscDSSetBdResidual()`
2309: @*/
2310: PetscErrorCode PetscDSGetBdResidual(PetscDS ds, PetscInt f, void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (**f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
2311: {
2312: PetscBdPointFunc *tmp0, *tmp1;
2313: PetscInt n0, n1;
2315: PetscFunctionBegin;
2317: 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);
2318: PetscCall(PetscWeakFormGetBdResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1));
2319: *f0 = tmp0 ? tmp0[0] : NULL;
2320: *f1 = tmp1 ? tmp1[0] : NULL;
2321: PetscFunctionReturn(PETSC_SUCCESS);
2322: }
2324: /*@C
2325: PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field
2327: Not Collective
2329: Input Parameters:
2330: + ds - The `PetscDS`
2331: . f - The test field number
2332: . f0 - boundary integrand for the test function term
2333: - f1 - boundary integrand for the test function gradient term
2335: Calling sequence of `f0`:
2336: + dim - the spatial dimension
2337: . Nf - the number of fields
2338: . NfAux - the number of auxiliary fields
2339: . uOff - the offset into u[] and u_t[] for each field
2340: . uOff_x - the offset into u_x[] for each field
2341: . u - each field evaluated at the current point
2342: . u_t - the time derivative of each field evaluated at the current point
2343: . u_x - the gradient of each field evaluated at the current point
2344: . aOff - the offset into a[] and a_t[] for each auxiliary field
2345: . aOff_x - the offset into a_x[] for each auxiliary field
2346: . a - each auxiliary field evaluated at the current point
2347: . a_t - the time derivative of each auxiliary field evaluated at the current point
2348: . a_x - the gradient of auxiliary each field evaluated at the current point
2349: . t - current time
2350: . x - coordinates of the current point
2351: . n - unit normal at the current point
2352: . numConstants - number of constant parameters
2353: . constants - constant parameters
2354: - f0 - output values at the current point
2356: Level: intermediate
2358: Note:
2359: The calling sequence of `f1` is identical, and therefore omitted for brevity.
2361: We are using a first order FEM model for the weak form\:
2363: $$
2364: \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
2365: $$
2367: .seealso: `PetscDS`, `PetscDSGetBdResidual()`
2368: @*/
2369: PetscErrorCode PetscDSSetBdResidual(PetscDS ds, PetscInt f, void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (*f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
2370: {
2371: PetscFunctionBegin;
2373: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2374: PetscCall(PetscWeakFormSetIndexBdResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1));
2375: PetscFunctionReturn(PETSC_SUCCESS);
2376: }
2378: /*@
2379: PetscDSHasBdJacobian - Indicates that boundary Jacobian functions have been set
2381: Not Collective
2383: Input Parameter:
2384: . ds - The `PetscDS`
2386: Output Parameter:
2387: . hasBdJac - flag that pointwise function for the boundary Jacobian has been set
2389: Level: intermediate
2391: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2392: @*/
2393: PetscErrorCode PetscDSHasBdJacobian(PetscDS ds, PetscBool *hasBdJac)
2394: {
2395: PetscFunctionBegin;
2397: PetscAssertPointer(hasBdJac, 2);
2398: PetscCall(PetscWeakFormHasBdJacobian(ds->wf, hasBdJac));
2399: PetscFunctionReturn(PETSC_SUCCESS);
2400: }
2402: /*@C
2403: PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field
2405: Not Collective
2407: Input Parameters:
2408: + ds - The `PetscDS`
2409: . f - The test field number
2410: - g - The field number
2412: Output Parameters:
2413: + g0 - integrand for the test and basis function term
2414: . g1 - integrand for the test function and basis function gradient term
2415: . g2 - integrand for the test function gradient and basis function term
2416: - g3 - integrand for the test function gradient and basis function gradient term
2418: Calling sequence of `g0`:
2419: + dim - the spatial dimension
2420: . Nf - the number of fields
2421: . NfAux - the number of auxiliary fields
2422: . uOff - the offset into u[] and u_t[] for each field
2423: . uOff_x - the offset into u_x[] for each field
2424: . u - each field evaluated at the current point
2425: . u_t - the time derivative of each field evaluated at the current point
2426: . u_x - the gradient of each field evaluated at the current point
2427: . aOff - the offset into a[] and a_t[] for each auxiliary field
2428: . aOff_x - the offset into a_x[] for each auxiliary field
2429: . a - each auxiliary field evaluated at the current point
2430: . a_t - the time derivative of each auxiliary field evaluated at the current point
2431: . a_x - the gradient of auxiliary each field evaluated at the current point
2432: . t - current time
2433: . u_tShift - the multiplier a for dF/dU_t
2434: . x - coordinates of the current point
2435: . n - normal at the current point
2436: . numConstants - number of constant parameters
2437: . constants - constant parameters
2438: - g0 - output values at the current point
2440: Level: intermediate
2442: Note:
2443: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2445: We are using a first order FEM model for the weak form\:
2447: $$
2448: \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
2449: + \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
2450: $$
2452: .seealso: `PetscDS`, `PetscDSSetBdJacobian()`
2453: @*/
2454: PetscErrorCode PetscDSGetBdJacobian(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
2455: {
2456: PetscBdPointJac *tmp0, *tmp1, *tmp2, *tmp3;
2457: PetscInt n0, n1, n2, n3;
2459: PetscFunctionBegin;
2461: 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);
2462: 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);
2463: PetscCall(PetscWeakFormGetBdJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2464: *g0 = tmp0 ? tmp0[0] : NULL;
2465: *g1 = tmp1 ? tmp1[0] : NULL;
2466: *g2 = tmp2 ? tmp2[0] : NULL;
2467: *g3 = tmp3 ? tmp3[0] : NULL;
2468: PetscFunctionReturn(PETSC_SUCCESS);
2469: }
2471: /*@C
2472: PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field
2474: Not Collective
2476: Input Parameters:
2477: + ds - The PetscDS
2478: . f - The test field number
2479: . g - The field number
2480: . g0 - integrand for the test and basis function term
2481: . g1 - integrand for the test function and basis function gradient term
2482: . g2 - integrand for the test function gradient and basis function term
2483: - g3 - integrand for the test function gradient and basis function gradient term
2485: Calling sequence of `g0`:
2486: + dim - the spatial dimension
2487: . Nf - the number of fields
2488: . NfAux - the number of auxiliary fields
2489: . uOff - the offset into u[] and u_t[] for each field
2490: . uOff_x - the offset into u_x[] for each field
2491: . u - each field evaluated at the current point
2492: . u_t - the time derivative of each field evaluated at the current point
2493: . u_x - the gradient of each field evaluated at the current point
2494: . aOff - the offset into a[] and a_t[] for each auxiliary field
2495: . aOff_x - the offset into a_x[] for each auxiliary field
2496: . a - each auxiliary field evaluated at the current point
2497: . a_t - the time derivative of each auxiliary field evaluated at the current point
2498: . a_x - the gradient of auxiliary each field evaluated at the current point
2499: . t - current time
2500: . u_tShift - the multiplier a for dF/dU_t
2501: . x - coordinates of the current point
2502: . n - normal at the current point
2503: . numConstants - number of constant parameters
2504: . constants - constant parameters
2505: - g0 - output values at the current point
2507: Level: intermediate
2509: Note:
2510: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2512: We are using a first order FEM model for the weak form\:
2514: $$
2515: \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
2516: + \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
2517: $$
2519: .seealso: `PetscDS`, `PetscDSGetBdJacobian()`
2520: @*/
2521: PetscErrorCode PetscDSSetBdJacobian(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
2522: {
2523: PetscFunctionBegin;
2529: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2530: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2531: PetscCall(PetscWeakFormSetIndexBdJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2532: PetscFunctionReturn(PETSC_SUCCESS);
2533: }
2535: /*@
2536: PetscDSHasBdJacobianPreconditioner - Signals that boundary Jacobian preconditioner functions have been set
2538: Not Collective
2540: Input Parameter:
2541: . ds - The `PetscDS`
2543: Output Parameter:
2544: . hasBdJacPre - flag that pointwise function for the boundary Jacobian preconditioner has been set
2546: Level: intermediate
2548: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2549: @*/
2550: PetscErrorCode PetscDSHasBdJacobianPreconditioner(PetscDS ds, PetscBool *hasBdJacPre)
2551: {
2552: PetscFunctionBegin;
2554: PetscAssertPointer(hasBdJacPre, 2);
2555: PetscCall(PetscWeakFormHasBdJacobianPreconditioner(ds->wf, hasBdJacPre));
2556: PetscFunctionReturn(PETSC_SUCCESS);
2557: }
2559: /*@C
2560: PetscDSGetBdJacobianPreconditioner - Get the pointwise boundary Jacobian preconditioner function for given test and basis field
2562: Not Collective; No Fortran Support
2564: Input Parameters:
2565: + ds - The `PetscDS`
2566: . f - The test field number
2567: - g - The field number
2569: Output Parameters:
2570: + g0 - integrand for the test and basis function term
2571: . g1 - integrand for the test function and basis function gradient term
2572: . g2 - integrand for the test function gradient and basis function term
2573: - g3 - integrand for the test function gradient and basis function gradient term
2575: Calling sequence of `g0`:
2576: + dim - the spatial dimension
2577: . Nf - the number of fields
2578: . NfAux - the number of auxiliary fields
2579: . uOff - the offset into u[] and u_t[] for each field
2580: . uOff_x - the offset into u_x[] for each field
2581: . u - each field evaluated at the current point
2582: . u_t - the time derivative of each field evaluated at the current point
2583: . u_x - the gradient of each field evaluated at the current point
2584: . aOff - the offset into a[] and a_t[] for each auxiliary field
2585: . aOff_x - the offset into a_x[] for each auxiliary field
2586: . a - each auxiliary field evaluated at the current point
2587: . a_t - the time derivative of each auxiliary field evaluated at the current point
2588: . a_x - the gradient of auxiliary each field evaluated at the current point
2589: . t - current time
2590: . u_tShift - the multiplier a for dF/dU_t
2591: . x - coordinates of the current point
2592: . n - normal at the current point
2593: . numConstants - number of constant parameters
2594: . constants - constant parameters
2595: - g0 - output values at the current point
2597: Level: intermediate
2599: Note:
2600: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2602: We are using a first order FEM model for the weak form\:
2604: $$
2605: \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
2606: + \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
2607: $$
2609: .seealso: `PetscDS`, `PetscDSSetBdJacobianPreconditioner()`
2610: @*/
2611: PetscErrorCode PetscDSGetBdJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
2612: {
2613: PetscBdPointJac *tmp0, *tmp1, *tmp2, *tmp3;
2614: PetscInt n0, n1, n2, n3;
2616: PetscFunctionBegin;
2618: 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);
2619: 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);
2620: PetscCall(PetscWeakFormGetBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2621: *g0 = tmp0 ? tmp0[0] : NULL;
2622: *g1 = tmp1 ? tmp1[0] : NULL;
2623: *g2 = tmp2 ? tmp2[0] : NULL;
2624: *g3 = tmp3 ? tmp3[0] : NULL;
2625: PetscFunctionReturn(PETSC_SUCCESS);
2626: }
2628: /*@C
2629: PetscDSSetBdJacobianPreconditioner - Set the pointwise boundary Jacobian preconditioner function for given test and basis field
2631: Not Collective; No Fortran Support
2633: Input Parameters:
2634: + ds - The `PetscDS`
2635: . f - The test field number
2636: . g - The field number
2637: . g0 - integrand for the test and basis function term
2638: . g1 - integrand for the test function and basis function gradient term
2639: . g2 - integrand for the test function gradient and basis function term
2640: - g3 - integrand for the test function gradient and basis function gradient term
2642: Calling sequence of `g0':
2643: + dim - the spatial dimension
2644: . Nf - the number of fields
2645: . NfAux - the number of auxiliary fields
2646: . uOff - the offset into u[] and u_t[] for each field
2647: . uOff_x - the offset into u_x[] for each field
2648: . u - each field evaluated at the current point
2649: . u_t - the time derivative of each field evaluated at the current point
2650: . u_x - the gradient of each field evaluated at the current point
2651: . aOff - the offset into a[] and a_t[] for each auxiliary field
2652: . aOff_x - the offset into a_x[] for each auxiliary field
2653: . a - each auxiliary field evaluated at the current point
2654: . a_t - the time derivative of each auxiliary field evaluated at the current point
2655: . a_x - the gradient of auxiliary each field evaluated at the current point
2656: . t - current time
2657: . u_tShift - the multiplier a for dF/dU_t
2658: . x - coordinates of the current point
2659: . n - normal at the current point
2660: . numConstants - number of constant parameters
2661: . constants - constant parameters
2662: - g0 - output values at the current point
2664: Level: intermediate
2666: Note:
2667: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2669: We are using a first order FEM model for the weak form\:
2671: $$
2672: \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
2673: + \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
2674: $$
2676: .seealso: `PetscDS`, `PetscDSGetBdJacobianPreconditioner()`
2677: @*/
2678: PetscErrorCode PetscDSSetBdJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (*g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (*g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
2679: {
2680: PetscFunctionBegin;
2686: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2687: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2688: PetscCall(PetscWeakFormSetIndexBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2689: PetscFunctionReturn(PETSC_SUCCESS);
2690: }
2692: /*@C
2693: PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field
2695: Not Collective
2697: Input Parameters:
2698: + prob - The PetscDS
2699: - f - The test field number
2701: Output Parameters:
2702: + sol - exact solution for the test field
2703: - ctx - exact solution context
2705: Calling sequence of `exactSol`:
2706: + dim - the spatial dimension
2707: . t - current time
2708: . x - coordinates of the current point
2709: . Nc - the number of field components
2710: . u - the solution field evaluated at the current point
2711: - ctx - a user context
2713: Level: intermediate
2715: .seealso: `PetscDS`, `PetscDSSetExactSolution()`, `PetscDSGetExactSolutionTimeDerivative()`
2716: @*/
2717: PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx)
2718: {
2719: PetscFunctionBegin;
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: if (sol) {
2723: PetscAssertPointer(sol, 3);
2724: *sol = prob->exactSol[f];
2725: }
2726: if (ctx) {
2727: PetscAssertPointer(ctx, 4);
2728: *ctx = prob->exactCtx[f];
2729: }
2730: PetscFunctionReturn(PETSC_SUCCESS);
2731: }
2733: /*@C
2734: PetscDSSetExactSolution - Set the pointwise exact solution function for a given test field
2736: Not Collective
2738: Input Parameters:
2739: + prob - The `PetscDS`
2740: . f - The test field number
2741: . sol - solution function for the test fields
2742: - ctx - solution context or `NULL`
2744: Calling sequence of `sol`:
2745: + dim - the spatial dimension
2746: . t - current time
2747: . x - coordinates of the current point
2748: . Nc - the number of field components
2749: . u - the solution field evaluated at the current point
2750: - ctx - a user context
2752: Level: intermediate
2754: .seealso: `PetscDS`, `PetscDSGetExactSolution()`
2755: @*/
2756: PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx)
2757: {
2758: PetscFunctionBegin;
2760: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2761: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
2762: if (sol) {
2764: prob->exactSol[f] = sol;
2765: }
2766: if (ctx) {
2768: prob->exactCtx[f] = ctx;
2769: }
2770: PetscFunctionReturn(PETSC_SUCCESS);
2771: }
2773: /*@C
2774: PetscDSGetExactSolutionTimeDerivative - Get the pointwise time derivative of the exact solution function for a given test field
2776: Not Collective
2778: Input Parameters:
2779: + prob - The `PetscDS`
2780: - f - The test field number
2782: Output Parameters:
2783: + sol - time derivative of the exact solution for the test field
2784: - ctx - time derivative of the exact solution context
2786: Calling sequence of `exactSol`:
2787: + dim - the spatial dimension
2788: . t - current time
2789: . x - coordinates of the current point
2790: . Nc - the number of field components
2791: . u - the solution field evaluated at the current point
2792: - ctx - a user context
2794: Level: intermediate
2796: .seealso: `PetscDS`, `PetscDSSetExactSolutionTimeDerivative()`, `PetscDSGetExactSolution()`
2797: @*/
2798: PetscErrorCode PetscDSGetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx)
2799: {
2800: PetscFunctionBegin;
2802: 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);
2803: if (sol) {
2804: PetscAssertPointer(sol, 3);
2805: *sol = prob->exactSol_t[f];
2806: }
2807: if (ctx) {
2808: PetscAssertPointer(ctx, 4);
2809: *ctx = prob->exactCtx_t[f];
2810: }
2811: PetscFunctionReturn(PETSC_SUCCESS);
2812: }
2814: /*@C
2815: PetscDSSetExactSolutionTimeDerivative - Set the pointwise time derivative of the exact solution function for a given test field
2817: Not Collective
2819: Input Parameters:
2820: + prob - The `PetscDS`
2821: . f - The test field number
2822: . sol - time derivative of the solution function for the test fields
2823: - ctx - time derivative of the solution context or `NULL`
2825: Calling sequence of `sol`:
2826: + dim - the spatial dimension
2827: . t - current time
2828: . x - coordinates of the current point
2829: . Nc - the number of field components
2830: . u - the solution field evaluated at the current point
2831: - ctx - a user context
2833: Level: intermediate
2835: .seealso: `PetscDS`, `PetscDSGetExactSolutionTimeDerivative()`, `PetscDSSetExactSolution()`
2836: @*/
2837: PetscErrorCode PetscDSSetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx)
2838: {
2839: PetscFunctionBegin;
2841: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2842: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
2843: if (sol) {
2845: prob->exactSol_t[f] = sol;
2846: }
2847: if (ctx) {
2849: prob->exactCtx_t[f] = ctx;
2850: }
2851: PetscFunctionReturn(PETSC_SUCCESS);
2852: }
2854: /*@C
2855: PetscDSGetConstants - Returns the array of constants passed to point functions
2857: Not Collective
2859: Input Parameter:
2860: . ds - The `PetscDS` object
2862: Output Parameters:
2863: + numConstants - The number of constants
2864: - constants - The array of constants, NULL if there are none
2866: Level: intermediate
2868: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSCreate()`
2869: @*/
2870: PetscErrorCode PetscDSGetConstants(PetscDS ds, PetscInt *numConstants, const PetscScalar *constants[])
2871: {
2872: PetscFunctionBegin;
2874: if (numConstants) {
2875: PetscAssertPointer(numConstants, 2);
2876: *numConstants = ds->numConstants;
2877: }
2878: if (constants) {
2879: PetscAssertPointer(constants, 3);
2880: *constants = ds->constants;
2881: }
2882: PetscFunctionReturn(PETSC_SUCCESS);
2883: }
2885: /*@C
2886: PetscDSSetConstants - Set the array of constants passed to point functions
2888: Not Collective
2890: Input Parameters:
2891: + ds - The `PetscDS` object
2892: . numConstants - The number of constants
2893: - constants - The array of constants, `NULL` if there are none
2895: Level: intermediate
2897: .seealso: `PetscDS`, `PetscDSGetConstants()`, `PetscDSCreate()`
2898: @*/
2899: PetscErrorCode PetscDSSetConstants(PetscDS ds, PetscInt numConstants, PetscScalar constants[])
2900: {
2901: PetscFunctionBegin;
2903: if (numConstants != ds->numConstants) {
2904: PetscCall(PetscFree(ds->constants));
2905: ds->numConstants = numConstants;
2906: PetscCall(PetscMalloc1(ds->numConstants + ds->numFuncConstants, &ds->constants));
2907: }
2908: if (ds->numConstants) {
2909: PetscAssertPointer(constants, 3);
2910: PetscCall(PetscArraycpy(ds->constants, constants, ds->numConstants));
2911: }
2912: PetscFunctionReturn(PETSC_SUCCESS);
2913: }
2915: /*@C
2916: PetscDSSetIntegrationParameters - Set the parameters for a particular integration
2918: Not Collective
2920: Input Parameters:
2921: + ds - The `PetscDS` object
2922: . fieldI - The test field for a given point function, or PETSC_DETERMINE
2923: - fieldJ - The basis field for a given point function, or PETSC_DETERMINE
2925: Level: intermediate
2927: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSGetConstants()`, `PetscDSCreate()`
2928: @*/
2929: PetscErrorCode PetscDSSetIntegrationParameters(PetscDS ds, PetscInt fieldI, PetscInt fieldJ)
2930: {
2931: PetscFunctionBegin;
2933: ds->constants[ds->numConstants] = fieldI;
2934: ds->constants[ds->numConstants + 1] = fieldJ;
2935: PetscFunctionReturn(PETSC_SUCCESS);
2936: }
2938: /*@C
2939: PetscDSSetCellParameters - Set the parameters for a particular cell
2941: Not Collective
2943: Input Parameters:
2944: + ds - The `PetscDS` object
2945: - volume - The cell volume
2947: Level: intermediate
2949: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSGetConstants()`, `PetscDSCreate()`
2950: @*/
2951: PetscErrorCode PetscDSSetCellParameters(PetscDS ds, PetscReal volume)
2952: {
2953: PetscFunctionBegin;
2955: ds->constants[ds->numConstants + 2] = volume;
2956: PetscFunctionReturn(PETSC_SUCCESS);
2957: }
2959: /*@
2960: PetscDSGetFieldIndex - Returns the index of the given field
2962: Not Collective
2964: Input Parameters:
2965: + prob - The `PetscDS` object
2966: - disc - The discretization object
2968: Output Parameter:
2969: . f - The field number
2971: Level: beginner
2973: .seealso: `PetscDS`, `PetscGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2974: @*/
2975: PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f)
2976: {
2977: PetscInt g;
2979: PetscFunctionBegin;
2981: PetscAssertPointer(f, 3);
2982: *f = -1;
2983: for (g = 0; g < prob->Nf; ++g) {
2984: if (disc == prob->disc[g]) break;
2985: }
2986: PetscCheck(g != prob->Nf, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS.");
2987: *f = g;
2988: PetscFunctionReturn(PETSC_SUCCESS);
2989: }
2991: /*@
2992: PetscDSGetFieldSize - Returns the size of the given field in the full space basis
2994: Not Collective
2996: Input Parameters:
2997: + prob - The `PetscDS` object
2998: - f - The field number
3000: Output Parameter:
3001: . size - The size
3003: Level: beginner
3005: .seealso: `PetscDS`, `PetscDSGetFieldOffset()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3006: @*/
3007: PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size)
3008: {
3009: PetscFunctionBegin;
3011: PetscAssertPointer(size, 3);
3012: 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);
3013: PetscCall(PetscDSSetUp(prob));
3014: *size = prob->Nb[f];
3015: PetscFunctionReturn(PETSC_SUCCESS);
3016: }
3018: /*@
3019: PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis
3021: Not Collective
3023: Input Parameters:
3024: + prob - The `PetscDS` object
3025: - f - The field number
3027: Output Parameter:
3028: . off - The offset
3030: Level: beginner
3032: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3033: @*/
3034: PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off)
3035: {
3036: PetscInt size, g;
3038: PetscFunctionBegin;
3040: PetscAssertPointer(off, 3);
3041: 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);
3042: *off = 0;
3043: for (g = 0; g < f; ++g) {
3044: PetscCall(PetscDSGetFieldSize(prob, g, &size));
3045: *off += size;
3046: }
3047: PetscFunctionReturn(PETSC_SUCCESS);
3048: }
3050: /*@
3051: PetscDSGetFieldOffsetCohesive - Returns the offset of the given field in the full space basis on a cohesive cell
3053: Not Collective
3055: Input Parameters:
3056: + ds - The `PetscDS` object
3057: - f - The field number
3059: Output Parameter:
3060: . off - The offset
3062: Level: beginner
3064: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3065: @*/
3066: PetscErrorCode PetscDSGetFieldOffsetCohesive(PetscDS ds, PetscInt f, PetscInt *off)
3067: {
3068: PetscInt size, g;
3070: PetscFunctionBegin;
3072: PetscAssertPointer(off, 3);
3073: 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);
3074: *off = 0;
3075: for (g = 0; g < f; ++g) {
3076: PetscBool cohesive;
3078: PetscCall(PetscDSGetCohesive(ds, g, &cohesive));
3079: PetscCall(PetscDSGetFieldSize(ds, g, &size));
3080: *off += cohesive ? size : size * 2;
3081: }
3082: PetscFunctionReturn(PETSC_SUCCESS);
3083: }
3085: /*@
3086: PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point
3088: Not Collective
3090: Input Parameter:
3091: . prob - The `PetscDS` object
3093: Output Parameter:
3094: . dimensions - The number of dimensions
3096: Level: beginner
3098: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3099: @*/
3100: PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[])
3101: {
3102: PetscFunctionBegin;
3104: PetscCall(PetscDSSetUp(prob));
3105: PetscAssertPointer(dimensions, 2);
3106: *dimensions = prob->Nb;
3107: PetscFunctionReturn(PETSC_SUCCESS);
3108: }
3110: /*@
3111: PetscDSGetComponents - Returns the number of components for each field on an evaluation point
3113: Not Collective
3115: Input Parameter:
3116: . prob - The `PetscDS` object
3118: Output Parameter:
3119: . components - The number of components
3121: Level: beginner
3123: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3124: @*/
3125: PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[])
3126: {
3127: PetscFunctionBegin;
3129: PetscCall(PetscDSSetUp(prob));
3130: PetscAssertPointer(components, 2);
3131: *components = prob->Nc;
3132: PetscFunctionReturn(PETSC_SUCCESS);
3133: }
3135: /*@
3136: PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point
3138: Not Collective
3140: Input Parameters:
3141: + prob - The `PetscDS` object
3142: - f - The field number
3144: Output Parameter:
3145: . off - The offset
3147: Level: beginner
3149: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3150: @*/
3151: PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off)
3152: {
3153: PetscFunctionBegin;
3155: PetscAssertPointer(off, 3);
3156: 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);
3157: PetscCall(PetscDSSetUp(prob));
3158: *off = prob->off[f];
3159: PetscFunctionReturn(PETSC_SUCCESS);
3160: }
3162: /*@
3163: PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point
3165: Not Collective
3167: Input Parameter:
3168: . prob - The `PetscDS` object
3170: Output Parameter:
3171: . offsets - The offsets
3173: Level: beginner
3175: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3176: @*/
3177: PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[])
3178: {
3179: PetscFunctionBegin;
3181: PetscAssertPointer(offsets, 2);
3182: PetscCall(PetscDSSetUp(prob));
3183: *offsets = prob->off;
3184: PetscFunctionReturn(PETSC_SUCCESS);
3185: }
3187: /*@
3188: PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point
3190: Not Collective
3192: Input Parameter:
3193: . prob - The `PetscDS` object
3195: Output Parameter:
3196: . offsets - The offsets
3198: Level: beginner
3200: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3201: @*/
3202: PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[])
3203: {
3204: PetscFunctionBegin;
3206: PetscAssertPointer(offsets, 2);
3207: PetscCall(PetscDSSetUp(prob));
3208: *offsets = prob->offDer;
3209: PetscFunctionReturn(PETSC_SUCCESS);
3210: }
3212: /*@
3213: PetscDSGetComponentOffsetsCohesive - Returns the offset of each field on an evaluation point
3215: Not Collective
3217: Input Parameters:
3218: + ds - The `PetscDS` object
3219: - s - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive
3221: Output Parameter:
3222: . offsets - The offsets
3224: Level: beginner
3226: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3227: @*/
3228: PetscErrorCode PetscDSGetComponentOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
3229: {
3230: PetscFunctionBegin;
3232: PetscAssertPointer(offsets, 3);
3233: PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS");
3234: PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s);
3235: PetscCall(PetscDSSetUp(ds));
3236: *offsets = ds->offCohesive[s];
3237: PetscFunctionReturn(PETSC_SUCCESS);
3238: }
3240: /*@
3241: PetscDSGetComponentDerivativeOffsetsCohesive - Returns the offset of each field derivative on an evaluation point
3243: Not Collective
3245: Input Parameters:
3246: + ds - The `PetscDS` object
3247: - s - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive
3249: Output Parameter:
3250: . offsets - The offsets
3252: Level: beginner
3254: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3255: @*/
3256: PetscErrorCode PetscDSGetComponentDerivativeOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
3257: {
3258: PetscFunctionBegin;
3260: PetscAssertPointer(offsets, 3);
3261: PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS");
3262: PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s);
3263: PetscCall(PetscDSSetUp(ds));
3264: *offsets = ds->offDerCohesive[s];
3265: PetscFunctionReturn(PETSC_SUCCESS);
3266: }
3268: /*@C
3269: PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization
3271: Not Collective
3273: Input Parameter:
3274: . prob - The `PetscDS` object
3276: Output Parameter:
3277: . T - The basis function and derivatives tabulation at quadrature points for each field
3279: Level: intermediate
3281: .seealso: `PetscDS`, `PetscTabulation`, `PetscDSCreate()`
3282: @*/
3283: PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscTabulation *T[])
3284: {
3285: PetscFunctionBegin;
3287: PetscAssertPointer(T, 2);
3288: PetscCall(PetscDSSetUp(prob));
3289: *T = prob->T;
3290: PetscFunctionReturn(PETSC_SUCCESS);
3291: }
3293: /*@C
3294: PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces
3296: Not Collective
3298: Input Parameter:
3299: . prob - The `PetscDS` object
3301: Output Parameter:
3302: . Tf - The basis function and derivative tabulation on each local face at quadrature points for each and field
3304: Level: intermediate
3306: .seealso: `PetscTabulation`, `PetscDS`, `PetscDSGetTabulation()`, `PetscDSCreate()`
3307: @*/
3308: PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscTabulation *Tf[])
3309: {
3310: PetscFunctionBegin;
3312: PetscAssertPointer(Tf, 2);
3313: PetscCall(PetscDSSetUp(prob));
3314: *Tf = prob->Tf;
3315: PetscFunctionReturn(PETSC_SUCCESS);
3316: }
3318: PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x)
3319: {
3320: PetscFunctionBegin;
3322: PetscCall(PetscDSSetUp(prob));
3323: if (u) {
3324: PetscAssertPointer(u, 2);
3325: *u = prob->u;
3326: }
3327: if (u_t) {
3328: PetscAssertPointer(u_t, 3);
3329: *u_t = prob->u_t;
3330: }
3331: if (u_x) {
3332: PetscAssertPointer(u_x, 4);
3333: *u_x = prob->u_x;
3334: }
3335: PetscFunctionReturn(PETSC_SUCCESS);
3336: }
3338: PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3)
3339: {
3340: PetscFunctionBegin;
3342: PetscCall(PetscDSSetUp(prob));
3343: if (f0) {
3344: PetscAssertPointer(f0, 2);
3345: *f0 = prob->f0;
3346: }
3347: if (f1) {
3348: PetscAssertPointer(f1, 3);
3349: *f1 = prob->f1;
3350: }
3351: if (g0) {
3352: PetscAssertPointer(g0, 4);
3353: *g0 = prob->g0;
3354: }
3355: if (g1) {
3356: PetscAssertPointer(g1, 5);
3357: *g1 = prob->g1;
3358: }
3359: if (g2) {
3360: PetscAssertPointer(g2, 6);
3361: *g2 = prob->g2;
3362: }
3363: if (g3) {
3364: PetscAssertPointer(g3, 7);
3365: *g3 = prob->g3;
3366: }
3367: PetscFunctionReturn(PETSC_SUCCESS);
3368: }
3370: PetscErrorCode PetscDSGetWorkspace(PetscDS prob, PetscReal **x, PetscScalar **basisReal, PetscScalar **basisDerReal, PetscScalar **testReal, PetscScalar **testDerReal)
3371: {
3372: PetscFunctionBegin;
3374: PetscCall(PetscDSSetUp(prob));
3375: if (x) {
3376: PetscAssertPointer(x, 2);
3377: *x = prob->x;
3378: }
3379: if (basisReal) {
3380: PetscAssertPointer(basisReal, 3);
3381: *basisReal = prob->basisReal;
3382: }
3383: if (basisDerReal) {
3384: PetscAssertPointer(basisDerReal, 4);
3385: *basisDerReal = prob->basisDerReal;
3386: }
3387: if (testReal) {
3388: PetscAssertPointer(testReal, 5);
3389: *testReal = prob->testReal;
3390: }
3391: if (testDerReal) {
3392: PetscAssertPointer(testDerReal, 6);
3393: *testDerReal = prob->testDerReal;
3394: }
3395: PetscFunctionReturn(PETSC_SUCCESS);
3396: }
3398: /*@C
3399: PetscDSAddBoundary - Add a boundary condition to the model.
3401: Collective
3403: Input Parameters:
3404: + ds - The PetscDS object
3405: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3406: . name - The BC name
3407: . label - The label defining constrained points
3408: . Nv - The number of `DMLabel` values for constrained points
3409: . values - An array of label values for constrained points
3410: . field - The field to constrain
3411: . Nc - The number of constrained field components (0 will constrain all fields)
3412: . comps - An array of constrained component numbers
3413: . bcFunc - A pointwise function giving boundary values
3414: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3415: - ctx - An optional user context for bcFunc
3417: Output Parameter:
3418: . bd - The boundary number
3420: Options Database Keys:
3421: + -bc_<boundary name> <num> - Overrides the boundary ids
3422: - -bc_<boundary name>_comp <num> - Overrides the boundary components
3424: Level: developer
3426: Note:
3427: Both `bcFunc` and `bcFunc_t` will depend on the boundary condition type. If the type if `DM_BC_ESSENTIAL`, then the calling sequence is\:
3429: $ void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])
3431: If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value, then the calling sequence is\:
3432: .vb
3433: void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3434: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3435: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3436: PetscReal time, const PetscReal x[], PetscScalar bcval[])
3437: .ve
3438: + dim - the spatial dimension
3439: . Nf - the number of fields
3440: . uOff - the offset into u[] and u_t[] for each field
3441: . uOff_x - the offset into u_x[] for each field
3442: . u - each field evaluated at the current point
3443: . u_t - the time derivative of each field evaluated at the current point
3444: . u_x - the gradient of each field evaluated at the current point
3445: . aOff - the offset into a[] and a_t[] for each auxiliary field
3446: . aOff_x - the offset into a_x[] for each auxiliary field
3447: . a - each auxiliary field evaluated at the current point
3448: . a_t - the time derivative of each auxiliary field evaluated at the current point
3449: . a_x - the gradient of auxiliary each field evaluated at the current point
3450: . t - current time
3451: . x - coordinates of the current point
3452: . numConstants - number of constant parameters
3453: . constants - constant parameters
3454: - bcval - output values at the current point
3456: Notes:
3457: The pointwise functions are used to provide boundary values for essential boundary
3458: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3459: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3460: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3462: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundaryByName()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3463: @*/
3464: PetscErrorCode PetscDSAddBoundary(PetscDS ds, DMBoundaryConditionType type, const char name[], DMLabel label, PetscInt Nv, const PetscInt values[], PetscInt field, PetscInt Nc, const PetscInt comps[], void (*bcFunc)(void), void (*bcFunc_t)(void), void *ctx, PetscInt *bd)
3465: {
3466: DSBoundary head = ds->boundary, b;
3467: PetscInt n = 0;
3468: const char *lname;
3470: PetscFunctionBegin;
3473: PetscAssertPointer(name, 3);
3478: 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);
3479: if (Nc > 0) {
3480: PetscInt *fcomps;
3481: PetscInt c;
3483: PetscCall(PetscDSGetComponents(ds, &fcomps));
3484: 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);
3485: for (c = 0; c < Nc; ++c) {
3486: 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);
3487: }
3488: }
3489: PetscCall(PetscNew(&b));
3490: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3491: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf));
3492: PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf));
3493: PetscCall(PetscMalloc1(Nv, &b->values));
3494: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3495: PetscCall(PetscMalloc1(Nc, &b->comps));
3496: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3497: PetscCall(PetscObjectGetName((PetscObject)label, &lname));
3498: PetscCall(PetscStrallocpy(lname, (char **)&b->lname));
3499: b->type = type;
3500: b->label = label;
3501: b->Nv = Nv;
3502: b->field = field;
3503: b->Nc = Nc;
3504: b->func = bcFunc;
3505: b->func_t = bcFunc_t;
3506: b->ctx = ctx;
3507: b->next = NULL;
3508: /* Append to linked list so that we can preserve the order */
3509: if (!head) ds->boundary = b;
3510: while (head) {
3511: if (!head->next) {
3512: head->next = b;
3513: head = b;
3514: }
3515: head = head->next;
3516: ++n;
3517: }
3518: if (bd) {
3519: PetscAssertPointer(bd, 13);
3520: *bd = n;
3521: }
3522: PetscFunctionReturn(PETSC_SUCCESS);
3523: }
3525: // PetscClangLinter pragma ignore: -fdoc-section-header-unknown
3526: /*@C
3527: PetscDSAddBoundaryByName - Add a boundary condition to the model.
3529: Collective
3531: Input Parameters:
3532: + ds - The `PetscDS` object
3533: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3534: . name - The BC name
3535: . lname - The naem of the label defining constrained points
3536: . Nv - The number of `DMLabel` values for constrained points
3537: . values - An array of label values for constrained points
3538: . field - The field to constrain
3539: . Nc - The number of constrained field components (0 will constrain all fields)
3540: . comps - An array of constrained component numbers
3541: . bcFunc - A pointwise function giving boundary values
3542: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3543: - ctx - An optional user context for bcFunc
3545: Output Parameter:
3546: . bd - The boundary number
3548: Options Database Keys:
3549: + -bc_<boundary name> <num> - Overrides the boundary ids
3550: - -bc_<boundary name>_comp <num> - Overrides the boundary components
3552: Calling Sequence of `bcFunc` and `bcFunc_t`:
3553: If the type is `DM_BC_ESSENTIAL`
3554: .vb
3555: void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])
3556: .ve
3557: If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value,
3558: .vb
3559: void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3560: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3561: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3562: PetscReal time, const PetscReal x[], PetscScalar bcval[])
3563: .ve
3564: + dim - the spatial dimension
3565: . Nf - the number of fields
3566: . uOff - the offset into u[] and u_t[] for each field
3567: . uOff_x - the offset into u_x[] for each field
3568: . u - each field evaluated at the current point
3569: . u_t - the time derivative of each field evaluated at the current point
3570: . u_x - the gradient of each field evaluated at the current point
3571: . aOff - the offset into a[] and a_t[] for each auxiliary field
3572: . aOff_x - the offset into a_x[] for each auxiliary field
3573: . a - each auxiliary field evaluated at the current point
3574: . a_t - the time derivative of each auxiliary field evaluated at the current point
3575: . a_x - the gradient of auxiliary each field evaluated at the current point
3576: . t - current time
3577: . x - coordinates of the current point
3578: . numConstants - number of constant parameters
3579: . constants - constant parameters
3580: - bcval - output values at the current point
3582: Level: developer
3584: Notes:
3585: The pointwise functions are used to provide boundary values for essential boundary
3586: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3587: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3588: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3590: This function should only be used with `DMFOREST` currently, since labels cannot be defined before the underlying `DMPLEX` is built.
3592: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3593: @*/
3594: PetscErrorCode PetscDSAddBoundaryByName(PetscDS ds, DMBoundaryConditionType type, const char name[], const char lname[], PetscInt Nv, const PetscInt values[], PetscInt field, PetscInt Nc, const PetscInt comps[], void (*bcFunc)(void), void (*bcFunc_t)(void), void *ctx, PetscInt *bd)
3595: {
3596: DSBoundary head = ds->boundary, b;
3597: PetscInt n = 0;
3599: PetscFunctionBegin;
3602: PetscAssertPointer(name, 3);
3603: PetscAssertPointer(lname, 4);
3607: PetscCall(PetscNew(&b));
3608: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3609: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf));
3610: PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf));
3611: PetscCall(PetscMalloc1(Nv, &b->values));
3612: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3613: PetscCall(PetscMalloc1(Nc, &b->comps));
3614: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3615: PetscCall(PetscStrallocpy(lname, (char **)&b->lname));
3616: b->type = type;
3617: b->label = NULL;
3618: b->Nv = Nv;
3619: b->field = field;
3620: b->Nc = Nc;
3621: b->func = bcFunc;
3622: b->func_t = bcFunc_t;
3623: b->ctx = ctx;
3624: b->next = NULL;
3625: /* Append to linked list so that we can preserve the order */
3626: if (!head) ds->boundary = b;
3627: while (head) {
3628: if (!head->next) {
3629: head->next = b;
3630: head = b;
3631: }
3632: head = head->next;
3633: ++n;
3634: }
3635: if (bd) {
3636: PetscAssertPointer(bd, 13);
3637: *bd = n;
3638: }
3639: PetscFunctionReturn(PETSC_SUCCESS);
3640: }
3642: /*@C
3643: PetscDSUpdateBoundary - Change a boundary condition for the model.
3645: Input Parameters:
3646: + ds - The `PetscDS` object
3647: . bd - The boundary condition number
3648: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3649: . name - The BC name
3650: . label - The label defining constrained points
3651: . Nv - The number of `DMLabel` ids for constrained points
3652: . values - An array of ids for constrained points
3653: . field - The field to constrain
3654: . Nc - The number of constrained field components
3655: . comps - An array of constrained component numbers
3656: . bcFunc - A pointwise function giving boundary values
3657: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3658: - ctx - An optional user context for bcFunc
3660: Level: developer
3662: Notes:
3663: The pointwise functions are used to provide boundary values for essential boundary
3664: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3665: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3666: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3668: The boundary condition number is the order in which it was registered. The user can get the number of boundary conditions from `PetscDSGetNumBoundary()`.
3669: See `PetscDSAddBoundary()` for a description of the calling sequences for the callbacks.
3671: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSGetNumBoundary()`, `DMLabel`
3672: @*/
3673: PetscErrorCode PetscDSUpdateBoundary(PetscDS ds, PetscInt bd, DMBoundaryConditionType type, const char name[], DMLabel label, PetscInt Nv, const PetscInt values[], PetscInt field, PetscInt Nc, const PetscInt comps[], void (*bcFunc)(void), void (*bcFunc_t)(void), void *ctx)
3674: {
3675: DSBoundary b = ds->boundary;
3676: PetscInt n = 0;
3678: PetscFunctionBegin;
3680: while (b) {
3681: if (n == bd) break;
3682: b = b->next;
3683: ++n;
3684: }
3685: PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n);
3686: if (name) {
3687: PetscCall(PetscFree(b->name));
3688: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3689: }
3690: b->type = type;
3691: if (label) {
3692: const char *name;
3694: b->label = label;
3695: PetscCall(PetscFree(b->lname));
3696: PetscCall(PetscObjectGetName((PetscObject)label, &name));
3697: PetscCall(PetscStrallocpy(name, (char **)&b->lname));
3698: }
3699: if (Nv >= 0) {
3700: b->Nv = Nv;
3701: PetscCall(PetscFree(b->values));
3702: PetscCall(PetscMalloc1(Nv, &b->values));
3703: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3704: }
3705: if (field >= 0) b->field = field;
3706: if (Nc >= 0) {
3707: b->Nc = Nc;
3708: PetscCall(PetscFree(b->comps));
3709: PetscCall(PetscMalloc1(Nc, &b->comps));
3710: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3711: }
3712: if (bcFunc) b->func = bcFunc;
3713: if (bcFunc_t) b->func_t = bcFunc_t;
3714: if (ctx) b->ctx = ctx;
3715: PetscFunctionReturn(PETSC_SUCCESS);
3716: }
3718: /*@
3719: PetscDSGetNumBoundary - Get the number of registered BC
3721: Input Parameter:
3722: . ds - The `PetscDS` object
3724: Output Parameter:
3725: . numBd - The number of BC
3727: Level: intermediate
3729: .seealso: `PetscDS`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`
3730: @*/
3731: PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd)
3732: {
3733: DSBoundary b = ds->boundary;
3735: PetscFunctionBegin;
3737: PetscAssertPointer(numBd, 2);
3738: *numBd = 0;
3739: while (b) {
3740: ++(*numBd);
3741: b = b->next;
3742: }
3743: PetscFunctionReturn(PETSC_SUCCESS);
3744: }
3746: /*@C
3747: PetscDSGetBoundary - Gets a boundary condition to the model
3749: Input Parameters:
3750: + ds - The `PetscDS` object
3751: - bd - The BC number
3753: Output Parameters:
3754: + wf - The `PetscWeakForm` holding the pointwise functions
3755: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3756: . name - The BC name
3757: . label - The label defining constrained points
3758: . Nv - The number of `DMLabel` ids for constrained points
3759: . values - An array of ids for constrained points
3760: . field - The field to constrain
3761: . Nc - The number of constrained field components
3762: . comps - An array of constrained component numbers
3763: . func - A pointwise function giving boundary values
3764: . func_t - A pointwise function giving the time derivative of the boundary values
3765: - ctx - An optional user context for bcFunc
3767: Options Database Keys:
3768: + -bc_<boundary name> <num> - Overrides the boundary ids
3769: - -bc_<boundary name>_comp <num> - Overrides the boundary components
3771: Level: developer
3773: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `DMLabel`
3774: @*/
3775: PetscErrorCode PetscDSGetBoundary(PetscDS ds, PetscInt bd, PetscWeakForm *wf, DMBoundaryConditionType *type, const char *name[], DMLabel *label, PetscInt *Nv, const PetscInt *values[], PetscInt *field, PetscInt *Nc, const PetscInt *comps[], void (**func)(void), void (**func_t)(void), void **ctx)
3776: {
3777: DSBoundary b = ds->boundary;
3778: PetscInt n = 0;
3780: PetscFunctionBegin;
3782: while (b) {
3783: if (n == bd) break;
3784: b = b->next;
3785: ++n;
3786: }
3787: PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n);
3788: if (wf) {
3789: PetscAssertPointer(wf, 3);
3790: *wf = b->wf;
3791: }
3792: if (type) {
3793: PetscAssertPointer(type, 4);
3794: *type = b->type;
3795: }
3796: if (name) {
3797: PetscAssertPointer(name, 5);
3798: *name = b->name;
3799: }
3800: if (label) {
3801: PetscAssertPointer(label, 6);
3802: *label = b->label;
3803: }
3804: if (Nv) {
3805: PetscAssertPointer(Nv, 7);
3806: *Nv = b->Nv;
3807: }
3808: if (values) {
3809: PetscAssertPointer(values, 8);
3810: *values = b->values;
3811: }
3812: if (field) {
3813: PetscAssertPointer(field, 9);
3814: *field = b->field;
3815: }
3816: if (Nc) {
3817: PetscAssertPointer(Nc, 10);
3818: *Nc = b->Nc;
3819: }
3820: if (comps) {
3821: PetscAssertPointer(comps, 11);
3822: *comps = b->comps;
3823: }
3824: if (func) {
3825: PetscAssertPointer(func, 12);
3826: *func = b->func;
3827: }
3828: if (func_t) {
3829: PetscAssertPointer(func_t, 13);
3830: *func_t = b->func_t;
3831: }
3832: if (ctx) {
3833: PetscAssertPointer(ctx, 14);
3834: *ctx = b->ctx;
3835: }
3836: PetscFunctionReturn(PETSC_SUCCESS);
3837: }
3839: /*@
3840: PetscDSUpdateBoundaryLabels - Update `DMLabel` in each boundary condition using the label name and the input `DM`
3842: Not Collective
3844: Input Parameters:
3845: + ds - The source `PetscDS` object
3846: - dm - The `DM` holding labels
3848: Level: intermediate
3850: .seealso: `PetscDS`, `DMBoundary`, `DM`, `PetscDSCopyBoundary()`, `PetscDSCreate()`, `DMGetLabel()`
3851: @*/
3852: PetscErrorCode PetscDSUpdateBoundaryLabels(PetscDS ds, DM dm)
3853: {
3854: DSBoundary b;
3856: PetscFunctionBegin;
3859: for (b = ds->boundary; b; b = b->next) {
3860: if (b->lname) PetscCall(DMGetLabel(dm, b->lname, &b->label));
3861: }
3862: PetscFunctionReturn(PETSC_SUCCESS);
3863: }
3865: static PetscErrorCode DSBoundaryDuplicate_Internal(DSBoundary b, DSBoundary *bNew)
3866: {
3867: PetscFunctionBegin;
3868: PetscCall(PetscNew(bNew));
3869: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &(*bNew)->wf));
3870: PetscCall(PetscWeakFormCopy(b->wf, (*bNew)->wf));
3871: PetscCall(PetscStrallocpy(b->name, (char **)&((*bNew)->name)));
3872: PetscCall(PetscStrallocpy(b->lname, (char **)&((*bNew)->lname)));
3873: (*bNew)->type = b->type;
3874: (*bNew)->label = b->label;
3875: (*bNew)->Nv = b->Nv;
3876: PetscCall(PetscMalloc1(b->Nv, &(*bNew)->values));
3877: PetscCall(PetscArraycpy((*bNew)->values, b->values, b->Nv));
3878: (*bNew)->field = b->field;
3879: (*bNew)->Nc = b->Nc;
3880: PetscCall(PetscMalloc1(b->Nc, &(*bNew)->comps));
3881: PetscCall(PetscArraycpy((*bNew)->comps, b->comps, b->Nc));
3882: (*bNew)->func = b->func;
3883: (*bNew)->func_t = b->func_t;
3884: (*bNew)->ctx = b->ctx;
3885: PetscFunctionReturn(PETSC_SUCCESS);
3886: }
3888: /*@
3889: PetscDSCopyBoundary - Copy all boundary condition objects to the new problem
3891: Not Collective
3893: Input Parameters:
3894: + ds - The source `PetscDS` object
3895: . numFields - The number of selected fields, or `PETSC_DEFAULT` for all fields
3896: - fields - The selected fields, or NULL for all fields
3898: Output Parameter:
3899: . newds - The target `PetscDS`, now with a copy of the boundary conditions
3901: Level: intermediate
3903: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3904: @*/
3905: PetscErrorCode PetscDSCopyBoundary(PetscDS ds, PetscInt numFields, const PetscInt fields[], PetscDS newds)
3906: {
3907: DSBoundary b, *lastnext;
3909: PetscFunctionBegin;
3912: if (ds == newds) PetscFunctionReturn(PETSC_SUCCESS);
3913: PetscCall(PetscDSDestroyBoundary(newds));
3914: lastnext = &newds->boundary;
3915: for (b = ds->boundary; b; b = b->next) {
3916: DSBoundary bNew;
3917: PetscInt fieldNew = -1;
3919: if (numFields > 0 && fields) {
3920: PetscInt f;
3922: for (f = 0; f < numFields; ++f)
3923: if (b->field == fields[f]) break;
3924: if (f == numFields) continue;
3925: fieldNew = f;
3926: }
3927: PetscCall(DSBoundaryDuplicate_Internal(b, &bNew));
3928: bNew->field = fieldNew < 0 ? b->field : fieldNew;
3929: *lastnext = bNew;
3930: lastnext = &bNew->next;
3931: }
3932: PetscFunctionReturn(PETSC_SUCCESS);
3933: }
3935: /*@
3936: PetscDSDestroyBoundary - Remove all `DMBoundary` objects from the `PetscDS`
3938: Not Collective
3940: Input Parameter:
3941: . ds - The `PetscDS` object
3943: Level: intermediate
3945: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`
3946: @*/
3947: PetscErrorCode PetscDSDestroyBoundary(PetscDS ds)
3948: {
3949: DSBoundary next = ds->boundary;
3951: PetscFunctionBegin;
3952: while (next) {
3953: DSBoundary b = next;
3955: next = b->next;
3956: PetscCall(PetscWeakFormDestroy(&b->wf));
3957: PetscCall(PetscFree(b->name));
3958: PetscCall(PetscFree(b->lname));
3959: PetscCall(PetscFree(b->values));
3960: PetscCall(PetscFree(b->comps));
3961: PetscCall(PetscFree(b));
3962: }
3963: PetscFunctionReturn(PETSC_SUCCESS);
3964: }
3966: /*@
3967: PetscDSSelectDiscretizations - Copy discretizations to the new problem with different field layout
3969: Not Collective
3971: Input Parameters:
3972: + prob - The `PetscDS` object
3973: . numFields - Number of new fields
3974: . fields - Old field number for each new field
3975: . minDegree - Minimum degree for a discretization, or `PETSC_DETERMINE` for no limit
3976: - maxDegree - Maximum degree for a discretization, or `PETSC_DETERMINE` for no limit
3978: Output Parameter:
3979: . newprob - The `PetscDS` copy
3981: Level: intermediate
3983: .seealso: `PetscDS`, `PetscDSSelectEquations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3984: @*/
3985: PetscErrorCode PetscDSSelectDiscretizations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscInt minDegree, PetscInt maxDegree, PetscDS newprob)
3986: {
3987: PetscInt Nf, Nfn, fn;
3989: PetscFunctionBegin;
3991: if (fields) PetscAssertPointer(fields, 3);
3993: PetscCall(PetscDSGetNumFields(prob, &Nf));
3994: PetscCall(PetscDSGetNumFields(newprob, &Nfn));
3995: numFields = numFields < 0 ? Nf : numFields;
3996: for (fn = 0; fn < numFields; ++fn) {
3997: const PetscInt f = fields ? fields[fn] : fn;
3998: PetscObject disc;
3999: PetscClassId id;
4001: if (f >= Nf) continue;
4002: PetscCall(PetscDSGetDiscretization(prob, f, &disc));
4003: PetscCallContinue(PetscObjectGetClassId(disc, &id));
4004: if (id == PETSCFE_CLASSID) {
4005: PetscFE fe;
4007: PetscCall(PetscFELimitDegree((PetscFE)disc, minDegree, maxDegree, &fe));
4008: PetscCall(PetscDSSetDiscretization(newprob, fn, (PetscObject)fe));
4009: PetscCall(PetscFEDestroy(&fe));
4010: } else {
4011: PetscCall(PetscDSSetDiscretization(newprob, fn, disc));
4012: }
4013: }
4014: PetscFunctionReturn(PETSC_SUCCESS);
4015: }
4017: /*@
4018: PetscDSSelectEquations - Copy pointwise function pointers to the new problem with different field layout
4020: Not Collective
4022: Input Parameters:
4023: + prob - The `PetscDS` object
4024: . numFields - Number of new fields
4025: - fields - Old field number for each new field
4027: Output Parameter:
4028: . newprob - The `PetscDS` copy
4030: Level: intermediate
4032: .seealso: `PetscDS`, `PetscDSSelectDiscretizations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
4033: @*/
4034: PetscErrorCode PetscDSSelectEquations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
4035: {
4036: PetscInt Nf, Nfn, fn, gn;
4038: PetscFunctionBegin;
4040: if (fields) PetscAssertPointer(fields, 3);
4042: PetscCall(PetscDSGetNumFields(prob, &Nf));
4043: PetscCall(PetscDSGetNumFields(newprob, &Nfn));
4044: PetscCheck(numFields <= Nfn, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_SIZ, "Number of fields %" PetscInt_FMT " to transfer must not be greater then the total number of fields %" PetscInt_FMT, numFields, Nfn);
4045: for (fn = 0; fn < numFields; ++fn) {
4046: const PetscInt f = fields ? fields[fn] : fn;
4047: PetscPointFunc obj;
4048: PetscPointFunc f0, f1;
4049: PetscBdPointFunc f0Bd, f1Bd;
4050: PetscRiemannFunc r;
4052: if (f >= Nf) continue;
4053: PetscCall(PetscDSGetObjective(prob, f, &obj));
4054: PetscCall(PetscDSGetResidual(prob, f, &f0, &f1));
4055: PetscCall(PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd));
4056: PetscCall(PetscDSGetRiemannSolver(prob, f, &r));
4057: PetscCall(PetscDSSetObjective(newprob, fn, obj));
4058: PetscCall(PetscDSSetResidual(newprob, fn, f0, f1));
4059: PetscCall(PetscDSSetBdResidual(newprob, fn, f0Bd, f1Bd));
4060: PetscCall(PetscDSSetRiemannSolver(newprob, fn, r));
4061: for (gn = 0; gn < numFields; ++gn) {
4062: const PetscInt g = fields ? fields[gn] : gn;
4063: PetscPointJac g0, g1, g2, g3;
4064: PetscPointJac g0p, g1p, g2p, g3p;
4065: PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd;
4067: if (g >= Nf) continue;
4068: PetscCall(PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3));
4069: PetscCall(PetscDSGetJacobianPreconditioner(prob, f, g, &g0p, &g1p, &g2p, &g3p));
4070: PetscCall(PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd));
4071: PetscCall(PetscDSSetJacobian(newprob, fn, gn, g0, g1, g2, g3));
4072: PetscCall(PetscDSSetJacobianPreconditioner(newprob, fn, gn, g0p, g1p, g2p, g3p));
4073: PetscCall(PetscDSSetBdJacobian(newprob, fn, gn, g0Bd, g1Bd, g2Bd, g3Bd));
4074: }
4075: }
4076: PetscFunctionReturn(PETSC_SUCCESS);
4077: }
4079: /*@
4080: PetscDSCopyEquations - Copy all pointwise function pointers to another `PetscDS`
4082: Not Collective
4084: Input Parameter:
4085: . prob - The `PetscDS` object
4087: Output Parameter:
4088: . newprob - The `PetscDS` copy
4090: Level: intermediate
4092: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
4093: @*/
4094: PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob)
4095: {
4096: PetscWeakForm wf, newwf;
4097: PetscInt Nf, Ng;
4099: PetscFunctionBegin;
4102: PetscCall(PetscDSGetNumFields(prob, &Nf));
4103: PetscCall(PetscDSGetNumFields(newprob, &Ng));
4104: PetscCheck(Nf == Ng, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %" PetscInt_FMT " != %" PetscInt_FMT, Nf, Ng);
4105: PetscCall(PetscDSGetWeakForm(prob, &wf));
4106: PetscCall(PetscDSGetWeakForm(newprob, &newwf));
4107: PetscCall(PetscWeakFormCopy(wf, newwf));
4108: PetscFunctionReturn(PETSC_SUCCESS);
4109: }
4111: /*@
4112: PetscDSCopyConstants - Copy all constants to another `PetscDS`
4114: Not Collective
4116: Input Parameter:
4117: . prob - The `PetscDS` object
4119: Output Parameter:
4120: . newprob - The `PetscDS` copy
4122: Level: intermediate
4124: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
4125: @*/
4126: PetscErrorCode PetscDSCopyConstants(PetscDS prob, PetscDS newprob)
4127: {
4128: PetscInt Nc;
4129: const PetscScalar *constants;
4131: PetscFunctionBegin;
4134: PetscCall(PetscDSGetConstants(prob, &Nc, &constants));
4135: PetscCall(PetscDSSetConstants(newprob, Nc, (PetscScalar *)constants));
4136: PetscFunctionReturn(PETSC_SUCCESS);
4137: }
4139: /*@
4140: PetscDSCopyExactSolutions - Copy all exact solutions to another `PetscDS`
4142: Not Collective
4144: Input Parameter:
4145: . ds - The `PetscDS` object
4147: Output Parameter:
4148: . newds - The `PetscDS` copy
4150: Level: intermediate
4152: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
4153: @*/
4154: PetscErrorCode PetscDSCopyExactSolutions(PetscDS ds, PetscDS newds)
4155: {
4156: PetscSimplePointFn *sol;
4157: void *ctx;
4158: PetscInt Nf, f;
4160: PetscFunctionBegin;
4163: PetscCall(PetscDSGetNumFields(ds, &Nf));
4164: for (f = 0; f < Nf; ++f) {
4165: PetscCall(PetscDSGetExactSolution(ds, f, &sol, &ctx));
4166: PetscCall(PetscDSSetExactSolution(newds, f, sol, ctx));
4167: PetscCall(PetscDSGetExactSolutionTimeDerivative(ds, f, &sol, &ctx));
4168: PetscCall(PetscDSSetExactSolutionTimeDerivative(newds, f, sol, ctx));
4169: }
4170: PetscFunctionReturn(PETSC_SUCCESS);
4171: }
4173: PetscErrorCode PetscDSCopy(PetscDS ds, PetscInt minDegree, PetscInt maxDegree, DM dmNew, PetscDS dsNew)
4174: {
4175: DSBoundary b;
4176: PetscInt cdim, Nf, f, d;
4177: PetscBool isCohesive;
4178: void *ctx;
4180: PetscFunctionBegin;
4181: PetscCall(PetscDSCopyConstants(ds, dsNew));
4182: PetscCall(PetscDSCopyExactSolutions(ds, dsNew));
4183: PetscCall(PetscDSSelectDiscretizations(ds, PETSC_DETERMINE, NULL, minDegree, maxDegree, dsNew));
4184: PetscCall(PetscDSCopyEquations(ds, dsNew));
4185: PetscCall(PetscDSGetNumFields(ds, &Nf));
4186: for (f = 0; f < Nf; ++f) {
4187: PetscCall(PetscDSGetContext(ds, f, &ctx));
4188: PetscCall(PetscDSSetContext(dsNew, f, ctx));
4189: PetscCall(PetscDSGetCohesive(ds, f, &isCohesive));
4190: PetscCall(PetscDSSetCohesive(dsNew, f, isCohesive));
4191: PetscCall(PetscDSGetJetDegree(ds, f, &d));
4192: PetscCall(PetscDSSetJetDegree(dsNew, f, d));
4193: }
4194: if (Nf) {
4195: PetscCall(PetscDSGetCoordinateDimension(ds, &cdim));
4196: PetscCall(PetscDSSetCoordinateDimension(dsNew, cdim));
4197: }
4198: PetscCall(PetscDSCopyBoundary(ds, PETSC_DETERMINE, NULL, dsNew));
4199: for (b = dsNew->boundary; b; b = b->next) {
4200: PetscCall(DMGetLabel(dmNew, b->lname, &b->label));
4201: /* Do not check if label exists here, since p4est calls this for the reference tree which does not have the labels */
4202: //PetscCheck(b->label,PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Label %s missing in new DM", name);
4203: }
4204: PetscFunctionReturn(PETSC_SUCCESS);
4205: }
4207: PetscErrorCode PetscDSGetHeightSubspace(PetscDS prob, PetscInt height, PetscDS *subprob)
4208: {
4209: PetscInt dim, Nf, f;
4211: PetscFunctionBegin;
4213: PetscAssertPointer(subprob, 3);
4214: if (height == 0) {
4215: *subprob = prob;
4216: PetscFunctionReturn(PETSC_SUCCESS);
4217: }
4218: PetscCall(PetscDSGetNumFields(prob, &Nf));
4219: PetscCall(PetscDSGetSpatialDimension(prob, &dim));
4220: PetscCheck(height <= dim, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_OUTOFRANGE, "DS can only handle height in [0, %" PetscInt_FMT "], not %" PetscInt_FMT, dim, height);
4221: if (!prob->subprobs) PetscCall(PetscCalloc1(dim, &prob->subprobs));
4222: if (!prob->subprobs[height - 1]) {
4223: PetscInt cdim;
4225: PetscCall(PetscDSCreate(PetscObjectComm((PetscObject)prob), &prob->subprobs[height - 1]));
4226: PetscCall(PetscDSGetCoordinateDimension(prob, &cdim));
4227: PetscCall(PetscDSSetCoordinateDimension(prob->subprobs[height - 1], cdim));
4228: for (f = 0; f < Nf; ++f) {
4229: PetscFE subfe;
4230: PetscObject obj;
4231: PetscClassId id;
4233: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
4234: PetscCall(PetscObjectGetClassId(obj, &id));
4235: if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetHeightSubspace((PetscFE)obj, height, &subfe));
4236: else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unsupported discretization type for field %" PetscInt_FMT, f);
4237: PetscCall(PetscDSSetDiscretization(prob->subprobs[height - 1], f, (PetscObject)subfe));
4238: }
4239: }
4240: *subprob = prob->subprobs[height - 1];
4241: PetscFunctionReturn(PETSC_SUCCESS);
4242: }
4244: PetscErrorCode PetscDSPermuteQuadPoint(PetscDS ds, PetscInt ornt, PetscInt field, PetscInt q, PetscInt *qperm)
4245: {
4246: IS permIS;
4247: PetscQuadrature quad;
4248: DMPolytopeType ct;
4249: const PetscInt *perm;
4250: PetscInt Na, Nq;
4252: PetscFunctionBeginHot;
4253: PetscCall(PetscFEGetQuadrature((PetscFE)ds->disc[field], &quad));
4254: PetscCall(PetscQuadratureGetData(quad, NULL, NULL, &Nq, NULL, NULL));
4255: PetscCall(PetscQuadratureGetCellType(quad, &ct));
4256: PetscCheck(q >= 0 && q < Nq, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Quadrature point %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", q, Nq);
4257: Na = DMPolytopeTypeGetNumArrangements(ct) / 2;
4258: 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);
4259: if (!ds->quadPerm[(PetscInt)ct]) PetscCall(PetscQuadratureComputePermutations(quad, NULL, &ds->quadPerm[(PetscInt)ct]));
4260: permIS = ds->quadPerm[(PetscInt)ct][ornt + Na];
4261: PetscCall(ISGetIndices(permIS, &perm));
4262: *qperm = perm[q];
4263: PetscCall(ISRestoreIndices(permIS, &perm));
4264: PetscFunctionReturn(PETSC_SUCCESS);
4265: }
4267: PetscErrorCode PetscDSGetDiscType_Internal(PetscDS ds, PetscInt f, PetscDiscType *disctype)
4268: {
4269: PetscObject obj;
4270: PetscClassId id;
4271: PetscInt Nf;
4273: PetscFunctionBegin;
4275: PetscAssertPointer(disctype, 3);
4276: *disctype = PETSC_DISC_NONE;
4277: PetscCall(PetscDSGetNumFields(ds, &Nf));
4278: PetscCheck(f < Nf, PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_SIZ, "Field %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, Nf);
4279: PetscCall(PetscDSGetDiscretization(ds, f, &obj));
4280: if (obj) {
4281: PetscCall(PetscObjectGetClassId(obj, &id));
4282: if (id == PETSCFE_CLASSID) *disctype = PETSC_DISC_FE;
4283: else *disctype = PETSC_DISC_FV;
4284: }
4285: PetscFunctionReturn(PETSC_SUCCESS);
4286: }
4288: static PetscErrorCode PetscDSDestroy_Basic(PetscDS ds)
4289: {
4290: PetscFunctionBegin;
4291: PetscCall(PetscFree(ds->data));
4292: PetscFunctionReturn(PETSC_SUCCESS);
4293: }
4295: static PetscErrorCode PetscDSInitialize_Basic(PetscDS ds)
4296: {
4297: PetscFunctionBegin;
4298: ds->ops->setfromoptions = NULL;
4299: ds->ops->setup = NULL;
4300: ds->ops->view = NULL;
4301: ds->ops->destroy = PetscDSDestroy_Basic;
4302: PetscFunctionReturn(PETSC_SUCCESS);
4303: }
4305: /*MC
4306: PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions
4308: Level: intermediate
4310: .seealso: `PetscDSType`, `PetscDSCreate()`, `PetscDSSetType()`
4311: M*/
4313: PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS ds)
4314: {
4315: PetscDS_Basic *b;
4317: PetscFunctionBegin;
4319: PetscCall(PetscNew(&b));
4320: ds->data = b;
4322: PetscCall(PetscDSInitialize_Basic(ds));
4323: PetscFunctionReturn(PETSC_SUCCESS);
4324: }