Actual source code: dtds.c
1: #include <petsc/private/petscdsimpl.h>
3: PetscClassId PETSCDS_CLASSID = 0;
5: PetscFunctionList PetscDSList = NULL;
6: PetscBool PetscDSRegisterAllCalled = PETSC_FALSE;
8: /* A PetscDS (Discrete System) encodes a set of equations posed in a discrete space, which represents a set of
9: nonlinear continuum equations. The equations can have multiple fields, each field having a different
10: discretization. In addition, different pieces of the domain can have different field combinations and equations.
12: The DS provides the user a description of the approximation space on any given cell. It also gives pointwise
13: functions representing the equations.
15: Each field is associated with a label, marking the cells on which it is supported. Note that a field can be
16: supported on the closure of a cell not in the label due to overlap of the boundary of neighboring cells. The DM
17: then creates a DS for each set of cells with identical approximation spaces. When assembling, the user asks for
18: the space associated with a given cell. DMPlex uses the labels associated with each DS in the default integration loop.
19: */
21: /*@C
22: PetscDSRegister - Adds a new `PetscDS` implementation
24: Not Collective; No Fortran Support
26: Input Parameters:
27: + sname - The name of a new user-defined creation routine
28: - function - The creation routine itself
30: Example Usage:
31: .vb
32: PetscDSRegister("my_ds", MyPetscDSCreate);
33: .ve
35: Then, your PetscDS type can be chosen with the procedural interface via
36: .vb
37: PetscDSCreate(MPI_Comm, PetscDS *);
38: PetscDSSetType(PetscDS, "my_ds");
39: .ve
40: or at runtime via the option
41: .vb
42: -petscds_type my_ds
43: .ve
45: Level: advanced
47: Note:
48: `PetscDSRegister()` may be called multiple times to add several user-defined `PetscDSs`
50: .seealso: `PetscDSType`, `PetscDS`, `PetscDSRegisterAll()`, `PetscDSRegisterDestroy()`
51: @*/
52: PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS))
53: {
54: PetscFunctionBegin;
55: PetscCall(PetscFunctionListAdd(&PetscDSList, sname, function));
56: PetscFunctionReturn(PETSC_SUCCESS);
57: }
59: /*@C
60: PetscDSSetType - Builds a particular `PetscDS`
62: Collective; No Fortran Support
64: Input Parameters:
65: + prob - The `PetscDS` object
66: - name - The `PetscDSType`
68: Options Database Key:
69: . -petscds_type <type> - Sets the PetscDS type; use -help for a list of available types
71: Level: intermediate
73: .seealso: `PetscDSType`, `PetscDS`, `PetscDSGetType()`, `PetscDSCreate()`
74: @*/
75: PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name)
76: {
77: PetscErrorCode (*r)(PetscDS);
78: PetscBool match;
80: PetscFunctionBegin;
82: PetscCall(PetscObjectTypeCompare((PetscObject)prob, name, &match));
83: if (match) PetscFunctionReturn(PETSC_SUCCESS);
85: PetscCall(PetscDSRegisterAll());
86: PetscCall(PetscFunctionListFind(PetscDSList, name, &r));
87: PetscCheck(r, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscDS type: %s", name);
89: PetscTryTypeMethod(prob, destroy);
90: prob->ops->destroy = NULL;
92: PetscCall((*r)(prob));
93: PetscCall(PetscObjectChangeTypeName((PetscObject)prob, name));
94: PetscFunctionReturn(PETSC_SUCCESS);
95: }
97: /*@C
98: PetscDSGetType - Gets the `PetscDSType` name (as a string) from the `PetscDS`
100: Not Collective; No Fortran Support
102: Input Parameter:
103: . prob - The `PetscDS`
105: Output Parameter:
106: . name - The `PetscDSType` name
108: Level: intermediate
110: .seealso: `PetscDSType`, `PetscDS`, `PetscDSSetType()`, `PetscDSCreate()`
111: @*/
112: PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name)
113: {
114: PetscFunctionBegin;
116: PetscAssertPointer(name, 2);
117: PetscCall(PetscDSRegisterAll());
118: *name = ((PetscObject)prob)->type_name;
119: PetscFunctionReturn(PETSC_SUCCESS);
120: }
122: static PetscErrorCode PetscDSView_Ascii(PetscDS ds, PetscViewer viewer)
123: {
124: PetscViewerFormat format;
125: const PetscScalar *constants;
126: PetscInt Nf, numConstants, f;
128: PetscFunctionBegin;
129: PetscCall(PetscDSGetNumFields(ds, &Nf));
130: PetscCall(PetscViewerGetFormat(viewer, &format));
131: PetscCall(PetscViewerASCIIPrintf(viewer, "Discrete System with %" PetscInt_FMT " fields\n", Nf));
132: PetscCall(PetscViewerASCIIPushTab(viewer));
133: PetscCall(PetscViewerASCIIPrintf(viewer, " cell total dim %" PetscInt_FMT " total comp %" PetscInt_FMT "\n", ds->totDim, ds->totComp));
134: if (ds->isCohesive) PetscCall(PetscViewerASCIIPrintf(viewer, " cohesive cell\n"));
135: for (f = 0; f < Nf; ++f) {
136: DSBoundary b;
137: PetscObject obj;
138: PetscClassId id;
139: PetscQuadrature q;
140: const char *name;
141: PetscInt Nc, Nq, Nqc;
143: PetscCall(PetscDSGetDiscretization(ds, f, &obj));
144: PetscCall(PetscObjectGetClassId(obj, &id));
145: PetscCall(PetscObjectGetName(obj, &name));
146: PetscCall(PetscViewerASCIIPrintf(viewer, "Field %s", name ? name : "<unknown>"));
147: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
148: if (id == PETSCFE_CLASSID) {
149: PetscCall(PetscFEGetNumComponents((PetscFE)obj, &Nc));
150: PetscCall(PetscFEGetQuadrature((PetscFE)obj, &q));
151: PetscCall(PetscViewerASCIIPrintf(viewer, " FEM"));
152: } else if (id == PETSCFV_CLASSID) {
153: PetscCall(PetscFVGetNumComponents((PetscFV)obj, &Nc));
154: PetscCall(PetscFVGetQuadrature((PetscFV)obj, &q));
155: PetscCall(PetscViewerASCIIPrintf(viewer, " FVM"));
156: } else SETERRQ(PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f);
157: if (Nc > 1) PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " components", Nc));
158: else PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " component ", Nc));
159: if (ds->implicit[f]) PetscCall(PetscViewerASCIIPrintf(viewer, " (implicit)"));
160: else PetscCall(PetscViewerASCIIPrintf(viewer, " (explicit)"));
161: if (q) {
162: PetscCall(PetscQuadratureGetData(q, NULL, &Nqc, &Nq, NULL, NULL));
163: PetscCall(PetscViewerASCIIPrintf(viewer, " (Nq %" PetscInt_FMT " Nqc %" PetscInt_FMT ")", Nq, Nqc));
164: }
165: PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT "-jet", ds->jetDegree[f]));
166: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
167: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
168: PetscCall(PetscViewerASCIIPushTab(viewer));
169: if (id == PETSCFE_CLASSID) PetscCall(PetscFEView((PetscFE)obj, viewer));
170: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVView((PetscFV)obj, viewer));
171: PetscCall(PetscViewerASCIIPopTab(viewer));
173: for (b = ds->boundary; b; b = b->next) {
174: char *name;
175: PetscInt c, i;
177: if (b->field != f) continue;
178: PetscCall(PetscViewerASCIIPushTab(viewer));
179: PetscCall(PetscViewerASCIIPrintf(viewer, "Boundary %s (%s) %s\n", b->name, b->lname, DMBoundaryConditionTypes[b->type]));
180: if (!b->Nc) {
181: PetscCall(PetscViewerASCIIPrintf(viewer, " all components\n"));
182: } else {
183: PetscCall(PetscViewerASCIIPrintf(viewer, " components: "));
184: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
185: for (c = 0; c < b->Nc; ++c) {
186: if (c > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", "));
187: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->comps[c]));
188: }
189: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
190: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
191: }
192: PetscCall(PetscViewerASCIIPrintf(viewer, " values: "));
193: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
194: for (i = 0; i < b->Nv; ++i) {
195: if (i > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", "));
196: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->values[i]));
197: }
198: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
199: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
200: #if defined(__clang__)
201: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat-pedantic")
202: #elif defined(__GNUC__) || defined(__GNUG__)
203: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat")
204: #endif
205: if (b->func) {
206: PetscCall(PetscDLAddr(b->func, &name));
207: if (name) PetscCall(PetscViewerASCIIPrintf(viewer, " func: %s\n", name));
208: else PetscCall(PetscViewerASCIIPrintf(viewer, " func: %p\n", b->func));
209: PetscCall(PetscFree(name));
210: }
211: if (b->func_t) {
212: PetscCall(PetscDLAddr(b->func_t, &name));
213: if (name) PetscCall(PetscViewerASCIIPrintf(viewer, " func_t: %s\n", name));
214: else PetscCall(PetscViewerASCIIPrintf(viewer, " func_t: %p\n", b->func_t));
215: PetscCall(PetscFree(name));
216: }
217: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_END()
218: PetscCall(PetscWeakFormView(b->wf, viewer));
219: PetscCall(PetscViewerASCIIPopTab(viewer));
220: }
221: }
222: PetscCall(PetscDSGetConstants(ds, &numConstants, &constants));
223: if (numConstants) {
224: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " constants\n", numConstants));
225: PetscCall(PetscViewerASCIIPushTab(viewer));
226: for (f = 0; f < numConstants; ++f) PetscCall(PetscViewerASCIIPrintf(viewer, "%g\n", (double)PetscRealPart(constants[f])));
227: PetscCall(PetscViewerASCIIPopTab(viewer));
228: }
229: PetscCall(PetscWeakFormView(ds->wf, viewer));
230: PetscCall(PetscViewerASCIIPopTab(viewer));
231: PetscFunctionReturn(PETSC_SUCCESS);
232: }
234: /*@C
235: PetscDSViewFromOptions - View a `PetscDS` based on values in the options database
237: Collective
239: Input Parameters:
240: + A - the `PetscDS` object
241: . obj - Optional object that provides the options prefix used in the search
242: - name - command line option
244: Level: intermediate
246: .seealso: `PetscDSType`, `PetscDS`, `PetscDSView()`, `PetscObjectViewFromOptions()`, `PetscDSCreate()`
247: @*/
248: PetscErrorCode PetscDSViewFromOptions(PetscDS A, PetscObject obj, const char name[])
249: {
250: PetscFunctionBegin;
252: PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
253: PetscFunctionReturn(PETSC_SUCCESS);
254: }
256: /*@C
257: PetscDSView - Views a `PetscDS`
259: Collective
261: Input Parameters:
262: + prob - the `PetscDS` object to view
263: - v - the viewer
265: Level: developer
267: .seealso: `PetscDSType`, `PetscDS`, `PetscViewer`, `PetscDSDestroy()`, `PetscDSViewFromOptions()`
268: @*/
269: PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v)
270: {
271: PetscBool iascii;
273: PetscFunctionBegin;
275: if (!v) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)prob), &v));
277: PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &iascii));
278: if (iascii) PetscCall(PetscDSView_Ascii(prob, v));
279: PetscTryTypeMethod(prob, view, v);
280: PetscFunctionReturn(PETSC_SUCCESS);
281: }
283: /*@
284: PetscDSSetFromOptions - sets parameters in a `PetscDS` from the options database
286: Collective
288: Input Parameter:
289: . prob - the `PetscDS` object to set options for
291: Options Database Keys:
292: + -petscds_type <type> - Set the `PetscDS` type
293: . -petscds_view <view opt> - View the `PetscDS`
294: . -petscds_jac_pre - Turn formation of a separate Jacobian preconditioner on or off
295: . -bc_<name> <ids> - Specify a list of label ids for a boundary condition
296: - -bc_<name>_comp <comps> - Specify a list of field components to constrain for a boundary condition
298: Level: intermediate
300: .seealso: `PetscDS`, `PetscDSView()`
301: @*/
302: PetscErrorCode PetscDSSetFromOptions(PetscDS prob)
303: {
304: DSBoundary b;
305: const char *defaultType;
306: char name[256];
307: PetscBool flg;
309: PetscFunctionBegin;
311: if (!((PetscObject)prob)->type_name) {
312: defaultType = PETSCDSBASIC;
313: } else {
314: defaultType = ((PetscObject)prob)->type_name;
315: }
316: PetscCall(PetscDSRegisterAll());
318: PetscObjectOptionsBegin((PetscObject)prob);
319: for (b = prob->boundary; b; b = b->next) {
320: char optname[1024];
321: PetscInt ids[1024], len = 1024;
322: PetscBool flg;
324: PetscCall(PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name));
325: PetscCall(PetscMemzero(ids, sizeof(ids)));
326: PetscCall(PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg));
327: if (flg) {
328: b->Nv = len;
329: PetscCall(PetscFree(b->values));
330: PetscCall(PetscMalloc1(len, &b->values));
331: PetscCall(PetscArraycpy(b->values, ids, len));
332: PetscCall(PetscWeakFormRewriteKeys(b->wf, b->label, len, b->values));
333: }
334: len = 1024;
335: PetscCall(PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name));
336: PetscCall(PetscMemzero(ids, sizeof(ids)));
337: PetscCall(PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg));
338: if (flg) {
339: b->Nc = len;
340: PetscCall(PetscFree(b->comps));
341: PetscCall(PetscMalloc1(len, &b->comps));
342: PetscCall(PetscArraycpy(b->comps, ids, len));
343: }
344: }
345: PetscCall(PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg));
346: if (flg) {
347: PetscCall(PetscDSSetType(prob, name));
348: } else if (!((PetscObject)prob)->type_name) {
349: PetscCall(PetscDSSetType(prob, defaultType));
350: }
351: PetscCall(PetscOptionsBool("-petscds_jac_pre", "Discrete System", "PetscDSUseJacobianPreconditioner", prob->useJacPre, &prob->useJacPre, &flg));
352: PetscCall(PetscOptionsBool("-petscds_force_quad", "Discrete System", "PetscDSSetForceQuad", prob->forceQuad, &prob->forceQuad, &flg));
353: 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: /*@C
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: *ds = NULL;
685: PetscCall(PetscDSInitializePackage());
687: PetscCall(PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView));
689: p->Nf = 0;
690: p->setup = PETSC_FALSE;
691: p->numConstants = 0;
692: p->constants = NULL;
693: p->dimEmbed = -1;
694: p->useJacPre = PETSC_TRUE;
695: p->forceQuad = PETSC_TRUE;
696: PetscCall(PetscWeakFormCreate(comm, &p->wf));
697: PetscCall(PetscArrayzero(p->quadPerm, DM_NUM_POLYTOPES));
699: *ds = p;
700: PetscFunctionReturn(PETSC_SUCCESS);
701: }
703: /*@
704: PetscDSGetNumFields - Returns the number of fields in the `PetscDS`
706: Not Collective
708: Input Parameter:
709: . prob - The `PetscDS` object
711: Output Parameter:
712: . Nf - The number of fields
714: Level: beginner
716: .seealso: `PetscDS`, `PetscDSGetSpatialDimension()`, `PetscDSCreate()`
717: @*/
718: PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf)
719: {
720: PetscFunctionBegin;
722: PetscAssertPointer(Nf, 2);
723: *Nf = prob->Nf;
724: PetscFunctionReturn(PETSC_SUCCESS);
725: }
727: /*@
728: PetscDSGetSpatialDimension - Returns the spatial dimension of the `PetscDS`, meaning the topological dimension of the discretizations
730: Not Collective
732: Input Parameter:
733: . prob - The `PetscDS` object
735: Output Parameter:
736: . dim - The spatial dimension
738: Level: beginner
740: .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
741: @*/
742: PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim)
743: {
744: PetscFunctionBegin;
746: PetscAssertPointer(dim, 2);
747: *dim = 0;
748: if (prob->Nf) {
749: PetscObject obj;
750: PetscClassId id;
752: PetscCall(PetscDSGetDiscretization(prob, 0, &obj));
753: if (obj) {
754: PetscCall(PetscObjectGetClassId(obj, &id));
755: if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetSpatialDimension((PetscFE)obj, dim));
756: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetSpatialDimension((PetscFV)obj, dim));
757: else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
758: }
759: }
760: PetscFunctionReturn(PETSC_SUCCESS);
761: }
763: /*@
764: PetscDSGetCoordinateDimension - Returns the coordinate dimension of the `PetscDS`, meaning the dimension of the space into which the discretiaztions are embedded
766: Not Collective
768: Input Parameter:
769: . prob - The `PetscDS` object
771: Output Parameter:
772: . dimEmbed - The coordinate dimension
774: Level: beginner
776: .seealso: `PetscDS`, `PetscDSSetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
777: @*/
778: PetscErrorCode PetscDSGetCoordinateDimension(PetscDS prob, PetscInt *dimEmbed)
779: {
780: PetscFunctionBegin;
782: PetscAssertPointer(dimEmbed, 2);
783: PetscCheck(prob->dimEmbed >= 0, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONGSTATE, "No coordinate dimension set for this DS");
784: *dimEmbed = prob->dimEmbed;
785: PetscFunctionReturn(PETSC_SUCCESS);
786: }
788: /*@
789: PetscDSSetCoordinateDimension - Set the coordinate dimension of the `PetscDS`, meaning the dimension of the space into which the discretiaztions are embedded
791: Logically Collective
793: Input Parameters:
794: + prob - The `PetscDS` object
795: - dimEmbed - The coordinate dimension
797: Level: beginner
799: .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
800: @*/
801: PetscErrorCode PetscDSSetCoordinateDimension(PetscDS prob, PetscInt dimEmbed)
802: {
803: PetscFunctionBegin;
805: PetscCheck(dimEmbed >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Coordinate dimension must be non-negative, not %" PetscInt_FMT, dimEmbed);
806: prob->dimEmbed = dimEmbed;
807: PetscFunctionReturn(PETSC_SUCCESS);
808: }
810: /*@
811: PetscDSGetForceQuad - Returns the flag to force matching quadratures among the field discretizations
813: Not collective
815: Input Parameter:
816: . ds - The `PetscDS` object
818: Output Parameter:
819: . forceQuad - The flag
821: Level: intermediate
823: .seealso: `PetscDS`, `PetscDSSetForceQuad()`, `PetscDSGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
824: @*/
825: PetscErrorCode PetscDSGetForceQuad(PetscDS ds, PetscBool *forceQuad)
826: {
827: PetscFunctionBegin;
829: PetscAssertPointer(forceQuad, 2);
830: *forceQuad = ds->forceQuad;
831: PetscFunctionReturn(PETSC_SUCCESS);
832: }
834: /*@
835: PetscDSSetForceQuad - Set the flag to force matching quadratures among the field discretizations
837: Logically collective on ds
839: Input Parameters:
840: + ds - The `PetscDS` object
841: - forceQuad - The flag
843: Level: intermediate
845: .seealso: `PetscDS`, `PetscDSGetForceQuad()`, `PetscDSGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
846: @*/
847: PetscErrorCode PetscDSSetForceQuad(PetscDS ds, PetscBool forceQuad)
848: {
849: PetscFunctionBegin;
851: ds->forceQuad = forceQuad;
852: PetscFunctionReturn(PETSC_SUCCESS);
853: }
855: /*@
856: PetscDSIsCohesive - Returns the flag indicating that this `PetscDS` is for a cohesive cell
858: Not Collective
860: Input Parameter:
861: . ds - The `PetscDS` object
863: Output Parameter:
864: . isCohesive - The flag
866: Level: developer
868: .seealso: `PetscDS`, `PetscDSGetNumCohesive()`, `PetscDSGetCohesive()`, `PetscDSSetCohesive()`, `PetscDSCreate()`
869: @*/
870: PetscErrorCode PetscDSIsCohesive(PetscDS ds, PetscBool *isCohesive)
871: {
872: PetscFunctionBegin;
874: PetscAssertPointer(isCohesive, 2);
875: *isCohesive = ds->isCohesive;
876: PetscFunctionReturn(PETSC_SUCCESS);
877: }
879: /*@
880: PetscDSGetNumCohesive - Returns the number of cohesive fields, meaning those defined on the interior of a cohesive cell
882: Not Collective
884: Input Parameter:
885: . ds - The `PetscDS` object
887: Output Parameter:
888: . numCohesive - The number of cohesive fields
890: Level: developer
892: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSCreate()`
893: @*/
894: PetscErrorCode PetscDSGetNumCohesive(PetscDS ds, PetscInt *numCohesive)
895: {
896: PetscInt f;
898: PetscFunctionBegin;
900: PetscAssertPointer(numCohesive, 2);
901: *numCohesive = 0;
902: for (f = 0; f < ds->Nf; ++f) *numCohesive += ds->cohesive[f] ? 1 : 0;
903: PetscFunctionReturn(PETSC_SUCCESS);
904: }
906: /*@
907: PetscDSGetCohesive - Returns the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell
909: Not Collective
911: Input Parameters:
912: + ds - The `PetscDS` object
913: - f - The field index
915: Output Parameter:
916: . isCohesive - The flag
918: Level: developer
920: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
921: @*/
922: PetscErrorCode PetscDSGetCohesive(PetscDS ds, PetscInt f, PetscBool *isCohesive)
923: {
924: PetscFunctionBegin;
926: PetscAssertPointer(isCohesive, 3);
927: 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);
928: *isCohesive = ds->cohesive[f];
929: PetscFunctionReturn(PETSC_SUCCESS);
930: }
932: /*@
933: PetscDSSetCohesive - Set the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell
935: Not Collective
937: Input Parameters:
938: + ds - The `PetscDS` object
939: . f - The field index
940: - isCohesive - The flag for a cohesive field
942: Level: developer
944: .seealso: `PetscDS`, `PetscDSGetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
945: @*/
946: PetscErrorCode PetscDSSetCohesive(PetscDS ds, PetscInt f, PetscBool isCohesive)
947: {
948: PetscInt i;
950: PetscFunctionBegin;
952: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
953: ds->cohesive[f] = isCohesive;
954: ds->isCohesive = PETSC_FALSE;
955: for (i = 0; i < ds->Nf; ++i) ds->isCohesive = ds->isCohesive || ds->cohesive[f] ? PETSC_TRUE : PETSC_FALSE;
956: PetscFunctionReturn(PETSC_SUCCESS);
957: }
959: /*@
960: PetscDSGetTotalDimension - Returns the total size of the approximation space for this system
962: Not Collective
964: Input Parameter:
965: . prob - The `PetscDS` object
967: Output Parameter:
968: . dim - The total problem dimension
970: Level: beginner
972: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
973: @*/
974: PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim)
975: {
976: PetscFunctionBegin;
978: PetscCall(PetscDSSetUp(prob));
979: PetscAssertPointer(dim, 2);
980: *dim = prob->totDim;
981: PetscFunctionReturn(PETSC_SUCCESS);
982: }
984: /*@
985: PetscDSGetTotalComponents - Returns the total number of components in this system
987: Not Collective
989: Input Parameter:
990: . prob - The `PetscDS` object
992: Output Parameter:
993: . Nc - The total number of components
995: Level: beginner
997: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
998: @*/
999: PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc)
1000: {
1001: PetscFunctionBegin;
1003: PetscCall(PetscDSSetUp(prob));
1004: PetscAssertPointer(Nc, 2);
1005: *Nc = prob->totComp;
1006: PetscFunctionReturn(PETSC_SUCCESS);
1007: }
1009: /*@
1010: PetscDSGetDiscretization - Returns the discretization object for the given field
1012: Not Collective
1014: Input Parameters:
1015: + prob - The `PetscDS` object
1016: - f - The field number
1018: Output Parameter:
1019: . disc - The discretization object
1021: Level: beginner
1023: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1024: @*/
1025: PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc)
1026: {
1027: PetscFunctionBeginHot;
1029: PetscAssertPointer(disc, 3);
1030: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
1031: *disc = prob->disc[f];
1032: PetscFunctionReturn(PETSC_SUCCESS);
1033: }
1035: /*@
1036: PetscDSSetDiscretization - Sets the discretization object for the given field
1038: Not Collective
1040: Input Parameters:
1041: + prob - The `PetscDS` object
1042: . f - The field number
1043: - disc - The discretization object
1045: Level: beginner
1047: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSGetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1048: @*/
1049: PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc)
1050: {
1051: PetscFunctionBegin;
1053: if (disc) PetscAssertPointer(disc, 3);
1054: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1055: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
1056: PetscCall(PetscObjectDereference(prob->disc[f]));
1057: prob->disc[f] = disc;
1058: PetscCall(PetscObjectReference(disc));
1059: if (disc) {
1060: PetscClassId id;
1062: PetscCall(PetscObjectGetClassId(disc, &id));
1063: if (id == PETSCFE_CLASSID) {
1064: PetscCall(PetscDSSetImplicit(prob, f, PETSC_TRUE));
1065: } else if (id == PETSCFV_CLASSID) {
1066: PetscCall(PetscDSSetImplicit(prob, f, PETSC_FALSE));
1067: }
1068: PetscCall(PetscDSSetJetDegree(prob, f, 1));
1069: }
1070: PetscFunctionReturn(PETSC_SUCCESS);
1071: }
1073: /*@
1074: PetscDSGetWeakForm - Returns the weak form object
1076: Not Collective
1078: Input Parameter:
1079: . ds - The `PetscDS` object
1081: Output Parameter:
1082: . wf - The weak form object
1084: Level: beginner
1086: .seealso: `PetscWeakForm`, `PetscDSSetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1087: @*/
1088: PetscErrorCode PetscDSGetWeakForm(PetscDS ds, PetscWeakForm *wf)
1089: {
1090: PetscFunctionBegin;
1092: PetscAssertPointer(wf, 2);
1093: *wf = ds->wf;
1094: PetscFunctionReturn(PETSC_SUCCESS);
1095: }
1097: /*@
1098: PetscDSSetWeakForm - Sets the weak form object
1100: Not Collective
1102: Input Parameters:
1103: + ds - The `PetscDS` object
1104: - wf - The weak form object
1106: Level: beginner
1108: .seealso: `PetscWeakForm`, `PetscDSGetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1109: @*/
1110: PetscErrorCode PetscDSSetWeakForm(PetscDS ds, PetscWeakForm wf)
1111: {
1112: PetscFunctionBegin;
1115: PetscCall(PetscObjectDereference((PetscObject)ds->wf));
1116: ds->wf = wf;
1117: PetscCall(PetscObjectReference((PetscObject)wf));
1118: PetscCall(PetscWeakFormSetNumFields(wf, ds->Nf));
1119: PetscFunctionReturn(PETSC_SUCCESS);
1120: }
1122: /*@
1123: PetscDSAddDiscretization - Adds a discretization object
1125: Not Collective
1127: Input Parameters:
1128: + prob - The `PetscDS` object
1129: - disc - The boundary discretization object
1131: Level: beginner
1133: .seealso: `PetscWeakForm`, `PetscDSGetDiscretization()`, `PetscDSSetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1134: @*/
1135: PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc)
1136: {
1137: PetscFunctionBegin;
1138: PetscCall(PetscDSSetDiscretization(prob, prob->Nf, disc));
1139: PetscFunctionReturn(PETSC_SUCCESS);
1140: }
1142: /*@
1143: PetscDSGetQuadrature - Returns the quadrature, which must agree for all fields in the `PetscDS`
1145: Not Collective
1147: Input Parameter:
1148: . prob - The `PetscDS` object
1150: Output Parameter:
1151: . q - The quadrature object
1153: Level: intermediate
1155: .seealso: `PetscDS`, `PetscQuadrature`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1156: @*/
1157: PetscErrorCode PetscDSGetQuadrature(PetscDS prob, PetscQuadrature *q)
1158: {
1159: PetscObject obj;
1160: PetscClassId id;
1162: PetscFunctionBegin;
1163: *q = NULL;
1164: if (!prob->Nf) PetscFunctionReturn(PETSC_SUCCESS);
1165: PetscCall(PetscDSGetDiscretization(prob, 0, &obj));
1166: PetscCall(PetscObjectGetClassId(obj, &id));
1167: if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetQuadrature((PetscFE)obj, q));
1168: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetQuadrature((PetscFV)obj, q));
1169: else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
1170: PetscFunctionReturn(PETSC_SUCCESS);
1171: }
1173: /*@
1174: PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for `TSIMEX`
1176: Not Collective
1178: Input Parameters:
1179: + prob - The `PetscDS` object
1180: - f - The field number
1182: Output Parameter:
1183: . implicit - The flag indicating what kind of solve to use for this field
1185: Level: developer
1187: .seealso: `TSIMEX`, `PetscDS`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1188: @*/
1189: PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit)
1190: {
1191: PetscFunctionBegin;
1193: PetscAssertPointer(implicit, 3);
1194: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
1195: *implicit = prob->implicit[f];
1196: PetscFunctionReturn(PETSC_SUCCESS);
1197: }
1199: /*@
1200: PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for `TSIMEX`
1202: Not Collective
1204: Input Parameters:
1205: + prob - The `PetscDS` object
1206: . f - The field number
1207: - implicit - The flag indicating what kind of solve to use for this field
1209: Level: developer
1211: .seealso: `TSIMEX`, `PetscDSGetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1212: @*/
1213: PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit)
1214: {
1215: PetscFunctionBegin;
1217: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
1218: prob->implicit[f] = implicit;
1219: PetscFunctionReturn(PETSC_SUCCESS);
1220: }
1222: /*@
1223: PetscDSGetJetDegree - Returns the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate.
1225: Not Collective
1227: Input Parameters:
1228: + ds - The `PetscDS` object
1229: - f - The field number
1231: Output Parameter:
1232: . k - The highest derivative we need to tabulate
1234: Level: developer
1236: .seealso: `PetscDS`, `PetscDSSetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1237: @*/
1238: PetscErrorCode PetscDSGetJetDegree(PetscDS ds, PetscInt f, PetscInt *k)
1239: {
1240: PetscFunctionBegin;
1242: PetscAssertPointer(k, 3);
1243: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1244: *k = ds->jetDegree[f];
1245: PetscFunctionReturn(PETSC_SUCCESS);
1246: }
1248: /*@
1249: PetscDSSetJetDegree - Set the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate.
1251: Not Collective
1253: Input Parameters:
1254: + ds - The `PetscDS` object
1255: . f - The field number
1256: - k - The highest derivative we need to tabulate
1258: Level: developer
1260: .seealso: `PetscDS`, `PetscDSGetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1261: @*/
1262: PetscErrorCode PetscDSSetJetDegree(PetscDS ds, PetscInt f, PetscInt k)
1263: {
1264: PetscFunctionBegin;
1266: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1267: ds->jetDegree[f] = k;
1268: PetscFunctionReturn(PETSC_SUCCESS);
1269: }
1271: /*@C
1272: PetscDSGetObjective - Get the pointwise objective function for a given test field
1274: Not Collective
1276: Input Parameters:
1277: + ds - The `PetscDS`
1278: - f - The test field number
1280: Output Parameter:
1281: . obj - integrand for the test function term
1283: Calling sequence of `obj`:
1284: + dim - the spatial dimension
1285: . Nf - the number of fields
1286: . NfAux - the number of auxiliary fields
1287: . uOff - the offset into u[] and u_t[] for each field
1288: . uOff_x - the offset into u_x[] for each field
1289: . u - each field evaluated at the current point
1290: . u_t - the time derivative of each field evaluated at the current point
1291: . u_x - the gradient of each field evaluated at the current point
1292: . aOff - the offset into a[] and a_t[] for each auxiliary field
1293: . aOff_x - the offset into a_x[] for each auxiliary field
1294: . a - each auxiliary field evaluated at the current point
1295: . a_t - the time derivative of each auxiliary field evaluated at the current point
1296: . a_x - the gradient of auxiliary each field evaluated at the current point
1297: . t - current time
1298: . x - coordinates of the current point
1299: . numConstants - number of constant parameters
1300: . constants - constant parameters
1301: - obj - output values at the current point
1303: Level: intermediate
1305: Note:
1306: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi obj(u, u_t, \nabla u, x, t)$
1308: .seealso: `PetscDS`, `PetscDSSetObjective()`, `PetscDSGetResidual()`
1309: @*/
1310: 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[]))
1311: {
1312: PetscPointFunc *tmp;
1313: PetscInt n;
1315: PetscFunctionBegin;
1317: PetscAssertPointer(obj, 3);
1318: 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);
1319: PetscCall(PetscWeakFormGetObjective(ds->wf, NULL, 0, f, 0, &n, &tmp));
1320: *obj = tmp ? tmp[0] : NULL;
1321: PetscFunctionReturn(PETSC_SUCCESS);
1322: }
1324: /*@C
1325: PetscDSSetObjective - Set the pointwise objective function for a given test field
1327: Not Collective
1329: Input Parameters:
1330: + ds - The `PetscDS`
1331: . f - The test field number
1332: - obj - integrand for the test function term
1334: Calling sequence of `obj`:
1335: + dim - the spatial dimension
1336: . Nf - the number of fields
1337: . NfAux - the number of auxiliary fields
1338: . uOff - the offset into u[] and u_t[] for each field
1339: . uOff_x - the offset into u_x[] for each field
1340: . u - each field evaluated at the current point
1341: . u_t - the time derivative of each field evaluated at the current point
1342: . u_x - the gradient of each field evaluated at the current point
1343: . aOff - the offset into a[] and a_t[] for each auxiliary field
1344: . aOff_x - the offset into a_x[] for each auxiliary field
1345: . a - each auxiliary field evaluated at the current point
1346: . a_t - the time derivative of each auxiliary field evaluated at the current point
1347: . a_x - the gradient of auxiliary each field evaluated at the current point
1348: . t - current time
1349: . x - coordinates of the current point
1350: . numConstants - number of constant parameters
1351: . constants - constant parameters
1352: - obj - output values at the current point
1354: Level: intermediate
1356: Note:
1357: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi obj(u, u_t, \nabla u, x, t)$
1359: .seealso: `PetscDS`, `PetscDSGetObjective()`, `PetscDSSetResidual()`
1360: @*/
1361: 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[]))
1362: {
1363: PetscFunctionBegin;
1366: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1367: PetscCall(PetscWeakFormSetIndexObjective(ds->wf, NULL, 0, f, 0, 0, obj));
1368: PetscFunctionReturn(PETSC_SUCCESS);
1369: }
1371: /*@C
1372: PetscDSGetResidual - Get the pointwise residual function for a given test field
1374: Not Collective
1376: Input Parameters:
1377: + ds - The `PetscDS`
1378: - f - The test field number
1380: Output Parameters:
1381: + f0 - integrand for the test function term
1382: - f1 - integrand for the test function gradient term
1384: Calling sequence of `f0`:
1385: + dim - the spatial dimension
1386: . Nf - the number of fields
1387: . NfAux - the number of auxiliary fields
1388: . uOff - the offset into u[] and u_t[] for each field
1389: . uOff_x - the offset into u_x[] for each field
1390: . u - each field evaluated at the current point
1391: . u_t - the time derivative of each field evaluated at the current point
1392: . u_x - the gradient of each field evaluated at the current point
1393: . aOff - the offset into a[] and a_t[] for each auxiliary field
1394: . aOff_x - the offset into a_x[] for each auxiliary field
1395: . a - each auxiliary field evaluated at the current point
1396: . a_t - the time derivative of each auxiliary field evaluated at the current point
1397: . a_x - the gradient of auxiliary each field evaluated at the current point
1398: . t - current time
1399: . x - coordinates of the current point
1400: . numConstants - number of constant parameters
1401: . constants - constant parameters
1402: - f0 - output values at the current point
1404: Level: intermediate
1406: Note:
1407: `f1` has an identical form and is omitted for brevity.
1409: 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)$
1411: .seealso: `PetscDS`, `PetscDSSetResidual()`
1412: @*/
1413: 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[]))
1414: {
1415: PetscPointFunc *tmp0, *tmp1;
1416: PetscInt n0, n1;
1418: PetscFunctionBegin;
1420: 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);
1421: PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1));
1422: *f0 = tmp0 ? tmp0[0] : NULL;
1423: *f1 = tmp1 ? tmp1[0] : NULL;
1424: PetscFunctionReturn(PETSC_SUCCESS);
1425: }
1427: /*@C
1428: PetscDSSetResidual - Set the pointwise residual function for a given test field
1430: Not Collective
1432: Input Parameters:
1433: + ds - The `PetscDS`
1434: . f - The test field number
1435: . f0 - integrand for the test function term
1436: - f1 - integrand for the test function gradient term
1438: Calling sequence of `f0`:
1439: + dim - the spatial dimension
1440: . Nf - the number of fields
1441: . NfAux - the number of auxiliary fields
1442: . uOff - the offset into u[] and u_t[] for each field
1443: . uOff_x - the offset into u_x[] for each field
1444: . u - each field evaluated at the current point
1445: . u_t - the time derivative of each field evaluated at the current point
1446: . u_x - the gradient of each field evaluated at the current point
1447: . aOff - the offset into a[] and a_t[] for each auxiliary field
1448: . aOff_x - the offset into a_x[] for each auxiliary field
1449: . a - each auxiliary field evaluated at the current point
1450: . a_t - the time derivative of each auxiliary field evaluated at the current point
1451: . a_x - the gradient of auxiliary each field evaluated at the current point
1452: . t - current time
1453: . x - coordinates of the current point
1454: . numConstants - number of constant parameters
1455: . constants - constant parameters
1456: - f0 - output values at the current point
1458: Level: intermediate
1460: Note:
1461: `f1` has an identical form and is omitted for brevity.
1463: 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)$
1465: .seealso: `PetscDS`, `PetscDSGetResidual()`
1466: @*/
1467: 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[]))
1468: {
1469: PetscFunctionBegin;
1473: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1474: PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1));
1475: PetscFunctionReturn(PETSC_SUCCESS);
1476: }
1478: /*@C
1479: PetscDSGetRHSResidual - Get the pointwise RHS residual function for explicit timestepping for a given test field
1481: Not Collective
1483: Input Parameters:
1484: + ds - The `PetscDS`
1485: - f - The test field number
1487: Output Parameters:
1488: + f0 - integrand for the test function term
1489: - f1 - integrand for the test function gradient term
1491: Calling sequence of `f0`:
1492: + dim - the spatial dimension
1493: . Nf - the number of fields
1494: . NfAux - the number of auxiliary fields
1495: . uOff - the offset into u[] and u_t[] for each field
1496: . uOff_x - the offset into u_x[] for each field
1497: . u - each field evaluated at the current point
1498: . u_t - the time derivative of each field evaluated at the current point
1499: . u_x - the gradient of each field evaluated at the current point
1500: . aOff - the offset into a[] and a_t[] for each auxiliary field
1501: . aOff_x - the offset into a_x[] for each auxiliary field
1502: . a - each auxiliary field evaluated at the current point
1503: . a_t - the time derivative of each auxiliary field evaluated at the current point
1504: . a_x - the gradient of auxiliary each field evaluated at the current point
1505: . t - current time
1506: . x - coordinates of the current point
1507: . numConstants - number of constant parameters
1508: . constants - constant parameters
1509: - f0 - output values at the current point
1511: Level: intermediate
1513: Note:
1514: `f1` has an identical form and is omitted for brevity.
1516: 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)$
1518: .seealso: `PetscDS`, `PetscDSSetRHSResidual()`
1519: @*/
1520: 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[]))
1521: {
1522: PetscPointFunc *tmp0, *tmp1;
1523: PetscInt n0, n1;
1525: PetscFunctionBegin;
1527: 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);
1528: PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 100, &n0, &tmp0, &n1, &tmp1));
1529: *f0 = tmp0 ? tmp0[0] : NULL;
1530: *f1 = tmp1 ? tmp1[0] : NULL;
1531: PetscFunctionReturn(PETSC_SUCCESS);
1532: }
1534: /*@C
1535: PetscDSSetRHSResidual - Set the pointwise residual function for explicit timestepping for a given test field
1537: Not Collective
1539: Input Parameters:
1540: + ds - The `PetscDS`
1541: . f - The test field number
1542: . f0 - integrand for the test function term
1543: - f1 - integrand for the test function gradient term
1545: Calling sequence for the callbacks `f0`:
1546: + dim - the spatial dimension
1547: . Nf - the number of fields
1548: . NfAux - the number of auxiliary fields
1549: . uOff - the offset into u[] and u_t[] for each field
1550: . uOff_x - the offset into u_x[] for each field
1551: . u - each field evaluated at the current point
1552: . u_t - the time derivative of each field evaluated at the current point
1553: . u_x - the gradient of each field evaluated at the current point
1554: . aOff - the offset into a[] and a_t[] for each auxiliary field
1555: . aOff_x - the offset into a_x[] for each auxiliary field
1556: . a - each auxiliary field evaluated at the current point
1557: . a_t - the time derivative of each auxiliary field evaluated at the current point
1558: . a_x - the gradient of auxiliary each field evaluated at the current point
1559: . t - current time
1560: . x - coordinates of the current point
1561: . numConstants - number of constant parameters
1562: . constants - constant parameters
1563: - f0 - output values at the current point
1565: Level: intermediate
1567: Note:
1568: `f1` has an identical form and is omitted for brevity.
1570: 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)$
1572: .seealso: `PetscDS`, `PetscDSGetResidual()`
1573: @*/
1574: 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[]))
1575: {
1576: PetscFunctionBegin;
1580: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1581: PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 100, 0, f0, 0, f1));
1582: PetscFunctionReturn(PETSC_SUCCESS);
1583: }
1585: /*@C
1586: PetscDSHasJacobian - Checks that the Jacobian functions have been set
1588: Not Collective
1590: Input Parameter:
1591: . ds - The `PetscDS`
1593: Output Parameter:
1594: . hasJac - flag that pointwise function for the Jacobian has been set
1596: Level: intermediate
1598: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1599: @*/
1600: PetscErrorCode PetscDSHasJacobian(PetscDS ds, PetscBool *hasJac)
1601: {
1602: PetscFunctionBegin;
1604: PetscCall(PetscWeakFormHasJacobian(ds->wf, hasJac));
1605: PetscFunctionReturn(PETSC_SUCCESS);
1606: }
1608: /*@C
1609: PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field
1611: Not Collective
1613: Input Parameters:
1614: + ds - The `PetscDS`
1615: . f - The test field number
1616: - g - The field number
1618: Output Parameters:
1619: + g0 - integrand for the test and basis function term
1620: . g1 - integrand for the test function and basis function gradient term
1621: . g2 - integrand for the test function gradient and basis function term
1622: - g3 - integrand for the test function gradient and basis function gradient term
1624: Calling sequence of `g0`:
1625: + dim - the spatial dimension
1626: . Nf - the number of fields
1627: . NfAux - the number of auxiliary fields
1628: . uOff - the offset into u[] and u_t[] for each field
1629: . uOff_x - the offset into u_x[] for each field
1630: . u - each field evaluated at the current point
1631: . u_t - the time derivative of each field evaluated at the current point
1632: . u_x - the gradient of each field evaluated at the current point
1633: . aOff - the offset into a[] and a_t[] for each auxiliary field
1634: . aOff_x - the offset into a_x[] for each auxiliary field
1635: . a - each auxiliary field evaluated at the current point
1636: . a_t - the time derivative of each auxiliary field evaluated at the current point
1637: . a_x - the gradient of auxiliary each field evaluated at the current point
1638: . t - current time
1639: . u_tShift - the multiplier a for dF/dU_t
1640: . x - coordinates of the current point
1641: . numConstants - number of constant parameters
1642: . constants - constant parameters
1643: - g0 - output values at the current point
1645: Level: intermediate
1647: Note:
1648: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1650: We are using a first order FEM model for the weak form\:
1652: $$
1653: \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1654: $$
1656: .seealso: `PetscDS`, `PetscDSSetJacobian()`
1657: @*/
1658: 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[]))
1659: {
1660: PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1661: PetscInt n0, n1, n2, n3;
1663: PetscFunctionBegin;
1665: 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);
1666: 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);
1667: PetscCall(PetscWeakFormGetJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1668: *g0 = tmp0 ? tmp0[0] : NULL;
1669: *g1 = tmp1 ? tmp1[0] : NULL;
1670: *g2 = tmp2 ? tmp2[0] : NULL;
1671: *g3 = tmp3 ? tmp3[0] : NULL;
1672: PetscFunctionReturn(PETSC_SUCCESS);
1673: }
1675: /*@C
1676: PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields
1678: Not Collective
1680: Input Parameters:
1681: + ds - The `PetscDS`
1682: . f - The test field number
1683: . g - The field number
1684: . g0 - integrand for the test and basis function term
1685: . g1 - integrand for the test function and basis function gradient term
1686: . g2 - integrand for the test function gradient and basis function term
1687: - g3 - integrand for the test function gradient and basis function gradient term
1689: Calling sequence of `g0`:
1690: + dim - the spatial dimension
1691: . Nf - the number of fields
1692: . NfAux - the number of auxiliary fields
1693: . uOff - the offset into u[] and u_t[] for each field
1694: . uOff_x - the offset into u_x[] for each field
1695: . u - each field evaluated at the current point
1696: . u_t - the time derivative of each field evaluated at the current point
1697: . u_x - the gradient of each field evaluated at the current point
1698: . aOff - the offset into a[] and a_t[] for each auxiliary field
1699: . aOff_x - the offset into a_x[] for each auxiliary field
1700: . a - each auxiliary field evaluated at the current point
1701: . a_t - the time derivative of each auxiliary field evaluated at the current point
1702: . a_x - the gradient of auxiliary each field evaluated at the current point
1703: . t - current time
1704: . u_tShift - the multiplier a for dF/dU_t
1705: . x - coordinates of the current point
1706: . numConstants - number of constant parameters
1707: . constants - constant parameters
1708: - g0 - output values at the current point
1710: Level: intermediate
1712: Note:
1713: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1715: We are using a first order FEM model for the weak form\:
1716: \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1718: .seealso: `PetscDS`, `PetscDSGetJacobian()`
1719: @*/
1720: 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[]))
1721: {
1722: PetscFunctionBegin;
1728: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1729: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1730: PetscCall(PetscWeakFormSetIndexJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1731: PetscFunctionReturn(PETSC_SUCCESS);
1732: }
1734: /*@C
1735: PetscDSUseJacobianPreconditioner - Set whether to construct a Jacobian preconditioner
1737: Not Collective
1739: Input Parameters:
1740: + prob - The `PetscDS`
1741: - useJacPre - flag that enables construction of a Jacobian preconditioner
1743: Level: intermediate
1745: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1746: @*/
1747: PetscErrorCode PetscDSUseJacobianPreconditioner(PetscDS prob, PetscBool useJacPre)
1748: {
1749: PetscFunctionBegin;
1751: prob->useJacPre = useJacPre;
1752: PetscFunctionReturn(PETSC_SUCCESS);
1753: }
1755: /*@C
1756: PetscDSHasJacobianPreconditioner - Checks if a Jacobian preconditioner matrix has been set
1758: Not Collective
1760: Input Parameter:
1761: . ds - The `PetscDS`
1763: Output Parameter:
1764: . hasJacPre - flag that pointwise function for Jacobian preconditioner matrix has been set
1766: Level: intermediate
1768: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1769: @*/
1770: PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS ds, PetscBool *hasJacPre)
1771: {
1772: PetscFunctionBegin;
1774: *hasJacPre = PETSC_FALSE;
1775: if (!ds->useJacPre) PetscFunctionReturn(PETSC_SUCCESS);
1776: PetscCall(PetscWeakFormHasJacobianPreconditioner(ds->wf, hasJacPre));
1777: PetscFunctionReturn(PETSC_SUCCESS);
1778: }
1780: /*@C
1781: PetscDSGetJacobianPreconditioner - Get the pointwise Jacobian preconditioner function for given test and basis field. If this is missing,
1782: the system matrix is used to build the preconditioner.
1784: Not Collective
1786: Input Parameters:
1787: + ds - The `PetscDS`
1788: . f - The test field number
1789: - g - The field number
1791: Output Parameters:
1792: + g0 - integrand for the test and basis function term
1793: . g1 - integrand for the test function and basis function gradient term
1794: . g2 - integrand for the test function gradient and basis function term
1795: - g3 - integrand for the test function gradient and basis function gradient term
1797: Calling sequence of `g0`:
1798: + dim - the spatial dimension
1799: . Nf - the number of fields
1800: . NfAux - the number of auxiliary fields
1801: . uOff - the offset into u[] and u_t[] for each field
1802: . uOff_x - the offset into u_x[] for each field
1803: . u - each field evaluated at the current point
1804: . u_t - the time derivative of each field evaluated at the current point
1805: . u_x - the gradient of each field evaluated at the current point
1806: . aOff - the offset into a[] and a_t[] for each auxiliary field
1807: . aOff_x - the offset into a_x[] for each auxiliary field
1808: . a - each auxiliary field evaluated at the current point
1809: . a_t - the time derivative of each auxiliary field evaluated at the current point
1810: . a_x - the gradient of auxiliary each field evaluated at the current point
1811: . t - current time
1812: . u_tShift - the multiplier a for dF/dU_t
1813: . x - coordinates of the current point
1814: . numConstants - number of constant parameters
1815: . constants - constant parameters
1816: - g0 - output values at the current point
1818: Level: intermediate
1820: Note:
1821: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1822: We are using a first order FEM model for the weak form\:
1823: \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1825: .seealso: `PetscDS`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1826: @*/
1827: 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[]))
1828: {
1829: PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1830: PetscInt n0, n1, n2, n3;
1832: PetscFunctionBegin;
1834: 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);
1835: 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);
1836: PetscCall(PetscWeakFormGetJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1837: *g0 = tmp0 ? tmp0[0] : NULL;
1838: *g1 = tmp1 ? tmp1[0] : NULL;
1839: *g2 = tmp2 ? tmp2[0] : NULL;
1840: *g3 = tmp3 ? tmp3[0] : NULL;
1841: PetscFunctionReturn(PETSC_SUCCESS);
1842: }
1844: /*@C
1845: PetscDSSetJacobianPreconditioner - Set the pointwise Jacobian preconditioner function for given test and basis fields.
1846: If this is missing, the system matrix is used to build the preconditioner.
1848: Not Collective
1850: Input Parameters:
1851: + ds - The `PetscDS`
1852: . f - The test field number
1853: . g - The field number
1854: . g0 - integrand for the test and basis function term
1855: . g1 - integrand for the test function and basis function gradient term
1856: . g2 - integrand for the test function gradient and basis function term
1857: - g3 - integrand for the test function gradient and basis function gradient term
1859: Calling sequence of `g0`:
1860: + dim - the spatial dimension
1861: . Nf - the number of fields
1862: . NfAux - the number of auxiliary fields
1863: . uOff - the offset into u[] and u_t[] for each field
1864: . uOff_x - the offset into u_x[] for each field
1865: . u - each field evaluated at the current point
1866: . u_t - the time derivative of each field evaluated at the current point
1867: . u_x - the gradient of each field evaluated at the current point
1868: . aOff - the offset into a[] and a_t[] for each auxiliary field
1869: . aOff_x - the offset into a_x[] for each auxiliary field
1870: . a - each auxiliary field evaluated at the current point
1871: . a_t - the time derivative of each auxiliary field evaluated at the current point
1872: . a_x - the gradient of auxiliary each field evaluated at the current point
1873: . t - current time
1874: . u_tShift - the multiplier a for dF/dU_t
1875: . x - coordinates of the current point
1876: . numConstants - number of constant parameters
1877: . constants - constant parameters
1878: - g0 - output values at the current point
1880: Level: intermediate
1882: Note:
1883: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1885: We are using a first order FEM model for the weak form\:
1886: \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1888: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobian()`
1889: @*/
1890: 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[]))
1891: {
1892: PetscFunctionBegin;
1898: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1899: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1900: PetscCall(PetscWeakFormSetIndexJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1901: PetscFunctionReturn(PETSC_SUCCESS);
1902: }
1904: /*@C
1905: PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, dF/du_t, has been set
1907: Not Collective
1909: Input Parameter:
1910: . ds - The `PetscDS`
1912: Output Parameter:
1913: . hasDynJac - flag that pointwise function for dynamic Jacobian has been set
1915: Level: intermediate
1917: .seealso: `PetscDS`, `PetscDSGetDynamicJacobian()`, `PetscDSSetDynamicJacobian()`, `PetscDSGetJacobian()`
1918: @*/
1919: PetscErrorCode PetscDSHasDynamicJacobian(PetscDS ds, PetscBool *hasDynJac)
1920: {
1921: PetscFunctionBegin;
1923: PetscCall(PetscWeakFormHasDynamicJacobian(ds->wf, hasDynJac));
1924: PetscFunctionReturn(PETSC_SUCCESS);
1925: }
1927: /*@C
1928: PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, dF/du_t, function for given test and basis field
1930: Not Collective
1932: Input Parameters:
1933: + ds - The `PetscDS`
1934: . f - The test field number
1935: - g - The field number
1937: Output Parameters:
1938: + g0 - integrand for the test and basis function term
1939: . g1 - integrand for the test function and basis function gradient term
1940: . g2 - integrand for the test function gradient and basis function term
1941: - g3 - integrand for the test function gradient and basis function gradient term
1943: Calling sequence of `g0`:
1944: + dim - the spatial dimension
1945: . Nf - the number of fields
1946: . NfAux - the number of auxiliary fields
1947: . uOff - the offset into u[] and u_t[] for each field
1948: . uOff_x - the offset into u_x[] for each field
1949: . u - each field evaluated at the current point
1950: . u_t - the time derivative of each field evaluated at the current point
1951: . u_x - the gradient of each field evaluated at the current point
1952: . aOff - the offset into a[] and a_t[] for each auxiliary field
1953: . aOff_x - the offset into a_x[] for each auxiliary field
1954: . a - each auxiliary field evaluated at the current point
1955: . a_t - the time derivative of each auxiliary field evaluated at the current point
1956: . a_x - the gradient of auxiliary each field evaluated at the current point
1957: . t - current time
1958: . u_tShift - the multiplier a for dF/dU_t
1959: . x - coordinates of the current point
1960: . numConstants - number of constant parameters
1961: . constants - constant parameters
1962: - g0 - output values at the current point
1964: Level: intermediate
1966: Note:
1967: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1969: We are using a first order FEM model for the weak form\:
1970: \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1972: .seealso: `PetscDS`, `PetscDSSetJacobian()`
1973: @*/
1974: 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[]))
1975: {
1976: PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1977: PetscInt n0, n1, n2, n3;
1979: PetscFunctionBegin;
1981: 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);
1982: 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);
1983: PetscCall(PetscWeakFormGetDynamicJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1984: *g0 = tmp0 ? tmp0[0] : NULL;
1985: *g1 = tmp1 ? tmp1[0] : NULL;
1986: *g2 = tmp2 ? tmp2[0] : NULL;
1987: *g3 = tmp3 ? tmp3[0] : NULL;
1988: PetscFunctionReturn(PETSC_SUCCESS);
1989: }
1991: /*@C
1992: PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, dF/du_t, function for given test and basis fields
1994: Not Collective
1996: Input Parameters:
1997: + ds - The `PetscDS`
1998: . f - The test field number
1999: . g - The field number
2000: . g0 - integrand for the test and basis function term
2001: . g1 - integrand for the test function and basis function gradient term
2002: . g2 - integrand for the test function gradient and basis function term
2003: - g3 - integrand for the test function gradient and basis function gradient term
2005: Calling sequence of `g0`:
2006: + dim - the spatial dimension
2007: . Nf - the number of fields
2008: . NfAux - the number of auxiliary fields
2009: . uOff - the offset into u[] and u_t[] for each field
2010: . uOff_x - the offset into u_x[] for each field
2011: . u - each field evaluated at the current point
2012: . u_t - the time derivative of each field evaluated at the current point
2013: . u_x - the gradient of each field evaluated at the current point
2014: . aOff - the offset into a[] and a_t[] for each auxiliary field
2015: . aOff_x - the offset into a_x[] for each auxiliary field
2016: . a - each auxiliary field evaluated at the current point
2017: . a_t - the time derivative of each auxiliary field evaluated at the current point
2018: . a_x - the gradient of auxiliary each field evaluated at the current point
2019: . t - current time
2020: . u_tShift - the multiplier a for dF/dU_t
2021: . x - coordinates of the current point
2022: . numConstants - number of constant parameters
2023: . constants - constant parameters
2024: - g0 - output values at the current point
2026: Level: intermediate
2028: Note:
2029: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2031: We are using a first order FEM model for the weak form\:
2032: \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
2034: .seealso: `PetscDS`, `PetscDSGetJacobian()`
2035: @*/
2036: 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[]))
2037: {
2038: PetscFunctionBegin;
2044: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2045: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2046: PetscCall(PetscWeakFormSetIndexDynamicJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2047: PetscFunctionReturn(PETSC_SUCCESS);
2048: }
2050: /*@C
2051: PetscDSGetRiemannSolver - Returns the Riemann solver for the given field
2053: Not Collective
2055: Input Parameters:
2056: + ds - The `PetscDS` object
2057: - f - The field number
2059: Output Parameter:
2060: . r - Riemann solver
2062: Calling sequence of `r`:
2063: + dim - The spatial dimension
2064: . Nf - The number of fields
2065: . x - The coordinates at a point on the interface
2066: . n - The normal vector to the interface
2067: . uL - The state vector to the left of the interface
2068: . uR - The state vector to the right of the interface
2069: . flux - output array of flux through the interface
2070: . numConstants - number of constant parameters
2071: . constants - constant parameters
2072: - ctx - optional user context
2074: Level: intermediate
2076: .seealso: `PetscDS`, `PetscDSSetRiemannSolver()`
2077: @*/
2078: 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))
2079: {
2080: PetscRiemannFunc *tmp;
2081: PetscInt n;
2083: PetscFunctionBegin;
2085: PetscAssertPointer(r, 3);
2086: 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);
2087: PetscCall(PetscWeakFormGetRiemannSolver(ds->wf, NULL, 0, f, 0, &n, &tmp));
2088: *r = tmp ? tmp[0] : NULL;
2089: PetscFunctionReturn(PETSC_SUCCESS);
2090: }
2092: /*@C
2093: PetscDSSetRiemannSolver - Sets the Riemann solver for the given field
2095: Not Collective
2097: Input Parameters:
2098: + ds - The `PetscDS` object
2099: . f - The field number
2100: - r - Riemann solver
2102: Calling sequence of `r`:
2103: + dim - The spatial dimension
2104: . Nf - The number of fields
2105: . x - The coordinates at a point on the interface
2106: . n - The normal vector to the interface
2107: . uL - The state vector to the left of the interface
2108: . uR - The state vector to the right of the interface
2109: . flux - output array of flux through the interface
2110: . numConstants - number of constant parameters
2111: . constants - constant parameters
2112: - ctx - optional user context
2114: Level: intermediate
2116: .seealso: `PetscDS`, `PetscDSGetRiemannSolver()`
2117: @*/
2118: 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))
2119: {
2120: PetscFunctionBegin;
2123: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2124: PetscCall(PetscWeakFormSetIndexRiemannSolver(ds->wf, NULL, 0, f, 0, 0, r));
2125: PetscFunctionReturn(PETSC_SUCCESS);
2126: }
2128: /*@C
2129: PetscDSGetUpdate - Get the pointwise update function for a given field
2131: Not Collective
2133: Input Parameters:
2134: + ds - The `PetscDS`
2135: - f - The field number
2137: Output Parameter:
2138: . update - update function
2140: Calling sequence of `update`:
2141: + dim - the spatial dimension
2142: . Nf - the number of fields
2143: . NfAux - the number of auxiliary fields
2144: . uOff - the offset into u[] and u_t[] for each field
2145: . uOff_x - the offset into u_x[] for each field
2146: . u - each field evaluated at the current point
2147: . u_t - the time derivative of each field evaluated at the current point
2148: . u_x - the gradient of each field evaluated at the current point
2149: . aOff - the offset into a[] and a_t[] for each auxiliary field
2150: . aOff_x - the offset into a_x[] for each auxiliary field
2151: . a - each auxiliary field evaluated at the current point
2152: . a_t - the time derivative of each auxiliary field evaluated at the current point
2153: . a_x - the gradient of auxiliary each field evaluated at the current point
2154: . t - current time
2155: . x - coordinates of the current point
2156: . numConstants - number of constant parameters
2157: . constants - constant parameters
2158: - uNew - new value for field at the current point
2160: Level: intermediate
2162: .seealso: `PetscDS`, `PetscDSSetUpdate()`, `PetscDSSetResidual()`
2163: @*/
2164: 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[]))
2165: {
2166: PetscFunctionBegin;
2168: 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);
2169: if (update) {
2170: PetscAssertPointer(update, 3);
2171: *update = ds->update[f];
2172: }
2173: PetscFunctionReturn(PETSC_SUCCESS);
2174: }
2176: /*@C
2177: PetscDSSetUpdate - Set the pointwise update function for a given field
2179: Not Collective
2181: Input Parameters:
2182: + ds - The `PetscDS`
2183: . f - The field number
2184: - update - update function
2186: Calling sequence of `update`:
2187: + dim - the spatial dimension
2188: . Nf - the number of fields
2189: . NfAux - the number of auxiliary fields
2190: . uOff - the offset into u[] and u_t[] for each field
2191: . uOff_x - the offset into u_x[] for each field
2192: . u - each field evaluated at the current point
2193: . u_t - the time derivative of each field evaluated at the current point
2194: . u_x - the gradient of each field evaluated at the current point
2195: . aOff - the offset into a[] and a_t[] for each auxiliary field
2196: . aOff_x - the offset into a_x[] for each auxiliary field
2197: . a - each auxiliary field evaluated at the current point
2198: . a_t - the time derivative of each auxiliary field evaluated at the current point
2199: . a_x - the gradient of auxiliary each field evaluated at the current point
2200: . t - current time
2201: . x - coordinates of the current point
2202: . numConstants - number of constant parameters
2203: . constants - constant parameters
2204: - uNew - new field values at the current point
2206: Level: intermediate
2208: .seealso: `PetscDS`, `PetscDSGetResidual()`
2209: @*/
2210: 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[]))
2211: {
2212: PetscFunctionBegin;
2215: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2216: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
2217: ds->update[f] = update;
2218: PetscFunctionReturn(PETSC_SUCCESS);
2219: }
2221: PetscErrorCode PetscDSGetContext(PetscDS ds, PetscInt f, void *ctx)
2222: {
2223: PetscFunctionBegin;
2225: 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);
2226: PetscAssertPointer(ctx, 3);
2227: *(void **)ctx = ds->ctx[f];
2228: PetscFunctionReturn(PETSC_SUCCESS);
2229: }
2231: PetscErrorCode PetscDSSetContext(PetscDS ds, PetscInt f, void *ctx)
2232: {
2233: PetscFunctionBegin;
2235: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2236: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
2237: ds->ctx[f] = ctx;
2238: PetscFunctionReturn(PETSC_SUCCESS);
2239: }
2241: /*@C
2242: PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field
2244: Not Collective
2246: Input Parameters:
2247: + ds - The PetscDS
2248: - f - The test field number
2250: Output Parameters:
2251: + f0 - boundary integrand for the test function term
2252: - f1 - boundary integrand for the test function gradient term
2254: Calling sequence of `f0`:
2255: + dim - the spatial dimension
2256: . Nf - the number of fields
2257: . NfAux - the number of auxiliary fields
2258: . uOff - the offset into u[] and u_t[] for each field
2259: . uOff_x - the offset into u_x[] for each field
2260: . u - each field evaluated at the current point
2261: . u_t - the time derivative of each field evaluated at the current point
2262: . u_x - the gradient of each field evaluated at the current point
2263: . aOff - the offset into a[] and a_t[] for each auxiliary field
2264: . aOff_x - the offset into a_x[] for each auxiliary field
2265: . a - each auxiliary field evaluated at the current point
2266: . a_t - the time derivative of each auxiliary field evaluated at the current point
2267: . a_x - the gradient of auxiliary each field evaluated at the current point
2268: . t - current time
2269: . x - coordinates of the current point
2270: . n - unit normal at the current point
2271: . numConstants - number of constant parameters
2272: . constants - constant parameters
2273: - f0 - output values at the current point
2275: Level: intermediate
2277: Note:
2278: The calling sequence of `f1` is identical, and therefore omitted for brevity.
2280: We are using a first order FEM model for the weak form\:
2281: \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
2283: .seealso: `PetscDS`, `PetscDSSetBdResidual()`
2284: @*/
2285: 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[]))
2286: {
2287: PetscBdPointFunc *tmp0, *tmp1;
2288: PetscInt n0, n1;
2290: PetscFunctionBegin;
2292: 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);
2293: PetscCall(PetscWeakFormGetBdResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1));
2294: *f0 = tmp0 ? tmp0[0] : NULL;
2295: *f1 = tmp1 ? tmp1[0] : NULL;
2296: PetscFunctionReturn(PETSC_SUCCESS);
2297: }
2299: /*@C
2300: PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field
2302: Not Collective
2304: Input Parameters:
2305: + ds - The `PetscDS`
2306: . f - The test field number
2307: . f0 - boundary integrand for the test function term
2308: - f1 - boundary integrand for the test function gradient term
2310: Calling sequence of `f0`:
2311: + dim - the spatial dimension
2312: . Nf - the number of fields
2313: . NfAux - the number of auxiliary fields
2314: . uOff - the offset into u[] and u_t[] for each field
2315: . uOff_x - the offset into u_x[] for each field
2316: . u - each field evaluated at the current point
2317: . u_t - the time derivative of each field evaluated at the current point
2318: . u_x - the gradient of each field evaluated at the current point
2319: . aOff - the offset into a[] and a_t[] for each auxiliary field
2320: . aOff_x - the offset into a_x[] for each auxiliary field
2321: . a - each auxiliary field evaluated at the current point
2322: . a_t - the time derivative of each auxiliary field evaluated at the current point
2323: . a_x - the gradient of auxiliary each field evaluated at the current point
2324: . t - current time
2325: . x - coordinates of the current point
2326: . n - unit normal at the current point
2327: . numConstants - number of constant parameters
2328: . constants - constant parameters
2329: - f0 - output values at the current point
2331: Level: intermediate
2333: Note:
2334: The calling sequence of `f1` is identical, and therefore omitted for brevity.
2336: We are using a first order FEM model for the weak form\:
2337: \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
2339: .seealso: `PetscDS`, `PetscDSGetBdResidual()`
2340: @*/
2341: 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[]))
2342: {
2343: PetscFunctionBegin;
2345: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2346: PetscCall(PetscWeakFormSetIndexBdResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1));
2347: PetscFunctionReturn(PETSC_SUCCESS);
2348: }
2350: /*@
2351: PetscDSHasBdJacobian - Indicates that boundary Jacobian functions have been set
2353: Not Collective
2355: Input Parameter:
2356: . ds - The `PetscDS`
2358: Output Parameter:
2359: . hasBdJac - flag that pointwise function for the boundary Jacobian has been set
2361: Level: intermediate
2363: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2364: @*/
2365: PetscErrorCode PetscDSHasBdJacobian(PetscDS ds, PetscBool *hasBdJac)
2366: {
2367: PetscFunctionBegin;
2369: PetscAssertPointer(hasBdJac, 2);
2370: PetscCall(PetscWeakFormHasBdJacobian(ds->wf, hasBdJac));
2371: PetscFunctionReturn(PETSC_SUCCESS);
2372: }
2374: /*@C
2375: PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field
2377: Not Collective
2379: Input Parameters:
2380: + ds - The `PetscDS`
2381: . f - The test field number
2382: - g - The field number
2384: Output Parameters:
2385: + g0 - integrand for the test and basis function term
2386: . g1 - integrand for the test function and basis function gradient term
2387: . g2 - integrand for the test function gradient and basis function term
2388: - g3 - integrand for the test function gradient and basis function gradient term
2390: Calling sequence of `g0`:
2391: + dim - the spatial dimension
2392: . Nf - the number of fields
2393: . NfAux - the number of auxiliary fields
2394: . uOff - the offset into u[] and u_t[] for each field
2395: . uOff_x - the offset into u_x[] for each field
2396: . u - each field evaluated at the current point
2397: . u_t - the time derivative of each field evaluated at the current point
2398: . u_x - the gradient of each field evaluated at the current point
2399: . aOff - the offset into a[] and a_t[] for each auxiliary field
2400: . aOff_x - the offset into a_x[] for each auxiliary field
2401: . a - each auxiliary field evaluated at the current point
2402: . a_t - the time derivative of each auxiliary field evaluated at the current point
2403: . a_x - the gradient of auxiliary each field evaluated at the current point
2404: . t - current time
2405: . u_tShift - the multiplier a for dF/dU_t
2406: . x - coordinates of the current point
2407: . n - normal at the current point
2408: . numConstants - number of constant parameters
2409: . constants - constant parameters
2410: - g0 - output values at the current point
2412: Level: intermediate
2414: Note:
2415: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2417: We are using a first order FEM model for the weak form\:
2418: \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi
2420: .seealso: `PetscDS`, `PetscDSSetBdJacobian()`
2421: @*/
2422: 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[]))
2423: {
2424: PetscBdPointJac *tmp0, *tmp1, *tmp2, *tmp3;
2425: PetscInt n0, n1, n2, n3;
2427: PetscFunctionBegin;
2429: 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);
2430: 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);
2431: PetscCall(PetscWeakFormGetBdJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2432: *g0 = tmp0 ? tmp0[0] : NULL;
2433: *g1 = tmp1 ? tmp1[0] : NULL;
2434: *g2 = tmp2 ? tmp2[0] : NULL;
2435: *g3 = tmp3 ? tmp3[0] : NULL;
2436: PetscFunctionReturn(PETSC_SUCCESS);
2437: }
2439: /*@C
2440: PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field
2442: Not Collective
2444: Input Parameters:
2445: + ds - The PetscDS
2446: . f - The test field number
2447: . g - The field number
2448: . g0 - integrand for the test and basis function term
2449: . g1 - integrand for the test function and basis function gradient term
2450: . g2 - integrand for the test function gradient and basis function term
2451: - g3 - integrand for the test function gradient and basis function gradient term
2453: Calling sequence of `g0`:
2454: + dim - the spatial dimension
2455: . Nf - the number of fields
2456: . NfAux - the number of auxiliary fields
2457: . uOff - the offset into u[] and u_t[] for each field
2458: . uOff_x - the offset into u_x[] for each field
2459: . u - each field evaluated at the current point
2460: . u_t - the time derivative of each field evaluated at the current point
2461: . u_x - the gradient of each field evaluated at the current point
2462: . aOff - the offset into a[] and a_t[] for each auxiliary field
2463: . aOff_x - the offset into a_x[] for each auxiliary field
2464: . a - each auxiliary field evaluated at the current point
2465: . a_t - the time derivative of each auxiliary field evaluated at the current point
2466: . a_x - the gradient of auxiliary each field evaluated at the current point
2467: . t - current time
2468: . u_tShift - the multiplier a for dF/dU_t
2469: . x - coordinates of the current point
2470: . n - normal at the current point
2471: . numConstants - number of constant parameters
2472: . constants - constant parameters
2473: - g0 - output values at the current point
2475: Level: intermediate
2477: Note:
2478: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2480: We are using a first order FEM model for the weak form\:
2481: \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi
2483: .seealso: `PetscDS`, `PetscDSGetBdJacobian()`
2484: @*/
2485: 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[]))
2486: {
2487: PetscFunctionBegin;
2493: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2494: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2495: PetscCall(PetscWeakFormSetIndexBdJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2496: PetscFunctionReturn(PETSC_SUCCESS);
2497: }
2499: /*@
2500: PetscDSHasBdJacobianPreconditioner - Signals that boundary Jacobian preconditioner functions have been set
2502: Not Collective
2504: Input Parameter:
2505: . ds - The `PetscDS`
2507: Output Parameter:
2508: . hasBdJacPre - flag that pointwise function for the boundary Jacobian preconditioner has been set
2510: Level: intermediate
2512: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2513: @*/
2514: PetscErrorCode PetscDSHasBdJacobianPreconditioner(PetscDS ds, PetscBool *hasBdJacPre)
2515: {
2516: PetscFunctionBegin;
2518: PetscAssertPointer(hasBdJacPre, 2);
2519: PetscCall(PetscWeakFormHasBdJacobianPreconditioner(ds->wf, hasBdJacPre));
2520: PetscFunctionReturn(PETSC_SUCCESS);
2521: }
2523: /*@C
2524: PetscDSGetBdJacobianPreconditioner - Get the pointwise boundary Jacobian preconditioner function for given test and basis field
2526: Not Collective; No Fortran Support
2528: Input Parameters:
2529: + ds - The `PetscDS`
2530: . f - The test field number
2531: - g - The field number
2533: Output Parameters:
2534: + g0 - integrand for the test and basis function term
2535: . g1 - integrand for the test function and basis function gradient term
2536: . g2 - integrand for the test function gradient and basis function term
2537: - g3 - integrand for the test function gradient and basis function gradient term
2539: Calling sequence of `g0`:
2540: + dim - the spatial dimension
2541: . Nf - the number of fields
2542: . NfAux - the number of auxiliary fields
2543: . uOff - the offset into u[] and u_t[] for each field
2544: . uOff_x - the offset into u_x[] for each field
2545: . u - each field evaluated at the current point
2546: . u_t - the time derivative of each field evaluated at the current point
2547: . u_x - the gradient of each field evaluated at the current point
2548: . aOff - the offset into a[] and a_t[] for each auxiliary field
2549: . aOff_x - the offset into a_x[] for each auxiliary field
2550: . a - each auxiliary field evaluated at the current point
2551: . a_t - the time derivative of each auxiliary field evaluated at the current point
2552: . a_x - the gradient of auxiliary each field evaluated at the current point
2553: . t - current time
2554: . u_tShift - the multiplier a for dF/dU_t
2555: . x - coordinates of the current point
2556: . n - normal at the current point
2557: . numConstants - number of constant parameters
2558: . constants - constant parameters
2559: - g0 - output values at the current point
2561: Level: intermediate
2563: Note:
2564: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2566: We are using a first order FEM model for the weak form\:
2567: \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi
2569: .seealso: `PetscDS`, `PetscDSSetBdJacobianPreconditioner()`
2570: @*/
2571: 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[]))
2572: {
2573: PetscBdPointJac *tmp0, *tmp1, *tmp2, *tmp3;
2574: PetscInt n0, n1, n2, n3;
2576: PetscFunctionBegin;
2578: 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);
2579: 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);
2580: PetscCall(PetscWeakFormGetBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2581: *g0 = tmp0 ? tmp0[0] : NULL;
2582: *g1 = tmp1 ? tmp1[0] : NULL;
2583: *g2 = tmp2 ? tmp2[0] : NULL;
2584: *g3 = tmp3 ? tmp3[0] : NULL;
2585: PetscFunctionReturn(PETSC_SUCCESS);
2586: }
2588: /*@C
2589: PetscDSSetBdJacobianPreconditioner - Set the pointwise boundary Jacobian preconditioner function for given test and basis field
2591: Not Collective; No Fortran Support
2593: Input Parameters:
2594: + ds - The `PetscDS`
2595: . f - The test field number
2596: . g - The field number
2597: . g0 - integrand for the test and basis function term
2598: . g1 - integrand for the test function and basis function gradient term
2599: . g2 - integrand for the test function gradient and basis function term
2600: - g3 - integrand for the test function gradient and basis function gradient term
2602: Calling sequence of `g0':
2603: + dim - the spatial dimension
2604: . Nf - the number of fields
2605: . NfAux - the number of auxiliary fields
2606: . uOff - the offset into u[] and u_t[] for each field
2607: . uOff_x - the offset into u_x[] for each field
2608: . u - each field evaluated at the current point
2609: . u_t - the time derivative of each field evaluated at the current point
2610: . u_x - the gradient of each field evaluated at the current point
2611: . aOff - the offset into a[] and a_t[] for each auxiliary field
2612: . aOff_x - the offset into a_x[] for each auxiliary field
2613: . a - each auxiliary field evaluated at the current point
2614: . a_t - the time derivative of each auxiliary field evaluated at the current point
2615: . a_x - the gradient of auxiliary each field evaluated at the current point
2616: . t - current time
2617: . u_tShift - the multiplier a for dF/dU_t
2618: . x - coordinates of the current point
2619: . n - normal at the current point
2620: . numConstants - number of constant parameters
2621: . constants - constant parameters
2622: - g0 - output values at the current point
2624: Level: intermediate
2626: Note:
2627: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2629: We are using a first order FEM model for the weak form\:
2630: \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi
2632: .seealso: `PetscDS`, `PetscDSGetBdJacobianPreconditioner()`
2633: @*/
2634: 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[]))
2635: {
2636: PetscFunctionBegin;
2642: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2643: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2644: PetscCall(PetscWeakFormSetIndexBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2645: PetscFunctionReturn(PETSC_SUCCESS);
2646: }
2648: /*@C
2649: PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field
2651: Not Collective
2653: Input Parameters:
2654: + prob - The PetscDS
2655: - f - The test field number
2657: Output Parameters:
2658: + sol - exact solution for the test field
2659: - ctx - exact solution context
2661: Calling sequence of `exactSol`:
2662: + dim - the spatial dimension
2663: . t - current time
2664: . x - coordinates of the current point
2665: . Nc - the number of field components
2666: . u - the solution field evaluated at the current point
2667: - ctx - a user context
2669: Level: intermediate
2671: .seealso: `PetscDS`, `PetscDSSetExactSolution()`, `PetscDSGetExactSolutionTimeDerivative()`
2672: @*/
2673: PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx)
2674: {
2675: PetscFunctionBegin;
2677: 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);
2678: if (sol) {
2679: PetscAssertPointer(sol, 3);
2680: *sol = prob->exactSol[f];
2681: }
2682: if (ctx) {
2683: PetscAssertPointer(ctx, 4);
2684: *ctx = prob->exactCtx[f];
2685: }
2686: PetscFunctionReturn(PETSC_SUCCESS);
2687: }
2689: /*@C
2690: PetscDSSetExactSolution - Set the pointwise exact solution function for a given test field
2692: Not Collective
2694: Input Parameters:
2695: + prob - The `PetscDS`
2696: . f - The test field number
2697: . sol - solution function for the test fields
2698: - ctx - solution context or `NULL`
2700: Calling sequence of `sol`:
2701: + dim - the spatial dimension
2702: . t - current time
2703: . x - coordinates of the current point
2704: . Nc - the number of field components
2705: . u - the solution field evaluated at the current point
2706: - ctx - a user context
2708: Level: intermediate
2710: .seealso: `PetscDS`, `PetscDSGetExactSolution()`
2711: @*/
2712: PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx)
2713: {
2714: PetscFunctionBegin;
2716: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2717: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
2718: if (sol) {
2720: prob->exactSol[f] = sol;
2721: }
2722: if (ctx) {
2724: prob->exactCtx[f] = ctx;
2725: }
2726: PetscFunctionReturn(PETSC_SUCCESS);
2727: }
2729: /*@C
2730: PetscDSGetExactSolutionTimeDerivative - Get the pointwise time derivative of the exact solution function for a given test field
2732: Not Collective
2734: Input Parameters:
2735: + prob - The `PetscDS`
2736: - f - The test field number
2738: Output Parameters:
2739: + sol - time derivative of the exact solution for the test field
2740: - ctx - time derivative of the exact solution context
2742: Calling sequence of `exactSol`:
2743: + dim - the spatial dimension
2744: . t - current time
2745: . x - coordinates of the current point
2746: . Nc - the number of field components
2747: . u - the solution field evaluated at the current point
2748: - ctx - a user context
2750: Level: intermediate
2752: .seealso: `PetscDS`, `PetscDSSetExactSolutionTimeDerivative()`, `PetscDSGetExactSolution()`
2753: @*/
2754: PetscErrorCode PetscDSGetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx)
2755: {
2756: PetscFunctionBegin;
2758: 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);
2759: if (sol) {
2760: PetscAssertPointer(sol, 3);
2761: *sol = prob->exactSol_t[f];
2762: }
2763: if (ctx) {
2764: PetscAssertPointer(ctx, 4);
2765: *ctx = prob->exactCtx_t[f];
2766: }
2767: PetscFunctionReturn(PETSC_SUCCESS);
2768: }
2770: /*@C
2771: PetscDSSetExactSolutionTimeDerivative - Set the pointwise time derivative of the exact solution function for a given test field
2773: Not Collective
2775: Input Parameters:
2776: + prob - The `PetscDS`
2777: . f - The test field number
2778: . sol - time derivative of the solution function for the test fields
2779: - ctx - time derivative of the solution context or `NULL`
2781: Calling sequence of `sol`:
2782: + dim - the spatial dimension
2783: . t - current time
2784: . x - coordinates of the current point
2785: . Nc - the number of field components
2786: . u - the solution field evaluated at the current point
2787: - ctx - a user context
2789: Level: intermediate
2791: .seealso: `PetscDS`, `PetscDSGetExactSolutionTimeDerivative()`, `PetscDSSetExactSolution()`
2792: @*/
2793: PetscErrorCode PetscDSSetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx)
2794: {
2795: PetscFunctionBegin;
2797: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2798: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
2799: if (sol) {
2801: prob->exactSol_t[f] = sol;
2802: }
2803: if (ctx) {
2805: prob->exactCtx_t[f] = ctx;
2806: }
2807: PetscFunctionReturn(PETSC_SUCCESS);
2808: }
2810: /*@C
2811: PetscDSGetConstants - Returns the array of constants passed to point functions
2813: Not Collective
2815: Input Parameter:
2816: . prob - The `PetscDS` object
2818: Output Parameters:
2819: + numConstants - The number of constants
2820: - constants - The array of constants, NULL if there are none
2822: Level: intermediate
2824: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSCreate()`
2825: @*/
2826: PetscErrorCode PetscDSGetConstants(PetscDS prob, PetscInt *numConstants, const PetscScalar *constants[])
2827: {
2828: PetscFunctionBegin;
2830: if (numConstants) {
2831: PetscAssertPointer(numConstants, 2);
2832: *numConstants = prob->numConstants;
2833: }
2834: if (constants) {
2835: PetscAssertPointer(constants, 3);
2836: *constants = prob->constants;
2837: }
2838: PetscFunctionReturn(PETSC_SUCCESS);
2839: }
2841: /*@C
2842: PetscDSSetConstants - Set the array of constants passed to point functions
2844: Not Collective
2846: Input Parameters:
2847: + prob - The `PetscDS` object
2848: . numConstants - The number of constants
2849: - constants - The array of constants, `NULL` if there are none
2851: Level: intermediate
2853: .seealso: `PetscDS`, `PetscDSGetConstants()`, `PetscDSCreate()`
2854: @*/
2855: PetscErrorCode PetscDSSetConstants(PetscDS prob, PetscInt numConstants, PetscScalar constants[])
2856: {
2857: PetscFunctionBegin;
2859: if (numConstants != prob->numConstants) {
2860: PetscCall(PetscFree(prob->constants));
2861: prob->numConstants = numConstants;
2862: if (prob->numConstants) {
2863: PetscCall(PetscMalloc1(prob->numConstants, &prob->constants));
2864: } else {
2865: prob->constants = NULL;
2866: }
2867: }
2868: if (prob->numConstants) {
2869: PetscAssertPointer(constants, 3);
2870: PetscCall(PetscArraycpy(prob->constants, constants, prob->numConstants));
2871: }
2872: PetscFunctionReturn(PETSC_SUCCESS);
2873: }
2875: /*@
2876: PetscDSGetFieldIndex - Returns the index of the given field
2878: Not Collective
2880: Input Parameters:
2881: + prob - The `PetscDS` object
2882: - disc - The discretization object
2884: Output Parameter:
2885: . f - The field number
2887: Level: beginner
2889: .seealso: `PetscDS`, `PetscGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2890: @*/
2891: PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f)
2892: {
2893: PetscInt g;
2895: PetscFunctionBegin;
2897: PetscAssertPointer(f, 3);
2898: *f = -1;
2899: for (g = 0; g < prob->Nf; ++g) {
2900: if (disc == prob->disc[g]) break;
2901: }
2902: PetscCheck(g != prob->Nf, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS.");
2903: *f = g;
2904: PetscFunctionReturn(PETSC_SUCCESS);
2905: }
2907: /*@
2908: PetscDSGetFieldSize - Returns the size of the given field in the full space basis
2910: Not Collective
2912: Input Parameters:
2913: + prob - The `PetscDS` object
2914: - f - The field number
2916: Output Parameter:
2917: . size - The size
2919: Level: beginner
2921: .seealso: `PetscDS`, `PetscDSGetFieldOffset()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2922: @*/
2923: PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size)
2924: {
2925: PetscFunctionBegin;
2927: PetscAssertPointer(size, 3);
2928: 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);
2929: PetscCall(PetscDSSetUp(prob));
2930: *size = prob->Nb[f];
2931: PetscFunctionReturn(PETSC_SUCCESS);
2932: }
2934: /*@
2935: PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis
2937: Not Collective
2939: Input Parameters:
2940: + prob - The `PetscDS` object
2941: - f - The field number
2943: Output Parameter:
2944: . off - The offset
2946: Level: beginner
2948: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2949: @*/
2950: PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off)
2951: {
2952: PetscInt size, g;
2954: PetscFunctionBegin;
2956: PetscAssertPointer(off, 3);
2957: 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);
2958: *off = 0;
2959: for (g = 0; g < f; ++g) {
2960: PetscCall(PetscDSGetFieldSize(prob, g, &size));
2961: *off += size;
2962: }
2963: PetscFunctionReturn(PETSC_SUCCESS);
2964: }
2966: /*@
2967: PetscDSGetFieldOffsetCohesive - Returns the offset of the given field in the full space basis on a cohesive cell
2969: Not Collective
2971: Input Parameters:
2972: + ds - The `PetscDS` object
2973: - f - The field number
2975: Output Parameter:
2976: . off - The offset
2978: Level: beginner
2980: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2981: @*/
2982: PetscErrorCode PetscDSGetFieldOffsetCohesive(PetscDS ds, PetscInt f, PetscInt *off)
2983: {
2984: PetscInt size, g;
2986: PetscFunctionBegin;
2988: PetscAssertPointer(off, 3);
2989: 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);
2990: *off = 0;
2991: for (g = 0; g < f; ++g) {
2992: PetscBool cohesive;
2994: PetscCall(PetscDSGetCohesive(ds, g, &cohesive));
2995: PetscCall(PetscDSGetFieldSize(ds, g, &size));
2996: *off += cohesive ? size : size * 2;
2997: }
2998: PetscFunctionReturn(PETSC_SUCCESS);
2999: }
3001: /*@
3002: PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point
3004: Not Collective
3006: Input Parameter:
3007: . prob - The `PetscDS` object
3009: Output Parameter:
3010: . dimensions - The number of dimensions
3012: Level: beginner
3014: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3015: @*/
3016: PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[])
3017: {
3018: PetscFunctionBegin;
3020: PetscCall(PetscDSSetUp(prob));
3021: PetscAssertPointer(dimensions, 2);
3022: *dimensions = prob->Nb;
3023: PetscFunctionReturn(PETSC_SUCCESS);
3024: }
3026: /*@
3027: PetscDSGetComponents - Returns the number of components for each field on an evaluation point
3029: Not Collective
3031: Input Parameter:
3032: . prob - The `PetscDS` object
3034: Output Parameter:
3035: . components - The number of components
3037: Level: beginner
3039: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3040: @*/
3041: PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[])
3042: {
3043: PetscFunctionBegin;
3045: PetscCall(PetscDSSetUp(prob));
3046: PetscAssertPointer(components, 2);
3047: *components = prob->Nc;
3048: PetscFunctionReturn(PETSC_SUCCESS);
3049: }
3051: /*@
3052: PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point
3054: Not Collective
3056: Input Parameters:
3057: + prob - The `PetscDS` object
3058: - f - The field number
3060: Output Parameter:
3061: . off - The offset
3063: Level: beginner
3065: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3066: @*/
3067: PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off)
3068: {
3069: PetscFunctionBegin;
3071: PetscAssertPointer(off, 3);
3072: 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);
3073: PetscCall(PetscDSSetUp(prob));
3074: *off = prob->off[f];
3075: PetscFunctionReturn(PETSC_SUCCESS);
3076: }
3078: /*@
3079: PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point
3081: Not Collective
3083: Input Parameter:
3084: . prob - The `PetscDS` object
3086: Output Parameter:
3087: . offsets - The offsets
3089: Level: beginner
3091: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3092: @*/
3093: PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[])
3094: {
3095: PetscFunctionBegin;
3097: PetscAssertPointer(offsets, 2);
3098: PetscCall(PetscDSSetUp(prob));
3099: *offsets = prob->off;
3100: PetscFunctionReturn(PETSC_SUCCESS);
3101: }
3103: /*@
3104: PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point
3106: Not Collective
3108: Input Parameter:
3109: . prob - The `PetscDS` object
3111: Output Parameter:
3112: . offsets - The offsets
3114: Level: beginner
3116: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3117: @*/
3118: PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[])
3119: {
3120: PetscFunctionBegin;
3122: PetscAssertPointer(offsets, 2);
3123: PetscCall(PetscDSSetUp(prob));
3124: *offsets = prob->offDer;
3125: PetscFunctionReturn(PETSC_SUCCESS);
3126: }
3128: /*@
3129: PetscDSGetComponentOffsetsCohesive - Returns the offset of each field on an evaluation point
3131: Not Collective
3133: Input Parameters:
3134: + ds - The `PetscDS` object
3135: - s - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive
3137: Output Parameter:
3138: . offsets - The offsets
3140: Level: beginner
3142: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3143: @*/
3144: PetscErrorCode PetscDSGetComponentOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
3145: {
3146: PetscFunctionBegin;
3148: PetscAssertPointer(offsets, 3);
3149: PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS");
3150: PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s);
3151: PetscCall(PetscDSSetUp(ds));
3152: *offsets = ds->offCohesive[s];
3153: PetscFunctionReturn(PETSC_SUCCESS);
3154: }
3156: /*@
3157: PetscDSGetComponentDerivativeOffsetsCohesive - Returns the offset of each field derivative on an evaluation point
3159: Not Collective
3161: Input Parameters:
3162: + ds - The `PetscDS` object
3163: - s - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive
3165: Output Parameter:
3166: . offsets - The offsets
3168: Level: beginner
3170: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3171: @*/
3172: PetscErrorCode PetscDSGetComponentDerivativeOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
3173: {
3174: PetscFunctionBegin;
3176: PetscAssertPointer(offsets, 3);
3177: PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS");
3178: PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s);
3179: PetscCall(PetscDSSetUp(ds));
3180: *offsets = ds->offDerCohesive[s];
3181: PetscFunctionReturn(PETSC_SUCCESS);
3182: }
3184: /*@C
3185: PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization
3187: Not Collective
3189: Input Parameter:
3190: . prob - The `PetscDS` object
3192: Output Parameter:
3193: . T - The basis function and derivatives tabulation at quadrature points for each field
3195: Level: intermediate
3197: .seealso: `PetscDS`, `PetscTabulation`, `PetscDSCreate()`
3198: @*/
3199: PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscTabulation *T[])
3200: {
3201: PetscFunctionBegin;
3203: PetscAssertPointer(T, 2);
3204: PetscCall(PetscDSSetUp(prob));
3205: *T = prob->T;
3206: PetscFunctionReturn(PETSC_SUCCESS);
3207: }
3209: /*@C
3210: PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces
3212: Not Collective
3214: Input Parameter:
3215: . prob - The `PetscDS` object
3217: Output Parameter:
3218: . Tf - The basis function and derivative tabulation on each local face at quadrature points for each and field
3220: Level: intermediate
3222: .seealso: `PetscTabulation`, `PetscDS`, `PetscDSGetTabulation()`, `PetscDSCreate()`
3223: @*/
3224: PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscTabulation *Tf[])
3225: {
3226: PetscFunctionBegin;
3228: PetscAssertPointer(Tf, 2);
3229: PetscCall(PetscDSSetUp(prob));
3230: *Tf = prob->Tf;
3231: PetscFunctionReturn(PETSC_SUCCESS);
3232: }
3234: PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x)
3235: {
3236: PetscFunctionBegin;
3238: PetscCall(PetscDSSetUp(prob));
3239: if (u) {
3240: PetscAssertPointer(u, 2);
3241: *u = prob->u;
3242: }
3243: if (u_t) {
3244: PetscAssertPointer(u_t, 3);
3245: *u_t = prob->u_t;
3246: }
3247: if (u_x) {
3248: PetscAssertPointer(u_x, 4);
3249: *u_x = prob->u_x;
3250: }
3251: PetscFunctionReturn(PETSC_SUCCESS);
3252: }
3254: PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3)
3255: {
3256: PetscFunctionBegin;
3258: PetscCall(PetscDSSetUp(prob));
3259: if (f0) {
3260: PetscAssertPointer(f0, 2);
3261: *f0 = prob->f0;
3262: }
3263: if (f1) {
3264: PetscAssertPointer(f1, 3);
3265: *f1 = prob->f1;
3266: }
3267: if (g0) {
3268: PetscAssertPointer(g0, 4);
3269: *g0 = prob->g0;
3270: }
3271: if (g1) {
3272: PetscAssertPointer(g1, 5);
3273: *g1 = prob->g1;
3274: }
3275: if (g2) {
3276: PetscAssertPointer(g2, 6);
3277: *g2 = prob->g2;
3278: }
3279: if (g3) {
3280: PetscAssertPointer(g3, 7);
3281: *g3 = prob->g3;
3282: }
3283: PetscFunctionReturn(PETSC_SUCCESS);
3284: }
3286: PetscErrorCode PetscDSGetWorkspace(PetscDS prob, PetscReal **x, PetscScalar **basisReal, PetscScalar **basisDerReal, PetscScalar **testReal, PetscScalar **testDerReal)
3287: {
3288: PetscFunctionBegin;
3290: PetscCall(PetscDSSetUp(prob));
3291: if (x) {
3292: PetscAssertPointer(x, 2);
3293: *x = prob->x;
3294: }
3295: if (basisReal) {
3296: PetscAssertPointer(basisReal, 3);
3297: *basisReal = prob->basisReal;
3298: }
3299: if (basisDerReal) {
3300: PetscAssertPointer(basisDerReal, 4);
3301: *basisDerReal = prob->basisDerReal;
3302: }
3303: if (testReal) {
3304: PetscAssertPointer(testReal, 5);
3305: *testReal = prob->testReal;
3306: }
3307: if (testDerReal) {
3308: PetscAssertPointer(testDerReal, 6);
3309: *testDerReal = prob->testDerReal;
3310: }
3311: PetscFunctionReturn(PETSC_SUCCESS);
3312: }
3314: /*@C
3315: PetscDSAddBoundary - Add a boundary condition to the model.
3317: Collective
3319: Input Parameters:
3320: + ds - The PetscDS object
3321: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3322: . name - The BC name
3323: . label - The label defining constrained points
3324: . Nv - The number of `DMLabel` values for constrained points
3325: . values - An array of label values for constrained points
3326: . field - The field to constrain
3327: . Nc - The number of constrained field components (0 will constrain all fields)
3328: . comps - An array of constrained component numbers
3329: . bcFunc - A pointwise function giving boundary values
3330: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3331: - ctx - An optional user context for bcFunc
3333: Output Parameter:
3334: . bd - The boundary number
3336: Options Database Keys:
3337: + -bc_<boundary name> <num> - Overrides the boundary ids
3338: - -bc_<boundary name>_comp <num> - Overrides the boundary components
3340: Level: developer
3342: Note:
3343: Both `bcFunc` and `bcFunc_t` will depend on the boundary condition type. If the type if `DM_BC_ESSENTIAL`, then the calling sequence is\:
3345: $ void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])
3347: If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value, then the calling sequence is\:
3348: .vb
3349: void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3350: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3351: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3352: PetscReal time, const PetscReal x[], PetscScalar bcval[])
3353: .ve
3354: + dim - the spatial dimension
3355: . Nf - the number of fields
3356: . uOff - the offset into u[] and u_t[] for each field
3357: . uOff_x - the offset into u_x[] for each field
3358: . u - each field evaluated at the current point
3359: . u_t - the time derivative of each field evaluated at the current point
3360: . u_x - the gradient of each field evaluated at the current point
3361: . aOff - the offset into a[] and a_t[] for each auxiliary field
3362: . aOff_x - the offset into a_x[] for each auxiliary field
3363: . a - each auxiliary field evaluated at the current point
3364: . a_t - the time derivative of each auxiliary field evaluated at the current point
3365: . a_x - the gradient of auxiliary each field evaluated at the current point
3366: . t - current time
3367: . x - coordinates of the current point
3368: . numConstants - number of constant parameters
3369: . constants - constant parameters
3370: - bcval - output values at the current point
3372: Notes:
3373: The pointwise functions are used to provide boundary values for essential boundary
3374: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3375: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3376: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3378: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundaryByName()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3379: @*/
3380: 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)
3381: {
3382: DSBoundary head = ds->boundary, b;
3383: PetscInt n = 0;
3384: const char *lname;
3386: PetscFunctionBegin;
3389: PetscAssertPointer(name, 3);
3394: 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);
3395: if (Nc > 0) {
3396: PetscInt *fcomps;
3397: PetscInt c;
3399: PetscCall(PetscDSGetComponents(ds, &fcomps));
3400: 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);
3401: for (c = 0; c < Nc; ++c) {
3402: 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);
3403: }
3404: }
3405: PetscCall(PetscNew(&b));
3406: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3407: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf));
3408: PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf));
3409: PetscCall(PetscMalloc1(Nv, &b->values));
3410: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3411: PetscCall(PetscMalloc1(Nc, &b->comps));
3412: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3413: PetscCall(PetscObjectGetName((PetscObject)label, &lname));
3414: PetscCall(PetscStrallocpy(lname, (char **)&b->lname));
3415: b->type = type;
3416: b->label = label;
3417: b->Nv = Nv;
3418: b->field = field;
3419: b->Nc = Nc;
3420: b->func = bcFunc;
3421: b->func_t = bcFunc_t;
3422: b->ctx = ctx;
3423: b->next = NULL;
3424: /* Append to linked list so that we can preserve the order */
3425: if (!head) ds->boundary = b;
3426: while (head) {
3427: if (!head->next) {
3428: head->next = b;
3429: head = b;
3430: }
3431: head = head->next;
3432: ++n;
3433: }
3434: if (bd) {
3435: PetscAssertPointer(bd, 13);
3436: *bd = n;
3437: }
3438: PetscFunctionReturn(PETSC_SUCCESS);
3439: }
3441: // PetscClangLinter pragma ignore: -fdoc-section-header-unknown
3442: /*@C
3443: PetscDSAddBoundaryByName - Add a boundary condition to the model.
3445: Collective
3447: Input Parameters:
3448: + ds - The `PetscDS` object
3449: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3450: . name - The BC name
3451: . lname - The naem of the label defining constrained points
3452: . Nv - The number of `DMLabel` values for constrained points
3453: . values - An array of label values for constrained points
3454: . field - The field to constrain
3455: . Nc - The number of constrained field components (0 will constrain all fields)
3456: . comps - An array of constrained component numbers
3457: . bcFunc - A pointwise function giving boundary values
3458: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3459: - ctx - An optional user context for bcFunc
3461: Output Parameter:
3462: . bd - The boundary number
3464: Options Database Keys:
3465: + -bc_<boundary name> <num> - Overrides the boundary ids
3466: - -bc_<boundary name>_comp <num> - Overrides the boundary components
3468: Calling Sequence of `bcFunc` and `bcFunc_t`:
3469: If the type is `DM_BC_ESSENTIAL`
3470: .vb
3471: void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])
3472: .ve
3473: If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value,
3474: .vb
3475: void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3476: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3477: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3478: PetscReal time, const PetscReal x[], PetscScalar bcval[])
3479: .ve
3480: + dim - the spatial dimension
3481: . Nf - the number of fields
3482: . uOff - the offset into u[] and u_t[] for each field
3483: . uOff_x - the offset into u_x[] for each field
3484: . u - each field evaluated at the current point
3485: . u_t - the time derivative of each field evaluated at the current point
3486: . u_x - the gradient of each field evaluated at the current point
3487: . aOff - the offset into a[] and a_t[] for each auxiliary field
3488: . aOff_x - the offset into a_x[] for each auxiliary field
3489: . a - each auxiliary field evaluated at the current point
3490: . a_t - the time derivative of each auxiliary field evaluated at the current point
3491: . a_x - the gradient of auxiliary each field evaluated at the current point
3492: . t - current time
3493: . x - coordinates of the current point
3494: . numConstants - number of constant parameters
3495: . constants - constant parameters
3496: - bcval - output values at the current point
3498: Level: developer
3500: Notes:
3501: The pointwise functions are used to provide boundary values for essential boundary
3502: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3503: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3504: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3506: This function should only be used with `DMFOREST` currently, since labels cannot be defined before the underlying `DMPLEX` is built.
3508: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3509: @*/
3510: 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)
3511: {
3512: DSBoundary head = ds->boundary, b;
3513: PetscInt n = 0;
3515: PetscFunctionBegin;
3518: PetscAssertPointer(name, 3);
3519: PetscAssertPointer(lname, 4);
3523: PetscCall(PetscNew(&b));
3524: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3525: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf));
3526: PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf));
3527: PetscCall(PetscMalloc1(Nv, &b->values));
3528: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3529: PetscCall(PetscMalloc1(Nc, &b->comps));
3530: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3531: PetscCall(PetscStrallocpy(lname, (char **)&b->lname));
3532: b->type = type;
3533: b->label = NULL;
3534: b->Nv = Nv;
3535: b->field = field;
3536: b->Nc = Nc;
3537: b->func = bcFunc;
3538: b->func_t = bcFunc_t;
3539: b->ctx = ctx;
3540: b->next = NULL;
3541: /* Append to linked list so that we can preserve the order */
3542: if (!head) ds->boundary = b;
3543: while (head) {
3544: if (!head->next) {
3545: head->next = b;
3546: head = b;
3547: }
3548: head = head->next;
3549: ++n;
3550: }
3551: if (bd) {
3552: PetscAssertPointer(bd, 13);
3553: *bd = n;
3554: }
3555: PetscFunctionReturn(PETSC_SUCCESS);
3556: }
3558: /*@C
3559: PetscDSUpdateBoundary - Change a boundary condition for the model.
3561: Input Parameters:
3562: + ds - The `PetscDS` object
3563: . bd - The boundary condition number
3564: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3565: . name - The BC name
3566: . label - The label defining constrained points
3567: . Nv - The number of `DMLabel` ids for constrained points
3568: . values - An array of ids for constrained points
3569: . field - The field to constrain
3570: . Nc - The number of constrained field components
3571: . comps - An array of constrained component numbers
3572: . bcFunc - A pointwise function giving boundary values
3573: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3574: - ctx - An optional user context for bcFunc
3576: Level: developer
3578: Notes:
3579: The pointwise functions are used to provide boundary values for essential boundary
3580: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3581: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3582: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3584: The boundary condition number is the order in which it was registered. The user can get the number of boundary conditions from `PetscDSGetNumBoundary()`.
3585: See `PetscDSAddBoundary()` for a description of the calling sequences for the callbacks.
3587: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSGetNumBoundary()`, `DMLabel`
3588: @*/
3589: 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)
3590: {
3591: DSBoundary b = ds->boundary;
3592: PetscInt n = 0;
3594: PetscFunctionBegin;
3596: while (b) {
3597: if (n == bd) break;
3598: b = b->next;
3599: ++n;
3600: }
3601: PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n);
3602: if (name) {
3603: PetscCall(PetscFree(b->name));
3604: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3605: }
3606: b->type = type;
3607: if (label) {
3608: const char *name;
3610: b->label = label;
3611: PetscCall(PetscFree(b->lname));
3612: PetscCall(PetscObjectGetName((PetscObject)label, &name));
3613: PetscCall(PetscStrallocpy(name, (char **)&b->lname));
3614: }
3615: if (Nv >= 0) {
3616: b->Nv = Nv;
3617: PetscCall(PetscFree(b->values));
3618: PetscCall(PetscMalloc1(Nv, &b->values));
3619: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3620: }
3621: if (field >= 0) b->field = field;
3622: if (Nc >= 0) {
3623: b->Nc = Nc;
3624: PetscCall(PetscFree(b->comps));
3625: PetscCall(PetscMalloc1(Nc, &b->comps));
3626: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3627: }
3628: if (bcFunc) b->func = bcFunc;
3629: if (bcFunc_t) b->func_t = bcFunc_t;
3630: if (ctx) b->ctx = ctx;
3631: PetscFunctionReturn(PETSC_SUCCESS);
3632: }
3634: /*@
3635: PetscDSGetNumBoundary - Get the number of registered BC
3637: Input Parameter:
3638: . ds - The `PetscDS` object
3640: Output Parameter:
3641: . numBd - The number of BC
3643: Level: intermediate
3645: .seealso: `PetscDS`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`
3646: @*/
3647: PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd)
3648: {
3649: DSBoundary b = ds->boundary;
3651: PetscFunctionBegin;
3653: PetscAssertPointer(numBd, 2);
3654: *numBd = 0;
3655: while (b) {
3656: ++(*numBd);
3657: b = b->next;
3658: }
3659: PetscFunctionReturn(PETSC_SUCCESS);
3660: }
3662: /*@C
3663: PetscDSGetBoundary - Gets a boundary condition to the model
3665: Input Parameters:
3666: + ds - The `PetscDS` object
3667: - bd - The BC number
3669: Output Parameters:
3670: + wf - The `PetscWeakForm` holding the pointwise functions
3671: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3672: . name - The BC name
3673: . label - The label defining constrained points
3674: . Nv - The number of `DMLabel` ids for constrained points
3675: . values - An array of ids for constrained points
3676: . field - The field to constrain
3677: . Nc - The number of constrained field components
3678: . comps - An array of constrained component numbers
3679: . func - A pointwise function giving boundary values
3680: . func_t - A pointwise function giving the time derivative of the boundary values
3681: - ctx - An optional user context for bcFunc
3683: Options Database Keys:
3684: + -bc_<boundary name> <num> - Overrides the boundary ids
3685: - -bc_<boundary name>_comp <num> - Overrides the boundary components
3687: Level: developer
3689: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `DMLabel`
3690: @*/
3691: 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)
3692: {
3693: DSBoundary b = ds->boundary;
3694: PetscInt n = 0;
3696: PetscFunctionBegin;
3698: while (b) {
3699: if (n == bd) break;
3700: b = b->next;
3701: ++n;
3702: }
3703: PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n);
3704: if (wf) {
3705: PetscAssertPointer(wf, 3);
3706: *wf = b->wf;
3707: }
3708: if (type) {
3709: PetscAssertPointer(type, 4);
3710: *type = b->type;
3711: }
3712: if (name) {
3713: PetscAssertPointer(name, 5);
3714: *name = b->name;
3715: }
3716: if (label) {
3717: PetscAssertPointer(label, 6);
3718: *label = b->label;
3719: }
3720: if (Nv) {
3721: PetscAssertPointer(Nv, 7);
3722: *Nv = b->Nv;
3723: }
3724: if (values) {
3725: PetscAssertPointer(values, 8);
3726: *values = b->values;
3727: }
3728: if (field) {
3729: PetscAssertPointer(field, 9);
3730: *field = b->field;
3731: }
3732: if (Nc) {
3733: PetscAssertPointer(Nc, 10);
3734: *Nc = b->Nc;
3735: }
3736: if (comps) {
3737: PetscAssertPointer(comps, 11);
3738: *comps = b->comps;
3739: }
3740: if (func) {
3741: PetscAssertPointer(func, 12);
3742: *func = b->func;
3743: }
3744: if (func_t) {
3745: PetscAssertPointer(func_t, 13);
3746: *func_t = b->func_t;
3747: }
3748: if (ctx) {
3749: PetscAssertPointer(ctx, 14);
3750: *ctx = b->ctx;
3751: }
3752: PetscFunctionReturn(PETSC_SUCCESS);
3753: }
3755: /*@
3756: PetscDSUpdateBoundaryLabels - Update `DMLabel` in each boundary condition using the label name and the input `DM`
3758: Not Collective
3760: Input Parameters:
3761: + ds - The source `PetscDS` object
3762: - dm - The `DM` holding labels
3764: Level: intermediate
3766: .seealso: `PetscDS`, `DMBoundary`, `DM`, `PetscDSCopyBoundary()`, `PetscDSCreate()`, `DMGetLabel()`
3767: @*/
3768: PetscErrorCode PetscDSUpdateBoundaryLabels(PetscDS ds, DM dm)
3769: {
3770: DSBoundary b;
3772: PetscFunctionBegin;
3775: for (b = ds->boundary; b; b = b->next) {
3776: if (b->lname) PetscCall(DMGetLabel(dm, b->lname, &b->label));
3777: }
3778: PetscFunctionReturn(PETSC_SUCCESS);
3779: }
3781: static PetscErrorCode DSBoundaryDuplicate_Internal(DSBoundary b, DSBoundary *bNew)
3782: {
3783: PetscFunctionBegin;
3784: PetscCall(PetscNew(bNew));
3785: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &(*bNew)->wf));
3786: PetscCall(PetscWeakFormCopy(b->wf, (*bNew)->wf));
3787: PetscCall(PetscStrallocpy(b->name, (char **)&((*bNew)->name)));
3788: PetscCall(PetscStrallocpy(b->lname, (char **)&((*bNew)->lname)));
3789: (*bNew)->type = b->type;
3790: (*bNew)->label = b->label;
3791: (*bNew)->Nv = b->Nv;
3792: PetscCall(PetscMalloc1(b->Nv, &(*bNew)->values));
3793: PetscCall(PetscArraycpy((*bNew)->values, b->values, b->Nv));
3794: (*bNew)->field = b->field;
3795: (*bNew)->Nc = b->Nc;
3796: PetscCall(PetscMalloc1(b->Nc, &(*bNew)->comps));
3797: PetscCall(PetscArraycpy((*bNew)->comps, b->comps, b->Nc));
3798: (*bNew)->func = b->func;
3799: (*bNew)->func_t = b->func_t;
3800: (*bNew)->ctx = b->ctx;
3801: PetscFunctionReturn(PETSC_SUCCESS);
3802: }
3804: /*@
3805: PetscDSCopyBoundary - Copy all boundary condition objects to the new problem
3807: Not Collective
3809: Input Parameters:
3810: + ds - The source `PetscDS` object
3811: . numFields - The number of selected fields, or `PETSC_DEFAULT` for all fields
3812: - fields - The selected fields, or NULL for all fields
3814: Output Parameter:
3815: . newds - The target `PetscDS`, now with a copy of the boundary conditions
3817: Level: intermediate
3819: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3820: @*/
3821: PetscErrorCode PetscDSCopyBoundary(PetscDS ds, PetscInt numFields, const PetscInt fields[], PetscDS newds)
3822: {
3823: DSBoundary b, *lastnext;
3825: PetscFunctionBegin;
3828: if (ds == newds) PetscFunctionReturn(PETSC_SUCCESS);
3829: PetscCall(PetscDSDestroyBoundary(newds));
3830: lastnext = &newds->boundary;
3831: for (b = ds->boundary; b; b = b->next) {
3832: DSBoundary bNew;
3833: PetscInt fieldNew = -1;
3835: if (numFields > 0 && fields) {
3836: PetscInt f;
3838: for (f = 0; f < numFields; ++f)
3839: if (b->field == fields[f]) break;
3840: if (f == numFields) continue;
3841: fieldNew = f;
3842: }
3843: PetscCall(DSBoundaryDuplicate_Internal(b, &bNew));
3844: bNew->field = fieldNew < 0 ? b->field : fieldNew;
3845: *lastnext = bNew;
3846: lastnext = &bNew->next;
3847: }
3848: PetscFunctionReturn(PETSC_SUCCESS);
3849: }
3851: /*@
3852: PetscDSDestroyBoundary - Remove all `DMBoundary` objects from the `PetscDS`
3854: Not Collective
3856: Input Parameter:
3857: . ds - The `PetscDS` object
3859: Level: intermediate
3861: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`
3862: @*/
3863: PetscErrorCode PetscDSDestroyBoundary(PetscDS ds)
3864: {
3865: DSBoundary next = ds->boundary;
3867: PetscFunctionBegin;
3868: while (next) {
3869: DSBoundary b = next;
3871: next = b->next;
3872: PetscCall(PetscWeakFormDestroy(&b->wf));
3873: PetscCall(PetscFree(b->name));
3874: PetscCall(PetscFree(b->lname));
3875: PetscCall(PetscFree(b->values));
3876: PetscCall(PetscFree(b->comps));
3877: PetscCall(PetscFree(b));
3878: }
3879: PetscFunctionReturn(PETSC_SUCCESS);
3880: }
3882: /*@
3883: PetscDSSelectDiscretizations - Copy discretizations to the new problem with different field layout
3885: Not Collective
3887: Input Parameters:
3888: + prob - The `PetscDS` object
3889: . numFields - Number of new fields
3890: - fields - Old field number for each new field
3892: Output Parameter:
3893: . newprob - The `PetscDS` copy
3895: Level: intermediate
3897: .seealso: `PetscDS`, `PetscDSSelectEquations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3898: @*/
3899: PetscErrorCode PetscDSSelectDiscretizations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
3900: {
3901: PetscInt Nf, Nfn, fn;
3903: PetscFunctionBegin;
3905: if (fields) PetscAssertPointer(fields, 3);
3907: PetscCall(PetscDSGetNumFields(prob, &Nf));
3908: PetscCall(PetscDSGetNumFields(newprob, &Nfn));
3909: numFields = numFields < 0 ? Nf : numFields;
3910: for (fn = 0; fn < numFields; ++fn) {
3911: const PetscInt f = fields ? fields[fn] : fn;
3912: PetscObject disc;
3914: if (f >= Nf) continue;
3915: PetscCall(PetscDSGetDiscretization(prob, f, &disc));
3916: PetscCall(PetscDSSetDiscretization(newprob, fn, disc));
3917: }
3918: PetscFunctionReturn(PETSC_SUCCESS);
3919: }
3921: /*@
3922: PetscDSSelectEquations - Copy pointwise function pointers to the new problem with different field layout
3924: Not Collective
3926: Input Parameters:
3927: + prob - The `PetscDS` object
3928: . numFields - Number of new fields
3929: - fields - Old field number for each new field
3931: Output Parameter:
3932: . newprob - The `PetscDS` copy
3934: Level: intermediate
3936: .seealso: `PetscDS`, `PetscDSSelectDiscretizations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3937: @*/
3938: PetscErrorCode PetscDSSelectEquations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
3939: {
3940: PetscInt Nf, Nfn, fn, gn;
3942: PetscFunctionBegin;
3944: if (fields) PetscAssertPointer(fields, 3);
3946: PetscCall(PetscDSGetNumFields(prob, &Nf));
3947: PetscCall(PetscDSGetNumFields(newprob, &Nfn));
3948: 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);
3949: for (fn = 0; fn < numFields; ++fn) {
3950: const PetscInt f = fields ? fields[fn] : fn;
3951: PetscPointFunc obj;
3952: PetscPointFunc f0, f1;
3953: PetscBdPointFunc f0Bd, f1Bd;
3954: PetscRiemannFunc r;
3956: if (f >= Nf) continue;
3957: PetscCall(PetscDSGetObjective(prob, f, &obj));
3958: PetscCall(PetscDSGetResidual(prob, f, &f0, &f1));
3959: PetscCall(PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd));
3960: PetscCall(PetscDSGetRiemannSolver(prob, f, &r));
3961: PetscCall(PetscDSSetObjective(newprob, fn, obj));
3962: PetscCall(PetscDSSetResidual(newprob, fn, f0, f1));
3963: PetscCall(PetscDSSetBdResidual(newprob, fn, f0Bd, f1Bd));
3964: PetscCall(PetscDSSetRiemannSolver(newprob, fn, r));
3965: for (gn = 0; gn < numFields; ++gn) {
3966: const PetscInt g = fields ? fields[gn] : gn;
3967: PetscPointJac g0, g1, g2, g3;
3968: PetscPointJac g0p, g1p, g2p, g3p;
3969: PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd;
3971: if (g >= Nf) continue;
3972: PetscCall(PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3));
3973: PetscCall(PetscDSGetJacobianPreconditioner(prob, f, g, &g0p, &g1p, &g2p, &g3p));
3974: PetscCall(PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd));
3975: PetscCall(PetscDSSetJacobian(newprob, fn, gn, g0, g1, g2, g3));
3976: PetscCall(PetscDSSetJacobianPreconditioner(newprob, fn, gn, g0p, g1p, g2p, g3p));
3977: PetscCall(PetscDSSetBdJacobian(newprob, fn, gn, g0Bd, g1Bd, g2Bd, g3Bd));
3978: }
3979: }
3980: PetscFunctionReturn(PETSC_SUCCESS);
3981: }
3983: /*@
3984: PetscDSCopyEquations - Copy all pointwise function pointers to another `PetscDS`
3986: Not Collective
3988: Input Parameter:
3989: . prob - The `PetscDS` object
3991: Output Parameter:
3992: . newprob - The `PetscDS` copy
3994: Level: intermediate
3996: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3997: @*/
3998: PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob)
3999: {
4000: PetscWeakForm wf, newwf;
4001: PetscInt Nf, Ng;
4003: PetscFunctionBegin;
4006: PetscCall(PetscDSGetNumFields(prob, &Nf));
4007: PetscCall(PetscDSGetNumFields(newprob, &Ng));
4008: PetscCheck(Nf == Ng, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %" PetscInt_FMT " != %" PetscInt_FMT, Nf, Ng);
4009: PetscCall(PetscDSGetWeakForm(prob, &wf));
4010: PetscCall(PetscDSGetWeakForm(newprob, &newwf));
4011: PetscCall(PetscWeakFormCopy(wf, newwf));
4012: PetscFunctionReturn(PETSC_SUCCESS);
4013: }
4015: /*@
4016: PetscDSCopyConstants - Copy all constants to another `PetscDS`
4018: Not Collective
4020: Input Parameter:
4021: . prob - The `PetscDS` object
4023: Output Parameter:
4024: . newprob - The `PetscDS` copy
4026: Level: intermediate
4028: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
4029: @*/
4030: PetscErrorCode PetscDSCopyConstants(PetscDS prob, PetscDS newprob)
4031: {
4032: PetscInt Nc;
4033: const PetscScalar *constants;
4035: PetscFunctionBegin;
4038: PetscCall(PetscDSGetConstants(prob, &Nc, &constants));
4039: PetscCall(PetscDSSetConstants(newprob, Nc, (PetscScalar *)constants));
4040: PetscFunctionReturn(PETSC_SUCCESS);
4041: }
4043: /*@
4044: PetscDSCopyExactSolutions - Copy all exact solutions to another `PetscDS`
4046: Not Collective
4048: Input Parameter:
4049: . ds - The `PetscDS` object
4051: Output Parameter:
4052: . newds - The `PetscDS` copy
4054: Level: intermediate
4056: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
4057: @*/
4058: PetscErrorCode PetscDSCopyExactSolutions(PetscDS ds, PetscDS newds)
4059: {
4060: PetscSimplePointFn *sol;
4061: void *ctx;
4062: PetscInt Nf, f;
4064: PetscFunctionBegin;
4067: PetscCall(PetscDSGetNumFields(ds, &Nf));
4068: for (f = 0; f < Nf; ++f) {
4069: PetscCall(PetscDSGetExactSolution(ds, f, &sol, &ctx));
4070: PetscCall(PetscDSSetExactSolution(newds, f, sol, ctx));
4071: PetscCall(PetscDSGetExactSolutionTimeDerivative(ds, f, &sol, &ctx));
4072: PetscCall(PetscDSSetExactSolutionTimeDerivative(newds, f, sol, ctx));
4073: }
4074: PetscFunctionReturn(PETSC_SUCCESS);
4075: }
4077: PetscErrorCode PetscDSCopy(PetscDS ds, DM dmNew, PetscDS dsNew)
4078: {
4079: DSBoundary b;
4080: PetscInt cdim, Nf, f, d;
4081: PetscBool isCohesive;
4082: void *ctx;
4084: PetscFunctionBegin;
4085: PetscCall(PetscDSCopyConstants(ds, dsNew));
4086: PetscCall(PetscDSCopyExactSolutions(ds, dsNew));
4087: PetscCall(PetscDSSelectDiscretizations(ds, PETSC_DETERMINE, NULL, dsNew));
4088: PetscCall(PetscDSCopyEquations(ds, dsNew));
4089: PetscCall(PetscDSGetNumFields(ds, &Nf));
4090: for (f = 0; f < Nf; ++f) {
4091: PetscCall(PetscDSGetContext(ds, f, &ctx));
4092: PetscCall(PetscDSSetContext(dsNew, f, ctx));
4093: PetscCall(PetscDSGetCohesive(ds, f, &isCohesive));
4094: PetscCall(PetscDSSetCohesive(dsNew, f, isCohesive));
4095: PetscCall(PetscDSGetJetDegree(ds, f, &d));
4096: PetscCall(PetscDSSetJetDegree(dsNew, f, d));
4097: }
4098: if (Nf) {
4099: PetscCall(PetscDSGetCoordinateDimension(ds, &cdim));
4100: PetscCall(PetscDSSetCoordinateDimension(dsNew, cdim));
4101: }
4102: PetscCall(PetscDSCopyBoundary(ds, PETSC_DETERMINE, NULL, dsNew));
4103: for (b = dsNew->boundary; b; b = b->next) {
4104: PetscCall(DMGetLabel(dmNew, b->lname, &b->label));
4105: /* Do not check if label exists here, since p4est calls this for the reference tree which does not have the labels */
4106: //PetscCheck(b->label,PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Label %s missing in new DM", name);
4107: }
4108: PetscFunctionReturn(PETSC_SUCCESS);
4109: }
4111: PetscErrorCode PetscDSGetHeightSubspace(PetscDS prob, PetscInt height, PetscDS *subprob)
4112: {
4113: PetscInt dim, Nf, f;
4115: PetscFunctionBegin;
4117: PetscAssertPointer(subprob, 3);
4118: if (height == 0) {
4119: *subprob = prob;
4120: PetscFunctionReturn(PETSC_SUCCESS);
4121: }
4122: PetscCall(PetscDSGetNumFields(prob, &Nf));
4123: PetscCall(PetscDSGetSpatialDimension(prob, &dim));
4124: PetscCheck(height <= dim, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_OUTOFRANGE, "DS can only handle height in [0, %" PetscInt_FMT "], not %" PetscInt_FMT, dim, height);
4125: if (!prob->subprobs) PetscCall(PetscCalloc1(dim, &prob->subprobs));
4126: if (!prob->subprobs[height - 1]) {
4127: PetscInt cdim;
4129: PetscCall(PetscDSCreate(PetscObjectComm((PetscObject)prob), &prob->subprobs[height - 1]));
4130: PetscCall(PetscDSGetCoordinateDimension(prob, &cdim));
4131: PetscCall(PetscDSSetCoordinateDimension(prob->subprobs[height - 1], cdim));
4132: for (f = 0; f < Nf; ++f) {
4133: PetscFE subfe;
4134: PetscObject obj;
4135: PetscClassId id;
4137: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
4138: PetscCall(PetscObjectGetClassId(obj, &id));
4139: if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetHeightSubspace((PetscFE)obj, height, &subfe));
4140: else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unsupported discretization type for field %" PetscInt_FMT, f);
4141: PetscCall(PetscDSSetDiscretization(prob->subprobs[height - 1], f, (PetscObject)subfe));
4142: }
4143: }
4144: *subprob = prob->subprobs[height - 1];
4145: PetscFunctionReturn(PETSC_SUCCESS);
4146: }
4148: PetscErrorCode PetscDSPermuteQuadPoint(PetscDS ds, PetscInt ornt, PetscInt field, PetscInt q, PetscInt *qperm)
4149: {
4150: IS permIS;
4151: PetscQuadrature quad;
4152: DMPolytopeType ct;
4153: const PetscInt *perm;
4154: PetscInt Na, Nq;
4156: PetscFunctionBeginHot;
4157: PetscCall(PetscFEGetQuadrature((PetscFE)ds->disc[field], &quad));
4158: PetscCall(PetscQuadratureGetData(quad, NULL, NULL, &Nq, NULL, NULL));
4159: PetscCall(PetscQuadratureGetCellType(quad, &ct));
4160: PetscCheck(q >= 0 && q < Nq, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Quadrature point %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", q, Nq);
4161: Na = DMPolytopeTypeGetNumArrangements(ct) / 2;
4162: 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);
4163: if (!ds->quadPerm[(PetscInt)ct]) PetscCall(PetscQuadratureComputePermutations(quad, NULL, &ds->quadPerm[(PetscInt)ct]));
4164: permIS = ds->quadPerm[(PetscInt)ct][ornt + Na];
4165: PetscCall(ISGetIndices(permIS, &perm));
4166: *qperm = perm[q];
4167: PetscCall(ISRestoreIndices(permIS, &perm));
4168: PetscFunctionReturn(PETSC_SUCCESS);
4169: }
4171: PetscErrorCode PetscDSGetDiscType_Internal(PetscDS ds, PetscInt f, PetscDiscType *disctype)
4172: {
4173: PetscObject obj;
4174: PetscClassId id;
4175: PetscInt Nf;
4177: PetscFunctionBegin;
4179: PetscAssertPointer(disctype, 3);
4180: *disctype = PETSC_DISC_NONE;
4181: PetscCall(PetscDSGetNumFields(ds, &Nf));
4182: PetscCheck(f < Nf, PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_SIZ, "Field %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, Nf);
4183: PetscCall(PetscDSGetDiscretization(ds, f, &obj));
4184: if (obj) {
4185: PetscCall(PetscObjectGetClassId(obj, &id));
4186: if (id == PETSCFE_CLASSID) *disctype = PETSC_DISC_FE;
4187: else *disctype = PETSC_DISC_FV;
4188: }
4189: PetscFunctionReturn(PETSC_SUCCESS);
4190: }
4192: static PetscErrorCode PetscDSDestroy_Basic(PetscDS ds)
4193: {
4194: PetscFunctionBegin;
4195: PetscCall(PetscFree(ds->data));
4196: PetscFunctionReturn(PETSC_SUCCESS);
4197: }
4199: static PetscErrorCode PetscDSInitialize_Basic(PetscDS ds)
4200: {
4201: PetscFunctionBegin;
4202: ds->ops->setfromoptions = NULL;
4203: ds->ops->setup = NULL;
4204: ds->ops->view = NULL;
4205: ds->ops->destroy = PetscDSDestroy_Basic;
4206: PetscFunctionReturn(PETSC_SUCCESS);
4207: }
4209: /*MC
4210: PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions
4212: Level: intermediate
4214: .seealso: `PetscDSType`, `PetscDSCreate()`, `PetscDSSetType()`
4215: M*/
4217: PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS ds)
4218: {
4219: PetscDS_Basic *b;
4221: PetscFunctionBegin;
4223: PetscCall(PetscNew(&b));
4224: ds->data = b;
4226: PetscCall(PetscDSInitialize_Basic(ds));
4227: PetscFunctionReturn(PETSC_SUCCESS);
4228: }