Actual source code: dtds.c
1: #include <petsc/private/petscdsimpl.h>
3: PetscClassId PETSCDS_CLASSID = 0;
5: PetscFunctionList PetscDSList = NULL;
6: PetscBool PetscDSRegisterAllCalled = PETSC_FALSE;
8: /* A PetscDS (Discrete System) encodes a set of equations posed in a discrete space, which represents a set of
9: nonlinear continuum equations. The equations can have multiple fields, each field having a different
10: discretization. In addition, different pieces of the domain can have different field combinations and equations.
12: The DS provides the user a description of the approximation space on any given cell. It also gives pointwise
13: functions representing the equations.
15: Each field is associated with a label, marking the cells on which it is supported. Note that a field can be
16: supported on the closure of a cell not in the label due to overlap of the boundary of neighboring cells. The DM
17: then creates a DS for each set of cells with identical approximation spaces. When assembling, the user asks for
18: the space associated with a given cell. DMPlex uses the labels associated with each DS in the default integration loop.
19: */
21: /*@C
22: PetscDSRegister - Adds a new `PetscDS` implementation
24: Not Collective; No Fortran Support
26: Input Parameters:
27: + sname - The name of a new user-defined creation routine
28: - function - The creation routine itself
30: Example Usage:
31: .vb
32: PetscDSRegister("my_ds", MyPetscDSCreate);
33: .ve
35: Then, your PetscDS type can be chosen with the procedural interface via
36: .vb
37: PetscDSCreate(MPI_Comm, PetscDS *);
38: PetscDSSetType(PetscDS, "my_ds");
39: .ve
40: or at runtime via the option
41: .vb
42: -petscds_type my_ds
43: .ve
45: Level: advanced
47: Note:
48: `PetscDSRegister()` may be called multiple times to add several user-defined `PetscDSs`
50: .seealso: `PetscDSType`, `PetscDS`, `PetscDSRegisterAll()`, `PetscDSRegisterDestroy()`
51: @*/
52: PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS))
53: {
54: PetscFunctionBegin;
55: PetscCall(PetscFunctionListAdd(&PetscDSList, sname, function));
56: PetscFunctionReturn(PETSC_SUCCESS);
57: }
59: /*@C
60: PetscDSSetType - Builds a particular `PetscDS`
62: Collective; No Fortran Support
64: Input Parameters:
65: + prob - The `PetscDS` object
66: - name - The `PetscDSType`
68: Options Database Key:
69: . -petscds_type <type> - Sets the PetscDS type; use -help for a list of available types
71: Level: intermediate
73: .seealso: `PetscDSType`, `PetscDS`, `PetscDSGetType()`, `PetscDSCreate()`
74: @*/
75: PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name)
76: {
77: PetscErrorCode (*r)(PetscDS);
78: PetscBool match;
80: PetscFunctionBegin;
82: PetscCall(PetscObjectTypeCompare((PetscObject)prob, name, &match));
83: if (match) PetscFunctionReturn(PETSC_SUCCESS);
85: PetscCall(PetscDSRegisterAll());
86: PetscCall(PetscFunctionListFind(PetscDSList, name, &r));
87: PetscCheck(r, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscDS type: %s", name);
89: PetscTryTypeMethod(prob, destroy);
90: prob->ops->destroy = NULL;
92: PetscCall((*r)(prob));
93: PetscCall(PetscObjectChangeTypeName((PetscObject)prob, name));
94: PetscFunctionReturn(PETSC_SUCCESS);
95: }
97: /*@C
98: PetscDSGetType - Gets the `PetscDSType` name (as a string) from the `PetscDS`
100: Not Collective; No Fortran Support
102: Input Parameter:
103: . prob - The `PetscDS`
105: Output Parameter:
106: . name - The `PetscDSType` name
108: Level: intermediate
110: .seealso: `PetscDSType`, `PetscDS`, `PetscDSSetType()`, `PetscDSCreate()`
111: @*/
112: PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name)
113: {
114: PetscFunctionBegin;
116: PetscAssertPointer(name, 2);
117: PetscCall(PetscDSRegisterAll());
118: *name = ((PetscObject)prob)->type_name;
119: PetscFunctionReturn(PETSC_SUCCESS);
120: }
122: static PetscErrorCode PetscDSView_Ascii(PetscDS ds, PetscViewer viewer)
123: {
124: PetscViewerFormat format;
125: const PetscScalar *constants;
126: PetscInt Nf, numConstants, f;
128: PetscFunctionBegin;
129: PetscCall(PetscDSGetNumFields(ds, &Nf));
130: PetscCall(PetscViewerGetFormat(viewer, &format));
131: PetscCall(PetscViewerASCIIPrintf(viewer, "Discrete System with %" PetscInt_FMT " fields\n", Nf));
132: PetscCall(PetscViewerASCIIPushTab(viewer));
133: PetscCall(PetscViewerASCIIPrintf(viewer, " cell total dim %" PetscInt_FMT " total comp %" PetscInt_FMT "\n", ds->totDim, ds->totComp));
134: if (ds->isCohesive) PetscCall(PetscViewerASCIIPrintf(viewer, " cohesive cell\n"));
135: for (f = 0; f < Nf; ++f) {
136: DSBoundary b;
137: PetscObject obj;
138: PetscClassId id;
139: PetscQuadrature q;
140: const char *name;
141: PetscInt Nc, Nq, Nqc;
143: PetscCall(PetscDSGetDiscretization(ds, f, &obj));
144: PetscCall(PetscObjectGetClassId(obj, &id));
145: PetscCall(PetscObjectGetName(obj, &name));
146: PetscCall(PetscViewerASCIIPrintf(viewer, "Field %s", name ? name : "<unknown>"));
147: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
148: if (id == PETSCFE_CLASSID) {
149: PetscCall(PetscFEGetNumComponents((PetscFE)obj, &Nc));
150: PetscCall(PetscFEGetQuadrature((PetscFE)obj, &q));
151: PetscCall(PetscViewerASCIIPrintf(viewer, " FEM"));
152: } else if (id == PETSCFV_CLASSID) {
153: PetscCall(PetscFVGetNumComponents((PetscFV)obj, &Nc));
154: PetscCall(PetscFVGetQuadrature((PetscFV)obj, &q));
155: PetscCall(PetscViewerASCIIPrintf(viewer, " FVM"));
156: } else SETERRQ(PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f);
157: if (Nc > 1) PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " components", Nc));
158: else PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " component ", Nc));
159: if (ds->implicit[f]) PetscCall(PetscViewerASCIIPrintf(viewer, " (implicit)"));
160: else PetscCall(PetscViewerASCIIPrintf(viewer, " (explicit)"));
161: if (q) {
162: PetscCall(PetscQuadratureGetData(q, NULL, &Nqc, &Nq, NULL, NULL));
163: PetscCall(PetscViewerASCIIPrintf(viewer, " (Nq %" PetscInt_FMT " Nqc %" PetscInt_FMT ")", Nq, Nqc));
164: }
165: PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT "-jet", ds->jetDegree[f]));
166: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
167: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
168: PetscCall(PetscViewerASCIIPushTab(viewer));
169: if (id == PETSCFE_CLASSID) PetscCall(PetscFEView((PetscFE)obj, viewer));
170: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVView((PetscFV)obj, viewer));
171: PetscCall(PetscViewerASCIIPopTab(viewer));
173: for (b = ds->boundary; b; b = b->next) {
174: char *name;
175: PetscInt c, i;
177: if (b->field != f) continue;
178: PetscCall(PetscViewerASCIIPushTab(viewer));
179: PetscCall(PetscViewerASCIIPrintf(viewer, "Boundary %s (%s) %s\n", b->name, b->lname, DMBoundaryConditionTypes[b->type]));
180: if (!b->Nc) {
181: PetscCall(PetscViewerASCIIPrintf(viewer, " all components\n"));
182: } else {
183: PetscCall(PetscViewerASCIIPrintf(viewer, " components: "));
184: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
185: for (c = 0; c < b->Nc; ++c) {
186: if (c > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", "));
187: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->comps[c]));
188: }
189: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
190: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
191: }
192: PetscCall(PetscViewerASCIIPrintf(viewer, " values: "));
193: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
194: for (i = 0; i < b->Nv; ++i) {
195: if (i > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", "));
196: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->values[i]));
197: }
198: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
199: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
200: #if defined(__clang__)
201: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat-pedantic")
202: #elif defined(__GNUC__) || defined(__GNUG__)
203: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat")
204: #endif
205: if (b->func) {
206: PetscCall(PetscDLAddr(b->func, &name));
207: if (name) PetscCall(PetscViewerASCIIPrintf(viewer, " func: %s\n", name));
208: else PetscCall(PetscViewerASCIIPrintf(viewer, " func: %p\n", b->func));
209: PetscCall(PetscFree(name));
210: }
211: if (b->func_t) {
212: PetscCall(PetscDLAddr(b->func_t, &name));
213: if (name) PetscCall(PetscViewerASCIIPrintf(viewer, " func_t: %s\n", name));
214: else PetscCall(PetscViewerASCIIPrintf(viewer, " func_t: %p\n", b->func_t));
215: PetscCall(PetscFree(name));
216: }
217: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_END()
218: PetscCall(PetscWeakFormView(b->wf, viewer));
219: PetscCall(PetscViewerASCIIPopTab(viewer));
220: }
221: }
222: PetscCall(PetscDSGetConstants(ds, &numConstants, &constants));
223: if (numConstants) {
224: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " constants\n", numConstants));
225: PetscCall(PetscViewerASCIIPushTab(viewer));
226: for (f = 0; f < numConstants; ++f) PetscCall(PetscViewerASCIIPrintf(viewer, "%g\n", (double)PetscRealPart(constants[f])));
227: PetscCall(PetscViewerASCIIPopTab(viewer));
228: }
229: PetscCall(PetscWeakFormView(ds->wf, viewer));
230: PetscCall(PetscViewerASCIIPopTab(viewer));
231: PetscFunctionReturn(PETSC_SUCCESS);
232: }
234: /*@C
235: PetscDSViewFromOptions - View a `PetscDS` based on values in the options database
237: Collective
239: Input Parameters:
240: + A - the `PetscDS` object
241: . obj - Optional object that provides the options prefix used in the search
242: - name - command line option
244: Level: intermediate
246: .seealso: `PetscDSType`, `PetscDS`, `PetscDSView()`, `PetscObjectViewFromOptions()`, `PetscDSCreate()`
247: @*/
248: PetscErrorCode PetscDSViewFromOptions(PetscDS A, PetscObject obj, const char name[])
249: {
250: PetscFunctionBegin;
252: PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
253: PetscFunctionReturn(PETSC_SUCCESS);
254: }
256: /*@C
257: PetscDSView - Views a `PetscDS`
259: Collective
261: Input Parameters:
262: + prob - the `PetscDS` object to view
263: - v - the viewer
265: Level: developer
267: .seealso: `PetscDSType`, `PetscDS`, `PetscViewer`, `PetscDSDestroy()`, `PetscDSViewFromOptions()`
268: @*/
269: PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v)
270: {
271: PetscBool iascii;
273: PetscFunctionBegin;
275: if (!v) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)prob), &v));
277: PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &iascii));
278: if (iascii) PetscCall(PetscDSView_Ascii(prob, v));
279: PetscTryTypeMethod(prob, view, v);
280: PetscFunctionReturn(PETSC_SUCCESS);
281: }
283: /*@
284: PetscDSSetFromOptions - sets parameters in a `PetscDS` from the options database
286: Collective
288: Input Parameter:
289: . prob - the `PetscDS` object to set options for
291: Options Database Keys:
292: + -petscds_type <type> - Set the `PetscDS` type
293: . -petscds_view <view opt> - View the `PetscDS`
294: . -petscds_jac_pre - Turn formation of a separate Jacobian preconditioner on or off
295: . -bc_<name> <ids> - Specify a list of label ids for a boundary condition
296: - -bc_<name>_comp <comps> - Specify a list of field components to constrain for a boundary condition
298: Level: intermediate
300: .seealso: `PetscDS`, `PetscDSView()`
301: @*/
302: PetscErrorCode PetscDSSetFromOptions(PetscDS prob)
303: {
304: DSBoundary b;
305: const char *defaultType;
306: char name[256];
307: PetscBool flg;
309: PetscFunctionBegin;
311: if (!((PetscObject)prob)->type_name) {
312: defaultType = PETSCDSBASIC;
313: } else {
314: defaultType = ((PetscObject)prob)->type_name;
315: }
316: PetscCall(PetscDSRegisterAll());
318: PetscObjectOptionsBegin((PetscObject)prob);
319: for (b = prob->boundary; b; b = b->next) {
320: char optname[1024];
321: PetscInt ids[1024], len = 1024;
322: PetscBool flg;
324: PetscCall(PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name));
325: PetscCall(PetscMemzero(ids, sizeof(ids)));
326: PetscCall(PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg));
327: if (flg) {
328: b->Nv = len;
329: PetscCall(PetscFree(b->values));
330: PetscCall(PetscMalloc1(len, &b->values));
331: PetscCall(PetscArraycpy(b->values, ids, len));
332: PetscCall(PetscWeakFormRewriteKeys(b->wf, b->label, len, b->values));
333: }
334: len = 1024;
335: PetscCall(PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name));
336: PetscCall(PetscMemzero(ids, sizeof(ids)));
337: PetscCall(PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg));
338: if (flg) {
339: b->Nc = len;
340: PetscCall(PetscFree(b->comps));
341: PetscCall(PetscMalloc1(len, &b->comps));
342: PetscCall(PetscArraycpy(b->comps, ids, len));
343: }
344: }
345: PetscCall(PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg));
346: if (flg) {
347: PetscCall(PetscDSSetType(prob, name));
348: } else if (!((PetscObject)prob)->type_name) {
349: PetscCall(PetscDSSetType(prob, defaultType));
350: }
351: PetscCall(PetscOptionsBool("-petscds_jac_pre", "Discrete System", "PetscDSUseJacobianPreconditioner", prob->useJacPre, &prob->useJacPre, &flg));
352: PetscCall(PetscOptionsBool("-petscds_force_quad", "Discrete System", "PetscDSSetForceQuad", prob->forceQuad, &prob->forceQuad, &flg));
353: PetscTryTypeMethod(prob, setfromoptions);
354: /* process any options handlers added with PetscObjectAddOptionsHandler() */
355: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)prob, PetscOptionsObject));
356: PetscOptionsEnd();
357: if (prob->Nf) PetscCall(PetscDSViewFromOptions(prob, NULL, "-petscds_view"));
358: PetscFunctionReturn(PETSC_SUCCESS);
359: }
361: /*@C
362: PetscDSSetUp - Construct data structures for the `PetscDS`
364: Collective
366: Input Parameter:
367: . prob - the `PetscDS` object to setup
369: Level: developer
371: .seealso: `PetscDS`, `PetscDSView()`, `PetscDSDestroy()`
372: @*/
373: PetscErrorCode PetscDSSetUp(PetscDS prob)
374: {
375: const PetscInt Nf = prob->Nf;
376: PetscBool hasH = PETSC_FALSE;
377: PetscInt maxOrder[4] = {-1, -1, -1, -1};
378: PetscInt dim, dimEmbed, NbMax = 0, NcMax = 0, NqMax = 0, NsMax = 1, f;
380: PetscFunctionBegin;
382: if (prob->setup) PetscFunctionReturn(PETSC_SUCCESS);
383: /* Calculate sizes */
384: PetscCall(PetscDSGetSpatialDimension(prob, &dim));
385: PetscCall(PetscDSGetCoordinateDimension(prob, &dimEmbed));
386: prob->totDim = prob->totComp = 0;
387: PetscCall(PetscMalloc2(Nf, &prob->Nc, Nf, &prob->Nb));
388: PetscCall(PetscCalloc2(Nf + 1, &prob->off, Nf + 1, &prob->offDer));
389: PetscCall(PetscCalloc6(Nf + 1, &prob->offCohesive[0], Nf + 1, &prob->offCohesive[1], Nf + 1, &prob->offCohesive[2], Nf + 1, &prob->offDerCohesive[0], Nf + 1, &prob->offDerCohesive[1], Nf + 1, &prob->offDerCohesive[2]));
390: PetscCall(PetscMalloc2(Nf, &prob->T, Nf, &prob->Tf));
391: if (prob->forceQuad) {
392: // Note: This assumes we have one kind of cell at each dimension.
393: // We can fix this by having quadrature hold the celltype
394: PetscQuadrature maxQuad[4] = {NULL, NULL, NULL, NULL};
396: for (f = 0; f < Nf; ++f) {
397: PetscObject obj;
398: PetscClassId id;
399: PetscQuadrature q = NULL, fq = NULL;
400: PetscInt dim = -1, order = -1, forder = -1;
402: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
403: if (!obj) continue;
404: PetscCall(PetscObjectGetClassId(obj, &id));
405: if (id == PETSCFE_CLASSID) {
406: PetscFE fe = (PetscFE)obj;
408: PetscCall(PetscFEGetQuadrature(fe, &q));
409: PetscCall(PetscFEGetFaceQuadrature(fe, &fq));
410: } else if (id == PETSCFV_CLASSID) {
411: PetscFV fv = (PetscFV)obj;
413: PetscCall(PetscFVGetQuadrature(fv, &q));
414: }
415: if (q) {
416: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
417: PetscCall(PetscQuadratureGetOrder(q, &order));
418: if (order > maxOrder[dim]) {
419: maxOrder[dim] = order;
420: maxQuad[dim] = q;
421: }
422: }
423: if (fq) {
424: PetscCall(PetscQuadratureGetData(fq, &dim, NULL, NULL, NULL, NULL));
425: PetscCall(PetscQuadratureGetOrder(fq, &forder));
426: if (forder > maxOrder[dim]) {
427: maxOrder[dim] = forder;
428: maxQuad[dim] = fq;
429: }
430: }
431: }
432: for (f = 0; f < Nf; ++f) {
433: PetscObject obj;
434: PetscClassId id;
435: PetscQuadrature q;
436: PetscInt dim;
438: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
439: if (!obj) continue;
440: PetscCall(PetscObjectGetClassId(obj, &id));
441: if (id == PETSCFE_CLASSID) {
442: PetscFE fe = (PetscFE)obj;
444: PetscCall(PetscFEGetQuadrature(fe, &q));
445: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
446: PetscCall(PetscFESetQuadrature(fe, maxQuad[dim]));
447: PetscCall(PetscFESetFaceQuadrature(fe, dim ? maxQuad[dim - 1] : NULL));
448: } else if (id == PETSCFV_CLASSID) {
449: PetscFV fv = (PetscFV)obj;
451: PetscCall(PetscFVGetQuadrature(fv, &q));
452: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
453: PetscCall(PetscFVSetQuadrature(fv, maxQuad[dim]));
454: }
455: }
456: }
457: for (f = 0; f < Nf; ++f) {
458: PetscObject obj;
459: PetscClassId id;
460: PetscQuadrature q = NULL;
461: PetscInt Nq = 0, Nb, Nc;
463: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
464: if (prob->jetDegree[f] > 1) hasH = PETSC_TRUE;
465: if (!obj) {
466: /* Empty mesh */
467: Nb = Nc = 0;
468: prob->T[f] = prob->Tf[f] = NULL;
469: } else {
470: PetscCall(PetscObjectGetClassId(obj, &id));
471: if (id == PETSCFE_CLASSID) {
472: PetscFE fe = (PetscFE)obj;
474: PetscCall(PetscFEGetQuadrature(fe, &q));
475: {
476: PetscQuadrature fq;
477: PetscInt dim, order;
479: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
480: PetscCall(PetscQuadratureGetOrder(q, &order));
481: if (maxOrder[dim] < 0) maxOrder[dim] = order;
482: PetscCheck(order == maxOrder[dim], PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Field %" PetscInt_FMT " cell quadrature order %" PetscInt_FMT " != %" PetscInt_FMT " DS cell quadrature order", f, order, maxOrder[dim]);
483: PetscCall(PetscFEGetFaceQuadrature(fe, &fq));
484: if (fq) {
485: PetscCall(PetscQuadratureGetData(fq, &dim, NULL, NULL, NULL, NULL));
486: PetscCall(PetscQuadratureGetOrder(fq, &order));
487: if (maxOrder[dim] < 0) maxOrder[dim] = order;
488: PetscCheck(order == maxOrder[dim], PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Field %" PetscInt_FMT " face quadrature order %" PetscInt_FMT " != %" PetscInt_FMT " DS face quadrature order", f, order, maxOrder[dim]);
489: }
490: }
491: PetscCall(PetscFEGetDimension(fe, &Nb));
492: PetscCall(PetscFEGetNumComponents(fe, &Nc));
493: PetscCall(PetscFEGetCellTabulation(fe, prob->jetDegree[f], &prob->T[f]));
494: PetscCall(PetscFEGetFaceTabulation(fe, prob->jetDegree[f], &prob->Tf[f]));
495: } else if (id == PETSCFV_CLASSID) {
496: PetscFV fv = (PetscFV)obj;
498: PetscCall(PetscFVGetQuadrature(fv, &q));
499: PetscCall(PetscFVGetNumComponents(fv, &Nc));
500: Nb = Nc;
501: PetscCall(PetscFVGetCellTabulation(fv, &prob->T[f]));
502: /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */
503: } else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f);
504: }
505: prob->Nc[f] = Nc;
506: prob->Nb[f] = Nb;
507: prob->off[f + 1] = Nc + prob->off[f];
508: prob->offDer[f + 1] = Nc * dim + prob->offDer[f];
509: prob->offCohesive[0][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[0][f];
510: prob->offDerCohesive[0][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[0][f];
511: prob->offCohesive[1][f] = (prob->cohesive[f] ? 0 : Nc) + prob->offCohesive[0][f];
512: prob->offDerCohesive[1][f] = (prob->cohesive[f] ? 0 : Nc) * dimEmbed + prob->offDerCohesive[0][f];
513: prob->offCohesive[2][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[2][f];
514: prob->offDerCohesive[2][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[2][f];
515: if (q) PetscCall(PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL));
516: NqMax = PetscMax(NqMax, Nq);
517: NbMax = PetscMax(NbMax, Nb);
518: NcMax = PetscMax(NcMax, Nc);
519: prob->totDim += Nb;
520: prob->totComp += Nc;
521: /* There are two faces for all fields on a cohesive cell, except for cohesive fields */
522: if (prob->isCohesive && !prob->cohesive[f]) prob->totDim += Nb;
523: }
524: prob->offCohesive[1][Nf] = prob->offCohesive[0][Nf];
525: prob->offDerCohesive[1][Nf] = prob->offDerCohesive[0][Nf];
526: /* Allocate works space */
527: NsMax = 2; /* A non-cohesive discretizations can be used on a cohesive cell, so we need this extra workspace for all DS */
528: PetscCall(PetscMalloc3(NsMax * prob->totComp, &prob->u, NsMax * prob->totComp, &prob->u_t, NsMax * prob->totComp * dimEmbed + (hasH ? NsMax * prob->totComp * dimEmbed * dimEmbed : 0), &prob->u_x));
529: PetscCall(PetscMalloc5(dimEmbed, &prob->x, NbMax * NcMax, &prob->basisReal, NbMax * NcMax * dimEmbed, &prob->basisDerReal, NbMax * NcMax, &prob->testReal, NbMax * NcMax * dimEmbed, &prob->testDerReal));
530: PetscCall(PetscMalloc6(NsMax * NqMax * NcMax, &prob->f0, NsMax * NqMax * NcMax * dimEmbed, &prob->f1, NsMax * NsMax * NqMax * NcMax * NcMax, &prob->g0, NsMax * NsMax * NqMax * NcMax * NcMax * dimEmbed, &prob->g1, NsMax * NsMax * NqMax * NcMax * NcMax * dimEmbed,
531: &prob->g2, NsMax * NsMax * NqMax * NcMax * NcMax * dimEmbed * dimEmbed, &prob->g3));
532: PetscTryTypeMethod(prob, setup);
533: prob->setup = PETSC_TRUE;
534: PetscFunctionReturn(PETSC_SUCCESS);
535: }
537: static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob)
538: {
539: PetscFunctionBegin;
540: PetscCall(PetscFree2(prob->Nc, prob->Nb));
541: PetscCall(PetscFree2(prob->off, prob->offDer));
542: PetscCall(PetscFree6(prob->offCohesive[0], prob->offCohesive[1], prob->offCohesive[2], prob->offDerCohesive[0], prob->offDerCohesive[1], prob->offDerCohesive[2]));
543: PetscCall(PetscFree2(prob->T, prob->Tf));
544: PetscCall(PetscFree3(prob->u, prob->u_t, prob->u_x));
545: PetscCall(PetscFree5(prob->x, prob->basisReal, prob->basisDerReal, prob->testReal, prob->testDerReal));
546: PetscCall(PetscFree6(prob->f0, prob->f1, prob->g0, prob->g1, prob->g2, prob->g3));
547: PetscFunctionReturn(PETSC_SUCCESS);
548: }
550: static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew)
551: {
552: PetscObject *tmpd;
553: PetscBool *tmpi;
554: PetscInt *tmpk;
555: PetscBool *tmpc;
556: PetscPointFunc *tmpup;
557: PetscSimplePointFunc *tmpexactSol, *tmpexactSol_t;
558: void **tmpexactCtx, **tmpexactCtx_t;
559: void **tmpctx;
560: PetscInt Nf = prob->Nf, f;
562: PetscFunctionBegin;
563: if (Nf >= NfNew) PetscFunctionReturn(PETSC_SUCCESS);
564: prob->setup = PETSC_FALSE;
565: PetscCall(PetscDSDestroyStructs_Static(prob));
566: PetscCall(PetscMalloc4(NfNew, &tmpd, NfNew, &tmpi, NfNew, &tmpc, NfNew, &tmpk));
567: for (f = 0; f < Nf; ++f) {
568: tmpd[f] = prob->disc[f];
569: tmpi[f] = prob->implicit[f];
570: tmpc[f] = prob->cohesive[f];
571: tmpk[f] = prob->jetDegree[f];
572: }
573: for (f = Nf; f < NfNew; ++f) {
574: tmpd[f] = NULL;
575: tmpi[f] = PETSC_TRUE, tmpc[f] = PETSC_FALSE;
576: tmpk[f] = 1;
577: }
578: PetscCall(PetscFree4(prob->disc, prob->implicit, prob->cohesive, prob->jetDegree));
579: PetscCall(PetscWeakFormSetNumFields(prob->wf, NfNew));
580: prob->Nf = NfNew;
581: prob->disc = tmpd;
582: prob->implicit = tmpi;
583: prob->cohesive = tmpc;
584: prob->jetDegree = tmpk;
585: PetscCall(PetscCalloc2(NfNew, &tmpup, NfNew, &tmpctx));
586: for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f];
587: for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f];
588: for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL;
589: for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL;
590: PetscCall(PetscFree2(prob->update, prob->ctx));
591: prob->update = tmpup;
592: prob->ctx = tmpctx;
593: PetscCall(PetscCalloc4(NfNew, &tmpexactSol, NfNew, &tmpexactCtx, NfNew, &tmpexactSol_t, NfNew, &tmpexactCtx_t));
594: for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f];
595: for (f = 0; f < Nf; ++f) tmpexactCtx[f] = prob->exactCtx[f];
596: for (f = 0; f < Nf; ++f) tmpexactSol_t[f] = prob->exactSol_t[f];
597: for (f = 0; f < Nf; ++f) tmpexactCtx_t[f] = prob->exactCtx_t[f];
598: for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL;
599: for (f = Nf; f < NfNew; ++f) tmpexactCtx[f] = NULL;
600: for (f = Nf; f < NfNew; ++f) tmpexactSol_t[f] = NULL;
601: for (f = Nf; f < NfNew; ++f) tmpexactCtx_t[f] = NULL;
602: PetscCall(PetscFree4(prob->exactSol, prob->exactCtx, prob->exactSol_t, prob->exactCtx_t));
603: prob->exactSol = tmpexactSol;
604: prob->exactCtx = tmpexactCtx;
605: prob->exactSol_t = tmpexactSol_t;
606: prob->exactCtx_t = tmpexactCtx_t;
607: PetscFunctionReturn(PETSC_SUCCESS);
608: }
610: /*@
611: PetscDSDestroy - Destroys a `PetscDS` object
613: Collective
615: Input Parameter:
616: . ds - the `PetscDS` object to destroy
618: Level: developer
620: .seealso: `PetscDSView()`
621: @*/
622: PetscErrorCode PetscDSDestroy(PetscDS *ds)
623: {
624: PetscInt f;
626: PetscFunctionBegin;
627: if (!*ds) PetscFunctionReturn(PETSC_SUCCESS);
630: if (--((PetscObject)(*ds))->refct > 0) {
631: *ds = NULL;
632: PetscFunctionReturn(PETSC_SUCCESS);
633: }
634: ((PetscObject)(*ds))->refct = 0;
635: if ((*ds)->subprobs) {
636: PetscInt dim, d;
638: PetscCall(PetscDSGetSpatialDimension(*ds, &dim));
639: for (d = 0; d < dim; ++d) PetscCall(PetscDSDestroy(&(*ds)->subprobs[d]));
640: }
641: PetscCall(PetscFree((*ds)->subprobs));
642: PetscCall(PetscDSDestroyStructs_Static(*ds));
643: for (f = 0; f < (*ds)->Nf; ++f) PetscCall(PetscObjectDereference((*ds)->disc[f]));
644: PetscCall(PetscFree4((*ds)->disc, (*ds)->implicit, (*ds)->cohesive, (*ds)->jetDegree));
645: PetscCall(PetscWeakFormDestroy(&(*ds)->wf));
646: PetscCall(PetscFree2((*ds)->update, (*ds)->ctx));
647: PetscCall(PetscFree4((*ds)->exactSol, (*ds)->exactCtx, (*ds)->exactSol_t, (*ds)->exactCtx_t));
648: PetscTryTypeMethod((*ds), destroy);
649: PetscCall(PetscDSDestroyBoundary(*ds));
650: PetscCall(PetscFree((*ds)->constants));
651: for (PetscInt c = 0; c < DM_NUM_POLYTOPES; ++c) {
652: const PetscInt Na = DMPolytopeTypeGetNumArrangments((DMPolytopeType)c);
653: if ((*ds)->quadPerm[c])
654: for (PetscInt o = 0; o < Na; ++o) PetscCall(ISDestroy(&(*ds)->quadPerm[c][o]));
655: PetscCall(PetscFree((*ds)->quadPerm[c]));
656: (*ds)->quadPerm[c] = NULL;
657: }
658: PetscCall(PetscHeaderDestroy(ds));
659: PetscFunctionReturn(PETSC_SUCCESS);
660: }
662: /*@
663: PetscDSCreate - Creates an empty `PetscDS` object. The type can then be set with `PetscDSSetType()`.
665: Collective
667: Input Parameter:
668: . comm - The communicator for the `PetscDS` object
670: Output Parameter:
671: . ds - The `PetscDS` object
673: Level: beginner
675: .seealso: `PetscDS`, `PetscDSSetType()`, `PETSCDSBASIC`, `PetscDSType`
676: @*/
677: PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *ds)
678: {
679: PetscDS p;
681: PetscFunctionBegin;
682: PetscAssertPointer(ds, 2);
683: *ds = NULL;
684: PetscCall(PetscDSInitializePackage());
686: PetscCall(PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView));
688: p->Nf = 0;
689: p->setup = PETSC_FALSE;
690: p->numConstants = 0;
691: p->constants = NULL;
692: p->dimEmbed = -1;
693: p->useJacPre = PETSC_TRUE;
694: p->forceQuad = PETSC_TRUE;
695: PetscCall(PetscWeakFormCreate(comm, &p->wf));
696: PetscCall(PetscArrayzero(p->quadPerm, DM_NUM_POLYTOPES));
698: *ds = p;
699: PetscFunctionReturn(PETSC_SUCCESS);
700: }
702: /*@
703: PetscDSGetNumFields - Returns the number of fields in the `PetscDS`
705: Not Collective
707: Input Parameter:
708: . prob - The `PetscDS` object
710: Output Parameter:
711: . Nf - The number of fields
713: Level: beginner
715: .seealso: `PetscDS`, `PetscDSGetSpatialDimension()`, `PetscDSCreate()`
716: @*/
717: PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf)
718: {
719: PetscFunctionBegin;
721: PetscAssertPointer(Nf, 2);
722: *Nf = prob->Nf;
723: PetscFunctionReturn(PETSC_SUCCESS);
724: }
726: /*@
727: PetscDSGetSpatialDimension - Returns the spatial dimension of the `PetscDS`, meaning the topological dimension of the discretizations
729: Not Collective
731: Input Parameter:
732: . prob - The `PetscDS` object
734: Output Parameter:
735: . dim - The spatial dimension
737: Level: beginner
739: .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
740: @*/
741: PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim)
742: {
743: PetscFunctionBegin;
745: PetscAssertPointer(dim, 2);
746: *dim = 0;
747: if (prob->Nf) {
748: PetscObject obj;
749: PetscClassId id;
751: PetscCall(PetscDSGetDiscretization(prob, 0, &obj));
752: if (obj) {
753: PetscCall(PetscObjectGetClassId(obj, &id));
754: if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetSpatialDimension((PetscFE)obj, dim));
755: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetSpatialDimension((PetscFV)obj, dim));
756: else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
757: }
758: }
759: PetscFunctionReturn(PETSC_SUCCESS);
760: }
762: /*@
763: PetscDSGetCoordinateDimension - Returns the coordinate dimension of the `PetscDS`, meaning the dimension of the space into which the discretiaztions are embedded
765: Not Collective
767: Input Parameter:
768: . prob - The `PetscDS` object
770: Output Parameter:
771: . dimEmbed - The coordinate dimension
773: Level: beginner
775: .seealso: `PetscDS`, `PetscDSSetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
776: @*/
777: PetscErrorCode PetscDSGetCoordinateDimension(PetscDS prob, PetscInt *dimEmbed)
778: {
779: PetscFunctionBegin;
781: PetscAssertPointer(dimEmbed, 2);
782: PetscCheck(prob->dimEmbed >= 0, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONGSTATE, "No coordinate dimension set for this DS");
783: *dimEmbed = prob->dimEmbed;
784: PetscFunctionReturn(PETSC_SUCCESS);
785: }
787: /*@
788: PetscDSSetCoordinateDimension - Set the coordinate dimension of the `PetscDS`, meaning the dimension of the space into which the discretiaztions are embedded
790: Logically Collective
792: Input Parameters:
793: + prob - The `PetscDS` object
794: - dimEmbed - The coordinate dimension
796: Level: beginner
798: .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
799: @*/
800: PetscErrorCode PetscDSSetCoordinateDimension(PetscDS prob, PetscInt dimEmbed)
801: {
802: PetscFunctionBegin;
804: PetscCheck(dimEmbed >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Coordinate dimension must be non-negative, not %" PetscInt_FMT, dimEmbed);
805: prob->dimEmbed = dimEmbed;
806: PetscFunctionReturn(PETSC_SUCCESS);
807: }
809: /*@
810: PetscDSGetForceQuad - Returns the flag to force matching quadratures among the field discretizations
812: Not collective
814: Input Parameter:
815: . ds - The `PetscDS` object
817: Output Parameter:
818: . forceQuad - The flag
820: Level: intermediate
822: .seealso: `PetscDS`, `PetscDSSetForceQuad()`, `PetscDSGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
823: @*/
824: PetscErrorCode PetscDSGetForceQuad(PetscDS ds, PetscBool *forceQuad)
825: {
826: PetscFunctionBegin;
828: PetscAssertPointer(forceQuad, 2);
829: *forceQuad = ds->forceQuad;
830: PetscFunctionReturn(PETSC_SUCCESS);
831: }
833: /*@
834: PetscDSSetForceQuad - Set the flag to force matching quadratures among the field discretizations
836: Logically collective on ds
838: Input Parameters:
839: + ds - The `PetscDS` object
840: - forceQuad - The flag
842: Level: intermediate
844: .seealso: `PetscDS`, `PetscDSGetForceQuad()`, `PetscDSGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
845: @*/
846: PetscErrorCode PetscDSSetForceQuad(PetscDS ds, PetscBool forceQuad)
847: {
848: PetscFunctionBegin;
850: ds->forceQuad = forceQuad;
851: PetscFunctionReturn(PETSC_SUCCESS);
852: }
854: /*@
855: PetscDSIsCohesive - Returns the flag indicating that this `PetscDS` is for a cohesive cell
857: Not Collective
859: Input Parameter:
860: . ds - The `PetscDS` object
862: Output Parameter:
863: . isCohesive - The flag
865: Level: developer
867: .seealso: `PetscDS`, `PetscDSGetNumCohesive()`, `PetscDSGetCohesive()`, `PetscDSSetCohesive()`, `PetscDSCreate()`
868: @*/
869: PetscErrorCode PetscDSIsCohesive(PetscDS ds, PetscBool *isCohesive)
870: {
871: PetscFunctionBegin;
873: PetscAssertPointer(isCohesive, 2);
874: *isCohesive = ds->isCohesive;
875: PetscFunctionReturn(PETSC_SUCCESS);
876: }
878: /*@
879: PetscDSGetNumCohesive - Returns the number of cohesive fields, meaning those defined on the interior of a cohesive cell
881: Not Collective
883: Input Parameter:
884: . ds - The `PetscDS` object
886: Output Parameter:
887: . numCohesive - The number of cohesive fields
889: Level: developer
891: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSCreate()`
892: @*/
893: PetscErrorCode PetscDSGetNumCohesive(PetscDS ds, PetscInt *numCohesive)
894: {
895: PetscInt f;
897: PetscFunctionBegin;
899: PetscAssertPointer(numCohesive, 2);
900: *numCohesive = 0;
901: for (f = 0; f < ds->Nf; ++f) *numCohesive += ds->cohesive[f] ? 1 : 0;
902: PetscFunctionReturn(PETSC_SUCCESS);
903: }
905: /*@
906: PetscDSGetCohesive - Returns the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell
908: Not Collective
910: Input Parameters:
911: + ds - The `PetscDS` object
912: - f - The field index
914: Output Parameter:
915: . isCohesive - The flag
917: Level: developer
919: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
920: @*/
921: PetscErrorCode PetscDSGetCohesive(PetscDS ds, PetscInt f, PetscBool *isCohesive)
922: {
923: PetscFunctionBegin;
925: PetscAssertPointer(isCohesive, 3);
926: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
927: *isCohesive = ds->cohesive[f];
928: PetscFunctionReturn(PETSC_SUCCESS);
929: }
931: /*@
932: PetscDSSetCohesive - Set the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell
934: Not Collective
936: Input Parameters:
937: + ds - The `PetscDS` object
938: . f - The field index
939: - isCohesive - The flag for a cohesive field
941: Level: developer
943: .seealso: `PetscDS`, `PetscDSGetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
944: @*/
945: PetscErrorCode PetscDSSetCohesive(PetscDS ds, PetscInt f, PetscBool isCohesive)
946: {
947: PetscInt i;
949: PetscFunctionBegin;
951: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
952: ds->cohesive[f] = isCohesive;
953: ds->isCohesive = PETSC_FALSE;
954: for (i = 0; i < ds->Nf; ++i) ds->isCohesive = ds->isCohesive || ds->cohesive[f] ? PETSC_TRUE : PETSC_FALSE;
955: PetscFunctionReturn(PETSC_SUCCESS);
956: }
958: /*@
959: PetscDSGetTotalDimension - Returns the total size of the approximation space for this system
961: Not Collective
963: Input Parameter:
964: . prob - The `PetscDS` object
966: Output Parameter:
967: . dim - The total problem dimension
969: Level: beginner
971: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
972: @*/
973: PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim)
974: {
975: PetscFunctionBegin;
977: PetscCall(PetscDSSetUp(prob));
978: PetscAssertPointer(dim, 2);
979: *dim = prob->totDim;
980: PetscFunctionReturn(PETSC_SUCCESS);
981: }
983: /*@
984: PetscDSGetTotalComponents - Returns the total number of components in this system
986: Not Collective
988: Input Parameter:
989: . prob - The `PetscDS` object
991: Output Parameter:
992: . Nc - The total number of components
994: Level: beginner
996: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
997: @*/
998: PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc)
999: {
1000: PetscFunctionBegin;
1002: PetscCall(PetscDSSetUp(prob));
1003: PetscAssertPointer(Nc, 2);
1004: *Nc = prob->totComp;
1005: PetscFunctionReturn(PETSC_SUCCESS);
1006: }
1008: /*@
1009: PetscDSGetDiscretization - Returns the discretization object for the given field
1011: Not Collective
1013: Input Parameters:
1014: + prob - The `PetscDS` object
1015: - f - The field number
1017: Output Parameter:
1018: . disc - The discretization object
1020: Level: beginner
1022: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1023: @*/
1024: PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc)
1025: {
1026: PetscFunctionBeginHot;
1028: PetscAssertPointer(disc, 3);
1029: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
1030: *disc = prob->disc[f];
1031: PetscFunctionReturn(PETSC_SUCCESS);
1032: }
1034: /*@
1035: PetscDSSetDiscretization - Sets the discretization object for the given field
1037: Not Collective
1039: Input Parameters:
1040: + prob - The `PetscDS` object
1041: . f - The field number
1042: - disc - The discretization object
1044: Level: beginner
1046: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSGetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1047: @*/
1048: PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc)
1049: {
1050: PetscFunctionBegin;
1052: if (disc) PetscAssertPointer(disc, 3);
1053: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1054: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
1055: PetscCall(PetscObjectDereference(prob->disc[f]));
1056: prob->disc[f] = disc;
1057: PetscCall(PetscObjectReference(disc));
1058: if (disc) {
1059: PetscClassId id;
1061: PetscCall(PetscObjectGetClassId(disc, &id));
1062: if (id == PETSCFE_CLASSID) {
1063: PetscCall(PetscDSSetImplicit(prob, f, PETSC_TRUE));
1064: } else if (id == PETSCFV_CLASSID) {
1065: PetscCall(PetscDSSetImplicit(prob, f, PETSC_FALSE));
1066: }
1067: PetscCall(PetscDSSetJetDegree(prob, f, 1));
1068: }
1069: PetscFunctionReturn(PETSC_SUCCESS);
1070: }
1072: /*@
1073: PetscDSGetWeakForm - Returns the weak form object
1075: Not Collective
1077: Input Parameter:
1078: . ds - The `PetscDS` object
1080: Output Parameter:
1081: . wf - The weak form object
1083: Level: beginner
1085: .seealso: `PetscWeakForm`, `PetscDSSetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1086: @*/
1087: PetscErrorCode PetscDSGetWeakForm(PetscDS ds, PetscWeakForm *wf)
1088: {
1089: PetscFunctionBegin;
1091: PetscAssertPointer(wf, 2);
1092: *wf = ds->wf;
1093: PetscFunctionReturn(PETSC_SUCCESS);
1094: }
1096: /*@
1097: PetscDSSetWeakForm - Sets the weak form object
1099: Not Collective
1101: Input Parameters:
1102: + ds - The `PetscDS` object
1103: - wf - The weak form object
1105: Level: beginner
1107: .seealso: `PetscWeakForm`, `PetscDSGetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1108: @*/
1109: PetscErrorCode PetscDSSetWeakForm(PetscDS ds, PetscWeakForm wf)
1110: {
1111: PetscFunctionBegin;
1114: PetscCall(PetscObjectDereference((PetscObject)ds->wf));
1115: ds->wf = wf;
1116: PetscCall(PetscObjectReference((PetscObject)wf));
1117: PetscCall(PetscWeakFormSetNumFields(wf, ds->Nf));
1118: PetscFunctionReturn(PETSC_SUCCESS);
1119: }
1121: /*@
1122: PetscDSAddDiscretization - Adds a discretization object
1124: Not Collective
1126: Input Parameters:
1127: + prob - The `PetscDS` object
1128: - disc - The boundary discretization object
1130: Level: beginner
1132: .seealso: `PetscWeakForm`, `PetscDSGetDiscretization()`, `PetscDSSetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1133: @*/
1134: PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc)
1135: {
1136: PetscFunctionBegin;
1137: PetscCall(PetscDSSetDiscretization(prob, prob->Nf, disc));
1138: PetscFunctionReturn(PETSC_SUCCESS);
1139: }
1141: /*@
1142: PetscDSGetQuadrature - Returns the quadrature, which must agree for all fields in the `PetscDS`
1144: Not Collective
1146: Input Parameter:
1147: . prob - The `PetscDS` object
1149: Output Parameter:
1150: . q - The quadrature object
1152: Level: intermediate
1154: .seealso: `PetscDS`, `PetscQuadrature`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1155: @*/
1156: PetscErrorCode PetscDSGetQuadrature(PetscDS prob, PetscQuadrature *q)
1157: {
1158: PetscObject obj;
1159: PetscClassId id;
1161: PetscFunctionBegin;
1162: *q = NULL;
1163: if (!prob->Nf) PetscFunctionReturn(PETSC_SUCCESS);
1164: PetscCall(PetscDSGetDiscretization(prob, 0, &obj));
1165: PetscCall(PetscObjectGetClassId(obj, &id));
1166: if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetQuadrature((PetscFE)obj, q));
1167: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetQuadrature((PetscFV)obj, q));
1168: else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
1169: PetscFunctionReturn(PETSC_SUCCESS);
1170: }
1172: /*@
1173: PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for `TSIMEX`
1175: Not Collective
1177: Input Parameters:
1178: + prob - The `PetscDS` object
1179: - f - The field number
1181: Output Parameter:
1182: . implicit - The flag indicating what kind of solve to use for this field
1184: Level: developer
1186: .seealso: `TSIMEX`, `PetscDS`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1187: @*/
1188: PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit)
1189: {
1190: PetscFunctionBegin;
1192: PetscAssertPointer(implicit, 3);
1193: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
1194: *implicit = prob->implicit[f];
1195: PetscFunctionReturn(PETSC_SUCCESS);
1196: }
1198: /*@
1199: PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for `TSIMEX`
1201: Not Collective
1203: Input Parameters:
1204: + prob - The `PetscDS` object
1205: . f - The field number
1206: - implicit - The flag indicating what kind of solve to use for this field
1208: Level: developer
1210: .seealso: `TSIMEX`, `PetscDSGetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1211: @*/
1212: PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit)
1213: {
1214: PetscFunctionBegin;
1216: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
1217: prob->implicit[f] = implicit;
1218: PetscFunctionReturn(PETSC_SUCCESS);
1219: }
1221: /*@
1222: PetscDSGetJetDegree - Returns the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate.
1224: Not Collective
1226: Input Parameters:
1227: + ds - The `PetscDS` object
1228: - f - The field number
1230: Output Parameter:
1231: . k - The highest derivative we need to tabulate
1233: Level: developer
1235: .seealso: `PetscDS`, `PetscDSSetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1236: @*/
1237: PetscErrorCode PetscDSGetJetDegree(PetscDS ds, PetscInt f, PetscInt *k)
1238: {
1239: PetscFunctionBegin;
1241: PetscAssertPointer(k, 3);
1242: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1243: *k = ds->jetDegree[f];
1244: PetscFunctionReturn(PETSC_SUCCESS);
1245: }
1247: /*@
1248: PetscDSSetJetDegree - Set the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate.
1250: Not Collective
1252: Input Parameters:
1253: + ds - The `PetscDS` object
1254: . f - The field number
1255: - k - The highest derivative we need to tabulate
1257: Level: developer
1259: .seealso: ``PetscDS`, `PetscDSGetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1260: @*/
1261: PetscErrorCode PetscDSSetJetDegree(PetscDS ds, PetscInt f, PetscInt k)
1262: {
1263: PetscFunctionBegin;
1265: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1266: ds->jetDegree[f] = k;
1267: PetscFunctionReturn(PETSC_SUCCESS);
1268: }
1270: /*@C
1271: PetscDSGetObjective - Get the pointwise objective function for a given test field
1273: Not Collective
1275: Input Parameters:
1276: + ds - The `PetscDS`
1277: - f - The test field number
1279: Output Parameter:
1280: . obj - integrand for the test function term
1282: Calling sequence of `obj`:
1283: + dim - the spatial dimension
1284: . Nf - the number of fields
1285: . NfAux - the number of auxiliary fields
1286: . uOff - the offset into u[] and u_t[] for each field
1287: . uOff_x - the offset into u_x[] for each field
1288: . u - each field evaluated at the current point
1289: . u_t - the time derivative of each field evaluated at the current point
1290: . u_x - the gradient of each field evaluated at the current point
1291: . aOff - the offset into a[] and a_t[] for each auxiliary field
1292: . aOff_x - the offset into a_x[] for each auxiliary field
1293: . a - each auxiliary field evaluated at the current point
1294: . a_t - the time derivative of each auxiliary field evaluated at the current point
1295: . a_x - the gradient of auxiliary each field evaluated at the current point
1296: . t - current time
1297: . x - coordinates of the current point
1298: . numConstants - number of constant parameters
1299: . constants - constant parameters
1300: - obj - output values at the current point
1302: Level: intermediate
1304: Note:
1305: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi obj(u, u_t, \nabla u, x, t)$
1307: .seealso: `PetscDS`, `PetscDSSetObjective()`, `PetscDSGetResidual()`
1308: @*/
1309: PetscErrorCode PetscDSGetObjective(PetscDS ds, PetscInt f, void (**obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
1310: {
1311: PetscPointFunc *tmp;
1312: PetscInt n;
1314: PetscFunctionBegin;
1316: PetscAssertPointer(obj, 3);
1317: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1318: PetscCall(PetscWeakFormGetObjective(ds->wf, NULL, 0, f, 0, &n, &tmp));
1319: *obj = tmp ? tmp[0] : NULL;
1320: PetscFunctionReturn(PETSC_SUCCESS);
1321: }
1323: /*@C
1324: PetscDSSetObjective - Set the pointwise objective function for a given test field
1326: Not Collective
1328: Input Parameters:
1329: + ds - The `PetscDS`
1330: . f - The test field number
1331: - obj - integrand for the test function term
1333: Calling sequence of `obj`:
1334: + dim - the spatial dimension
1335: . Nf - the number of fields
1336: . NfAux - the number of auxiliary fields
1337: . uOff - the offset into u[] and u_t[] for each field
1338: . uOff_x - the offset into u_x[] for each field
1339: . u - each field evaluated at the current point
1340: . u_t - the time derivative of each field evaluated at the current point
1341: . u_x - the gradient of each field evaluated at the current point
1342: . aOff - the offset into a[] and a_t[] for each auxiliary field
1343: . aOff_x - the offset into a_x[] for each auxiliary field
1344: . a - each auxiliary field evaluated at the current point
1345: . a_t - the time derivative of each auxiliary field evaluated at the current point
1346: . a_x - the gradient of auxiliary each field evaluated at the current point
1347: . t - current time
1348: . x - coordinates of the current point
1349: . numConstants - number of constant parameters
1350: . constants - constant parameters
1351: - obj - output values at the current point
1353: Level: intermediate
1355: Note:
1356: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi obj(u, u_t, \nabla u, x, t)$
1358: .seealso: `PetscDS`, `PetscDSGetObjective()`, `PetscDSSetResidual()`
1359: @*/
1360: PetscErrorCode PetscDSSetObjective(PetscDS ds, PetscInt f, void (*obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
1361: {
1362: PetscFunctionBegin;
1365: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1366: PetscCall(PetscWeakFormSetIndexObjective(ds->wf, NULL, 0, f, 0, 0, obj));
1367: PetscFunctionReturn(PETSC_SUCCESS);
1368: }
1370: /*@C
1371: PetscDSGetResidual - Get the pointwise residual function for a given test field
1373: Not Collective
1375: Input Parameters:
1376: + ds - The `PetscDS`
1377: - f - The test field number
1379: Output Parameters:
1380: + f0 - integrand for the test function term
1381: - f1 - integrand for the test function gradient term
1383: Calling sequence of `f0`:
1384: + dim - the spatial dimension
1385: . Nf - the number of fields
1386: . NfAux - the number of auxiliary fields
1387: . uOff - the offset into u[] and u_t[] for each field
1388: . uOff_x - the offset into u_x[] for each field
1389: . u - each field evaluated at the current point
1390: . u_t - the time derivative of each field evaluated at the current point
1391: . u_x - the gradient of each field evaluated at the current point
1392: . aOff - the offset into a[] and a_t[] for each auxiliary field
1393: . aOff_x - the offset into a_x[] for each auxiliary field
1394: . a - each auxiliary field evaluated at the current point
1395: . a_t - the time derivative of each auxiliary field evaluated at the current point
1396: . a_x - the gradient of auxiliary each field evaluated at the current point
1397: . t - current time
1398: . x - coordinates of the current point
1399: . numConstants - number of constant parameters
1400: . constants - constant parameters
1401: - f0 - output values at the current point
1403: Level: intermediate
1405: Note:
1406: `f1` has an identical form and is omitted for brevity.
1408: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t)$
1410: .seealso: `PetscDS`, `PetscDSSetResidual()`
1411: @*/
1412: PetscErrorCode PetscDSGetResidual(PetscDS ds, PetscInt f, void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (**f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1413: {
1414: PetscPointFunc *tmp0, *tmp1;
1415: PetscInt n0, n1;
1417: PetscFunctionBegin;
1419: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1420: PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1));
1421: *f0 = tmp0 ? tmp0[0] : NULL;
1422: *f1 = tmp1 ? tmp1[0] : NULL;
1423: PetscFunctionReturn(PETSC_SUCCESS);
1424: }
1426: /*@C
1427: PetscDSSetResidual - Set the pointwise residual function for a given test field
1429: Not Collective
1431: Input Parameters:
1432: + ds - The `PetscDS`
1433: . f - The test field number
1434: . f0 - integrand for the test function term
1435: - f1 - integrand for the test function gradient term
1437: Calling sequence of `f0`:
1438: + dim - the spatial dimension
1439: . Nf - the number of fields
1440: . NfAux - the number of auxiliary fields
1441: . uOff - the offset into u[] and u_t[] for each field
1442: . uOff_x - the offset into u_x[] for each field
1443: . u - each field evaluated at the current point
1444: . u_t - the time derivative of each field evaluated at the current point
1445: . u_x - the gradient of each field evaluated at the current point
1446: . aOff - the offset into a[] and a_t[] for each auxiliary field
1447: . aOff_x - the offset into a_x[] for each auxiliary field
1448: . a - each auxiliary field evaluated at the current point
1449: . a_t - the time derivative of each auxiliary field evaluated at the current point
1450: . a_x - the gradient of auxiliary each field evaluated at the current point
1451: . t - current time
1452: . x - coordinates of the current point
1453: . numConstants - number of constant parameters
1454: . constants - constant parameters
1455: - f0 - output values at the current point
1457: Level: intermediate
1459: Note:
1460: `f1` has an identical form and is omitted for brevity.
1462: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t)$
1464: .seealso: `PetscDS`, `PetscDSGetResidual()`
1465: @*/
1466: PetscErrorCode PetscDSSetResidual(PetscDS ds, PetscInt f, void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (*f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1467: {
1468: PetscFunctionBegin;
1472: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1473: PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1));
1474: PetscFunctionReturn(PETSC_SUCCESS);
1475: }
1477: /*@C
1478: PetscDSGetRHSResidual - Get the pointwise RHS residual function for explicit timestepping for a given test field
1480: Not Collective
1482: Input Parameters:
1483: + ds - The `PetscDS`
1484: - f - The test field number
1486: Output Parameters:
1487: + f0 - integrand for the test function term
1488: - f1 - integrand for the test function gradient term
1490: Calling sequence of `f0`:
1491: + dim - the spatial dimension
1492: . Nf - the number of fields
1493: . NfAux - the number of auxiliary fields
1494: . uOff - the offset into u[] and u_t[] for each field
1495: . uOff_x - the offset into u_x[] for each field
1496: . u - each field evaluated at the current point
1497: . u_t - the time derivative of each field evaluated at the current point
1498: . u_x - the gradient of each field evaluated at the current point
1499: . aOff - the offset into a[] and a_t[] for each auxiliary field
1500: . aOff_x - the offset into a_x[] for each auxiliary field
1501: . a - each auxiliary field evaluated at the current point
1502: . a_t - the time derivative of each auxiliary field evaluated at the current point
1503: . a_x - the gradient of auxiliary each field evaluated at the current point
1504: . t - current time
1505: . x - coordinates of the current point
1506: . numConstants - number of constant parameters
1507: . constants - constant parameters
1508: - f0 - output values at the current point
1510: Level: intermediate
1512: Note:
1513: `f1` has an identical form and is omitted for brevity.
1515: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t)$
1517: .seealso: `PetscDS`, `PetscDSSetRHSResidual()`
1518: @*/
1519: PetscErrorCode PetscDSGetRHSResidual(PetscDS ds, PetscInt f, void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (**f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1520: {
1521: PetscPointFunc *tmp0, *tmp1;
1522: PetscInt n0, n1;
1524: PetscFunctionBegin;
1526: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1527: PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 100, &n0, &tmp0, &n1, &tmp1));
1528: *f0 = tmp0 ? tmp0[0] : NULL;
1529: *f1 = tmp1 ? tmp1[0] : NULL;
1530: PetscFunctionReturn(PETSC_SUCCESS);
1531: }
1533: /*@C
1534: PetscDSSetRHSResidual - Set the pointwise residual function for explicit timestepping for a given test field
1536: Not Collective
1538: Input Parameters:
1539: + ds - The `PetscDS`
1540: . f - The test field number
1541: . f0 - integrand for the test function term
1542: - f1 - integrand for the test function gradient term
1544: Calling sequence for the callbacks `f0`:
1545: + dim - the spatial dimension
1546: . Nf - the number of fields
1547: . NfAux - the number of auxiliary fields
1548: . uOff - the offset into u[] and u_t[] for each field
1549: . uOff_x - the offset into u_x[] for each field
1550: . u - each field evaluated at the current point
1551: . u_t - the time derivative of each field evaluated at the current point
1552: . u_x - the gradient of each field evaluated at the current point
1553: . aOff - the offset into a[] and a_t[] for each auxiliary field
1554: . aOff_x - the offset into a_x[] for each auxiliary field
1555: . a - each auxiliary field evaluated at the current point
1556: . a_t - the time derivative of each auxiliary field evaluated at the current point
1557: . a_x - the gradient of auxiliary each field evaluated at the current point
1558: . t - current time
1559: . x - coordinates of the current point
1560: . numConstants - number of constant parameters
1561: . constants - constant parameters
1562: - f0 - output values at the current point
1564: Level: intermediate
1566: Note:
1567: `f1` has an identical form and is omitted for brevity.
1569: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t)$
1571: .seealso: `PetscDS`, `PetscDSGetResidual()`
1572: @*/
1573: PetscErrorCode PetscDSSetRHSResidual(PetscDS ds, PetscInt f, void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]), void (*f1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1574: {
1575: PetscFunctionBegin;
1579: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1580: PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 100, 0, f0, 0, f1));
1581: PetscFunctionReturn(PETSC_SUCCESS);
1582: }
1584: /*@C
1585: PetscDSHasJacobian - Checks that the Jacobian functions have been set
1587: Not Collective
1589: Input Parameter:
1590: . ds - The `PetscDS`
1592: Output Parameter:
1593: . hasJac - flag that pointwise function for the Jacobian has been set
1595: Level: intermediate
1597: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1598: @*/
1599: PetscErrorCode PetscDSHasJacobian(PetscDS ds, PetscBool *hasJac)
1600: {
1601: PetscFunctionBegin;
1603: PetscCall(PetscWeakFormHasJacobian(ds->wf, hasJac));
1604: PetscFunctionReturn(PETSC_SUCCESS);
1605: }
1607: /*@C
1608: PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field
1610: Not Collective
1612: Input Parameters:
1613: + ds - The `PetscDS`
1614: . f - The test field number
1615: - g - The field number
1617: Output Parameters:
1618: + g0 - integrand for the test and basis function term
1619: . g1 - integrand for the test function and basis function gradient term
1620: . g2 - integrand for the test function gradient and basis function term
1621: - g3 - integrand for the test function gradient and basis function gradient term
1623: Calling sequence of `g0`:
1624: + dim - the spatial dimension
1625: . Nf - the number of fields
1626: . NfAux - the number of auxiliary fields
1627: . uOff - the offset into u[] and u_t[] for each field
1628: . uOff_x - the offset into u_x[] for each field
1629: . u - each field evaluated at the current point
1630: . u_t - the time derivative of each field evaluated at the current point
1631: . u_x - the gradient of each field evaluated at the current point
1632: . aOff - the offset into a[] and a_t[] for each auxiliary field
1633: . aOff_x - the offset into a_x[] for each auxiliary field
1634: . a - each auxiliary field evaluated at the current point
1635: . a_t - the time derivative of each auxiliary field evaluated at the current point
1636: . a_x - the gradient of auxiliary each field evaluated at the current point
1637: . t - current time
1638: . u_tShift - the multiplier a for dF/dU_t
1639: . x - coordinates of the current point
1640: . numConstants - number of constant parameters
1641: . constants - constant parameters
1642: - g0 - output values at the current point
1644: Level: intermediate
1646: Note:
1647: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1649: We are using a first order FEM model for the weak form\:
1651: $$
1652: \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
1653: $$
1655: .seealso: `PetscDS`, `PetscDSSetJacobian()`
1656: @*/
1657: PetscErrorCode PetscDSGetJacobian(PetscDS ds, PetscInt f, PetscInt g, void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux, const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[], const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[], PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]), void (**g1)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g2)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void (**g3)(PetscInt, PetscInt, PetscInt, const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[], PetscReal, PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]))
1658: {
1659: PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1660: PetscInt n0, n1, n2, n3;
1662: PetscFunctionBegin;
1664: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1665: PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf);
1666: PetscCall(PetscWeakFormGetJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1667: *g0 = tmp0 ? tmp0[0] : NULL;
1668: *g1 = tmp1 ? tmp1[0] : NULL;
1669: *g2 = tmp2 ? tmp2[0] : NULL;
1670: *g3 = tmp3 ? tmp3[0] : NULL;
1671: PetscFunctionReturn(PETSC_SUCCESS);
1672: }
1674: /*@C
1675: PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields
1677: Not Collective
1679: Input Parameters:
1680: + ds - The `PetscDS`
1681: . f - The test field number
1682: . g - The field number
1683: . g0 - integrand for the test and basis function term
1684: . g1 - integrand for the test function and basis function gradient term
1685: . g2 - integrand for the test function gradient and basis function term
1686: - g3 - integrand for the test function gradient and basis function gradient term
1688: Calling sequence of `g0`:
1689: + dim - the spatial dimension
1690: . Nf - the number of fields
1691: . NfAux - the number of auxiliary fields
1692: . uOff - the offset into u[] and u_t[] for each field
1693: . uOff_x - the offset into u_x[] for each field
1694: . u - each field evaluated at the current point
1695: . u_t - the time derivative of each field evaluated at the current point
1696: . u_x - the gradient of each field evaluated at the current point
1697: . aOff - the offset into a[] and a_t[] for each auxiliary field
1698: . aOff_x - the offset into a_x[] for each auxiliary field
1699: . a - each auxiliary field evaluated at the current point
1700: . a_t - the time derivative of each auxiliary field evaluated at the current point
1701: . a_x - the gradient of auxiliary each field evaluated at the current point
1702: . t - current time
1703: . u_tShift - the multiplier a for dF/dU_t
1704: . x - coordinates of the current point
1705: . numConstants - number of constant parameters
1706: . constants - constant parameters
1707: - g0 - output values at the current point
1709: Level: intermediate
1711: Note:
1712: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1714: We are using a first order FEM model for the weak form\:
1715: \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
1717: .seealso: `PetscDS`, `PetscDSGetJacobian()`
1718: @*/
1719: 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[]))
1720: {
1721: PetscFunctionBegin;
1727: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1728: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1729: PetscCall(PetscWeakFormSetIndexJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1730: PetscFunctionReturn(PETSC_SUCCESS);
1731: }
1733: /*@C
1734: PetscDSUseJacobianPreconditioner - Set whether to construct a Jacobian preconditioner
1736: Not Collective
1738: Input Parameters:
1739: + prob - The `PetscDS`
1740: - useJacPre - flag that enables construction of a Jacobian preconditioner
1742: Level: intermediate
1744: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1745: @*/
1746: PetscErrorCode PetscDSUseJacobianPreconditioner(PetscDS prob, PetscBool useJacPre)
1747: {
1748: PetscFunctionBegin;
1750: prob->useJacPre = useJacPre;
1751: PetscFunctionReturn(PETSC_SUCCESS);
1752: }
1754: /*@C
1755: PetscDSHasJacobianPreconditioner - Checks if a Jacobian preconditioner matrix has been set
1757: Not Collective
1759: Input Parameter:
1760: . ds - The `PetscDS`
1762: Output Parameter:
1763: . hasJacPre - flag that pointwise function for Jacobian preconditioner matrix has been set
1765: Level: intermediate
1767: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1768: @*/
1769: PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS ds, PetscBool *hasJacPre)
1770: {
1771: PetscFunctionBegin;
1773: *hasJacPre = PETSC_FALSE;
1774: if (!ds->useJacPre) PetscFunctionReturn(PETSC_SUCCESS);
1775: PetscCall(PetscWeakFormHasJacobianPreconditioner(ds->wf, hasJacPre));
1776: PetscFunctionReturn(PETSC_SUCCESS);
1777: }
1779: /*@C
1780: PetscDSGetJacobianPreconditioner - Get the pointwise Jacobian preconditioner function for given test and basis field. If this is missing,
1781: the system matrix is used to build the preconditioner.
1783: Not Collective
1785: Input Parameters:
1786: + ds - The `PetscDS`
1787: . f - The test field number
1788: - g - The field number
1790: Output Parameters:
1791: + g0 - integrand for the test and basis function term
1792: . g1 - integrand for the test function and basis function gradient term
1793: . g2 - integrand for the test function gradient and basis function term
1794: - g3 - integrand for the test function gradient and basis function gradient term
1796: Calling sequence of `g0`:
1797: + dim - the spatial dimension
1798: . Nf - the number of fields
1799: . NfAux - the number of auxiliary fields
1800: . uOff - the offset into u[] and u_t[] for each field
1801: . uOff_x - the offset into u_x[] for each field
1802: . u - each field evaluated at the current point
1803: . u_t - the time derivative of each field evaluated at the current point
1804: . u_x - the gradient of each field evaluated at the current point
1805: . aOff - the offset into a[] and a_t[] for each auxiliary field
1806: . aOff_x - the offset into a_x[] for each auxiliary field
1807: . a - each auxiliary field evaluated at the current point
1808: . a_t - the time derivative of each auxiliary field evaluated at the current point
1809: . a_x - the gradient of auxiliary each field evaluated at the current point
1810: . t - current time
1811: . u_tShift - the multiplier a for dF/dU_t
1812: . x - coordinates of the current point
1813: . numConstants - number of constant parameters
1814: . constants - constant parameters
1815: - g0 - output values at the current point
1817: Level: intermediate
1819: Note:
1820: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1821: We are using a first order FEM model for the weak form\:
1822: \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
1824: .seealso: `PetscDS`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1825: @*/
1826: 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[]))
1827: {
1828: PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1829: PetscInt n0, n1, n2, n3;
1831: PetscFunctionBegin;
1833: 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);
1834: 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);
1835: PetscCall(PetscWeakFormGetJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1836: *g0 = tmp0 ? tmp0[0] : NULL;
1837: *g1 = tmp1 ? tmp1[0] : NULL;
1838: *g2 = tmp2 ? tmp2[0] : NULL;
1839: *g3 = tmp3 ? tmp3[0] : NULL;
1840: PetscFunctionReturn(PETSC_SUCCESS);
1841: }
1843: /*@C
1844: PetscDSSetJacobianPreconditioner - Set the pointwise Jacobian preconditioner function for given test and basis fields.
1845: If this is missing, the system matrix is used to build the preconditioner.
1847: Not Collective
1849: Input Parameters:
1850: + ds - The `PetscDS`
1851: . f - The test field number
1852: . g - The field number
1853: . g0 - integrand for the test and basis function term
1854: . g1 - integrand for the test function and basis function gradient term
1855: . g2 - integrand for the test function gradient and basis function term
1856: - g3 - integrand for the test function gradient and basis function gradient term
1858: Calling sequence of `g0`:
1859: + dim - the spatial dimension
1860: . Nf - the number of fields
1861: . NfAux - the number of auxiliary fields
1862: . uOff - the offset into u[] and u_t[] for each field
1863: . uOff_x - the offset into u_x[] for each field
1864: . u - each field evaluated at the current point
1865: . u_t - the time derivative of each field evaluated at the current point
1866: . u_x - the gradient of each field evaluated at the current point
1867: . aOff - the offset into a[] and a_t[] for each auxiliary field
1868: . aOff_x - the offset into a_x[] for each auxiliary field
1869: . a - each auxiliary field evaluated at the current point
1870: . a_t - the time derivative of each auxiliary field evaluated at the current point
1871: . a_x - the gradient of auxiliary each field evaluated at the current point
1872: . t - current time
1873: . u_tShift - the multiplier a for dF/dU_t
1874: . x - coordinates of the current point
1875: . numConstants - number of constant parameters
1876: . constants - constant parameters
1877: - g0 - output values at the current point
1879: Level: intermediate
1881: Note:
1882: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1884: We are using a first order FEM model for the weak form\:
1885: \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
1887: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobian()`
1888: @*/
1889: 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[]))
1890: {
1891: PetscFunctionBegin;
1897: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1898: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1899: PetscCall(PetscWeakFormSetIndexJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1900: PetscFunctionReturn(PETSC_SUCCESS);
1901: }
1903: /*@C
1904: PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, dF/du_t, has been set
1906: Not Collective
1908: Input Parameter:
1909: . ds - The `PetscDS`
1911: Output Parameter:
1912: . hasDynJac - flag that pointwise function for dynamic Jacobian has been set
1914: Level: intermediate
1916: .seealso: `PetscDS`, `PetscDSGetDynamicJacobian()`, `PetscDSSetDynamicJacobian()`, `PetscDSGetJacobian()`
1917: @*/
1918: PetscErrorCode PetscDSHasDynamicJacobian(PetscDS ds, PetscBool *hasDynJac)
1919: {
1920: PetscFunctionBegin;
1922: PetscCall(PetscWeakFormHasDynamicJacobian(ds->wf, hasDynJac));
1923: PetscFunctionReturn(PETSC_SUCCESS);
1924: }
1926: /*@C
1927: PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, dF/du_t, function for given test and basis field
1929: Not Collective
1931: Input Parameters:
1932: + ds - The `PetscDS`
1933: . f - The test field number
1934: - g - The field number
1936: Output Parameters:
1937: + g0 - integrand for the test and basis function term
1938: . g1 - integrand for the test function and basis function gradient term
1939: . g2 - integrand for the test function gradient and basis function term
1940: - g3 - integrand for the test function gradient and basis function gradient term
1942: Calling sequence of `g0`:
1943: + dim - the spatial dimension
1944: . Nf - the number of fields
1945: . NfAux - the number of auxiliary fields
1946: . uOff - the offset into u[] and u_t[] for each field
1947: . uOff_x - the offset into u_x[] for each field
1948: . u - each field evaluated at the current point
1949: . u_t - the time derivative of each field evaluated at the current point
1950: . u_x - the gradient of each field evaluated at the current point
1951: . aOff - the offset into a[] and a_t[] for each auxiliary field
1952: . aOff_x - the offset into a_x[] for each auxiliary field
1953: . a - each auxiliary field evaluated at the current point
1954: . a_t - the time derivative of each auxiliary field evaluated at the current point
1955: . a_x - the gradient of auxiliary each field evaluated at the current point
1956: . t - current time
1957: . u_tShift - the multiplier a for dF/dU_t
1958: . x - coordinates of the current point
1959: . numConstants - number of constant parameters
1960: . constants - constant parameters
1961: - g0 - output values at the current point
1963: Level: intermediate
1965: Note:
1966: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
1968: We are using a first order FEM model for the weak form\:
1969: \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
1971: .seealso: `PetscDS`, `PetscDSSetJacobian()`
1972: @*/
1973: 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[]))
1974: {
1975: PetscPointJac *tmp0, *tmp1, *tmp2, *tmp3;
1976: PetscInt n0, n1, n2, n3;
1978: PetscFunctionBegin;
1980: 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);
1981: 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);
1982: PetscCall(PetscWeakFormGetDynamicJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1983: *g0 = tmp0 ? tmp0[0] : NULL;
1984: *g1 = tmp1 ? tmp1[0] : NULL;
1985: *g2 = tmp2 ? tmp2[0] : NULL;
1986: *g3 = tmp3 ? tmp3[0] : NULL;
1987: PetscFunctionReturn(PETSC_SUCCESS);
1988: }
1990: /*@C
1991: PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, dF/du_t, function for given test and basis fields
1993: Not Collective
1995: Input Parameters:
1996: + ds - The `PetscDS`
1997: . f - The test field number
1998: . g - The field number
1999: . g0 - integrand for the test and basis function term
2000: . g1 - integrand for the test function and basis function gradient term
2001: . g2 - integrand for the test function gradient and basis function term
2002: - g3 - integrand for the test function gradient and basis function gradient term
2004: Calling sequence of `g0`:
2005: + dim - the spatial dimension
2006: . Nf - the number of fields
2007: . NfAux - the number of auxiliary fields
2008: . uOff - the offset into u[] and u_t[] for each field
2009: . uOff_x - the offset into u_x[] for each field
2010: . u - each field evaluated at the current point
2011: . u_t - the time derivative of each field evaluated at the current point
2012: . u_x - the gradient of each field evaluated at the current point
2013: . aOff - the offset into a[] and a_t[] for each auxiliary field
2014: . aOff_x - the offset into a_x[] for each auxiliary field
2015: . a - each auxiliary field evaluated at the current point
2016: . a_t - the time derivative of each auxiliary field evaluated at the current point
2017: . a_x - the gradient of auxiliary each field evaluated at the current point
2018: . t - current time
2019: . u_tShift - the multiplier a for dF/dU_t
2020: . x - coordinates of the current point
2021: . numConstants - number of constant parameters
2022: . constants - constant parameters
2023: - g0 - output values at the current point
2025: Level: intermediate
2027: Note:
2028: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2030: We are using a first order FEM model for the weak form\:
2031: \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
2033: .seealso: `PetscDS`, `PetscDSGetJacobian()`
2034: @*/
2035: 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[]))
2036: {
2037: PetscFunctionBegin;
2043: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2044: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2045: PetscCall(PetscWeakFormSetIndexDynamicJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2046: PetscFunctionReturn(PETSC_SUCCESS);
2047: }
2049: /*@C
2050: PetscDSGetRiemannSolver - Returns the Riemann solver for the given field
2052: Not Collective
2054: Input Parameters:
2055: + ds - The `PetscDS` object
2056: - f - The field number
2058: Output Parameter:
2059: . r - Riemann solver
2061: Calling sequence of `r`:
2062: + dim - The spatial dimension
2063: . Nf - The number of fields
2064: . x - The coordinates at a point on the interface
2065: . n - The normal vector to the interface
2066: . uL - The state vector to the left of the interface
2067: . uR - The state vector to the right of the interface
2068: . flux - output array of flux through the interface
2069: . numConstants - number of constant parameters
2070: . constants - constant parameters
2071: - ctx - optional user context
2073: Level: intermediate
2075: .seealso: `PetscDS`, `PetscDSSetRiemannSolver()`
2076: @*/
2077: 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))
2078: {
2079: PetscRiemannFunc *tmp;
2080: PetscInt n;
2082: PetscFunctionBegin;
2084: PetscAssertPointer(r, 3);
2085: 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);
2086: PetscCall(PetscWeakFormGetRiemannSolver(ds->wf, NULL, 0, f, 0, &n, &tmp));
2087: *r = tmp ? tmp[0] : NULL;
2088: PetscFunctionReturn(PETSC_SUCCESS);
2089: }
2091: /*@C
2092: PetscDSSetRiemannSolver - Sets the Riemann solver for the given field
2094: Not Collective
2096: Input Parameters:
2097: + ds - The `PetscDS` object
2098: . f - The field number
2099: - r - Riemann solver
2101: Calling sequence of `r`:
2102: + dim - The spatial dimension
2103: . Nf - The number of fields
2104: . x - The coordinates at a point on the interface
2105: . n - The normal vector to the interface
2106: . uL - The state vector to the left of the interface
2107: . uR - The state vector to the right of the interface
2108: . flux - output array of flux through the interface
2109: . numConstants - number of constant parameters
2110: . constants - constant parameters
2111: - ctx - optional user context
2113: Level: intermediate
2115: .seealso: `PetscDS`, `PetscDSGetRiemannSolver()`
2116: @*/
2117: 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))
2118: {
2119: PetscFunctionBegin;
2122: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2123: PetscCall(PetscWeakFormSetIndexRiemannSolver(ds->wf, NULL, 0, f, 0, 0, r));
2124: PetscFunctionReturn(PETSC_SUCCESS);
2125: }
2127: /*@C
2128: PetscDSGetUpdate - Get the pointwise update function for a given field
2130: Not Collective
2132: Input Parameters:
2133: + ds - The `PetscDS`
2134: - f - The field number
2136: Output Parameter:
2137: . update - update function
2139: Calling sequence of `update`:
2140: + dim - the spatial dimension
2141: . Nf - the number of fields
2142: . NfAux - the number of auxiliary fields
2143: . uOff - the offset into u[] and u_t[] for each field
2144: . uOff_x - the offset into u_x[] for each field
2145: . u - each field evaluated at the current point
2146: . u_t - the time derivative of each field evaluated at the current point
2147: . u_x - the gradient of each field evaluated at the current point
2148: . aOff - the offset into a[] and a_t[] for each auxiliary field
2149: . aOff_x - the offset into a_x[] for each auxiliary field
2150: . a - each auxiliary field evaluated at the current point
2151: . a_t - the time derivative of each auxiliary field evaluated at the current point
2152: . a_x - the gradient of auxiliary each field evaluated at the current point
2153: . t - current time
2154: . x - coordinates of the current point
2155: . numConstants - number of constant parameters
2156: . constants - constant parameters
2157: - uNew - new value for field at the current point
2159: Level: intermediate
2161: .seealso: `PetscDS`, `PetscDSSetUpdate()`, `PetscDSSetResidual()`
2162: @*/
2163: 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[]))
2164: {
2165: PetscFunctionBegin;
2167: 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);
2168: if (update) {
2169: PetscAssertPointer(update, 3);
2170: *update = ds->update[f];
2171: }
2172: PetscFunctionReturn(PETSC_SUCCESS);
2173: }
2175: /*@C
2176: PetscDSSetUpdate - Set the pointwise update function for a given field
2178: Not Collective
2180: Input Parameters:
2181: + ds - The `PetscDS`
2182: . f - The field number
2183: - update - update function
2185: Calling sequence of `update`:
2186: + dim - the spatial dimension
2187: . Nf - the number of fields
2188: . NfAux - the number of auxiliary fields
2189: . uOff - the offset into u[] and u_t[] for each field
2190: . uOff_x - the offset into u_x[] for each field
2191: . u - each field evaluated at the current point
2192: . u_t - the time derivative of each field evaluated at the current point
2193: . u_x - the gradient of each field evaluated at the current point
2194: . aOff - the offset into a[] and a_t[] for each auxiliary field
2195: . aOff_x - the offset into a_x[] for each auxiliary field
2196: . a - each auxiliary field evaluated at the current point
2197: . a_t - the time derivative of each auxiliary field evaluated at the current point
2198: . a_x - the gradient of auxiliary each field evaluated at the current point
2199: . t - current time
2200: . x - coordinates of the current point
2201: . numConstants - number of constant parameters
2202: . constants - constant parameters
2203: - uNew - new field values at the current point
2205: Level: intermediate
2207: .seealso: `PetscDS`, `PetscDSGetResidual()`
2208: @*/
2209: 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[]))
2210: {
2211: PetscFunctionBegin;
2214: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2215: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
2216: ds->update[f] = update;
2217: PetscFunctionReturn(PETSC_SUCCESS);
2218: }
2220: PetscErrorCode PetscDSGetContext(PetscDS ds, PetscInt f, void *ctx)
2221: {
2222: PetscFunctionBegin;
2224: 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);
2225: PetscAssertPointer(ctx, 3);
2226: *(void **)ctx = ds->ctx[f];
2227: PetscFunctionReturn(PETSC_SUCCESS);
2228: }
2230: PetscErrorCode PetscDSSetContext(PetscDS ds, PetscInt f, void *ctx)
2231: {
2232: PetscFunctionBegin;
2234: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2235: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
2236: ds->ctx[f] = ctx;
2237: PetscFunctionReturn(PETSC_SUCCESS);
2238: }
2240: /*@C
2241: PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field
2243: Not Collective
2245: Input Parameters:
2246: + ds - The PetscDS
2247: - f - The test field number
2249: Output Parameters:
2250: + f0 - boundary integrand for the test function term
2251: - f1 - boundary integrand for the test function gradient term
2253: Calling sequence of `f0`:
2254: + dim - the spatial dimension
2255: . Nf - the number of fields
2256: . NfAux - the number of auxiliary fields
2257: . uOff - the offset into u[] and u_t[] for each field
2258: . uOff_x - the offset into u_x[] for each field
2259: . u - each field evaluated at the current point
2260: . u_t - the time derivative of each field evaluated at the current point
2261: . u_x - the gradient of each field evaluated at the current point
2262: . aOff - the offset into a[] and a_t[] for each auxiliary field
2263: . aOff_x - the offset into a_x[] for each auxiliary field
2264: . a - each auxiliary field evaluated at the current point
2265: . a_t - the time derivative of each auxiliary field evaluated at the current point
2266: . a_x - the gradient of auxiliary each field evaluated at the current point
2267: . t - current time
2268: . x - coordinates of the current point
2269: . n - unit normal at the current point
2270: . numConstants - number of constant parameters
2271: . constants - constant parameters
2272: - f0 - output values at the current point
2274: Level: intermediate
2276: Note:
2277: The calling sequence of `f1` is identical, and therefore omitted for brevity.
2279: We are using a first order FEM model for the weak form\:
2280: \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
2282: .seealso: `PetscDS`, `PetscDSSetBdResidual()`
2283: @*/
2284: 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[]))
2285: {
2286: PetscBdPointFunc *tmp0, *tmp1;
2287: PetscInt n0, n1;
2289: PetscFunctionBegin;
2291: 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);
2292: PetscCall(PetscWeakFormGetBdResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1));
2293: *f0 = tmp0 ? tmp0[0] : NULL;
2294: *f1 = tmp1 ? tmp1[0] : NULL;
2295: PetscFunctionReturn(PETSC_SUCCESS);
2296: }
2298: /*@C
2299: PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field
2301: Not Collective
2303: Input Parameters:
2304: + ds - The `PetscDS`
2305: . f - The test field number
2306: . f0 - boundary integrand for the test function term
2307: - f1 - boundary integrand for the test function gradient term
2309: Calling sequence of `f0`:
2310: + dim - the spatial dimension
2311: . Nf - the number of fields
2312: . NfAux - the number of auxiliary fields
2313: . uOff - the offset into u[] and u_t[] for each field
2314: . uOff_x - the offset into u_x[] for each field
2315: . u - each field evaluated at the current point
2316: . u_t - the time derivative of each field evaluated at the current point
2317: . u_x - the gradient of each field evaluated at the current point
2318: . aOff - the offset into a[] and a_t[] for each auxiliary field
2319: . aOff_x - the offset into a_x[] for each auxiliary field
2320: . a - each auxiliary field evaluated at the current point
2321: . a_t - the time derivative of each auxiliary field evaluated at the current point
2322: . a_x - the gradient of auxiliary each field evaluated at the current point
2323: . t - current time
2324: . x - coordinates of the current point
2325: . n - unit normal at the current point
2326: . numConstants - number of constant parameters
2327: . constants - constant parameters
2328: - f0 - output values at the current point
2330: Level: intermediate
2332: Note:
2333: The calling sequence of `f1` is identical, and therefore omitted for brevity.
2335: We are using a first order FEM model for the weak form\:
2336: \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
2338: .seealso: `PetscDS`, `PetscDSGetBdResidual()`
2339: @*/
2340: 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[]))
2341: {
2342: PetscFunctionBegin;
2344: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2345: PetscCall(PetscWeakFormSetIndexBdResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1));
2346: PetscFunctionReturn(PETSC_SUCCESS);
2347: }
2349: /*@
2350: PetscDSHasBdJacobian - Indicates that boundary Jacobian functions have been set
2352: Not Collective
2354: Input Parameter:
2355: . ds - The `PetscDS`
2357: Output Parameter:
2358: . hasBdJac - flag that pointwise function for the boundary Jacobian has been set
2360: Level: intermediate
2362: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2363: @*/
2364: PetscErrorCode PetscDSHasBdJacobian(PetscDS ds, PetscBool *hasBdJac)
2365: {
2366: PetscFunctionBegin;
2368: PetscAssertPointer(hasBdJac, 2);
2369: PetscCall(PetscWeakFormHasBdJacobian(ds->wf, hasBdJac));
2370: PetscFunctionReturn(PETSC_SUCCESS);
2371: }
2373: /*@C
2374: PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field
2376: Not Collective
2378: Input Parameters:
2379: + ds - The `PetscDS`
2380: . f - The test field number
2381: - g - The field number
2383: Output Parameters:
2384: + g0 - integrand for the test and basis function term
2385: . g1 - integrand for the test function and basis function gradient term
2386: . g2 - integrand for the test function gradient and basis function term
2387: - g3 - integrand for the test function gradient and basis function gradient term
2389: Calling sequence of `g0`:
2390: + dim - the spatial dimension
2391: . Nf - the number of fields
2392: . NfAux - the number of auxiliary fields
2393: . uOff - the offset into u[] and u_t[] for each field
2394: . uOff_x - the offset into u_x[] for each field
2395: . u - each field evaluated at the current point
2396: . u_t - the time derivative of each field evaluated at the current point
2397: . u_x - the gradient of each field evaluated at the current point
2398: . aOff - the offset into a[] and a_t[] for each auxiliary field
2399: . aOff_x - the offset into a_x[] for each auxiliary field
2400: . a - each auxiliary field evaluated at the current point
2401: . a_t - the time derivative of each auxiliary field evaluated at the current point
2402: . a_x - the gradient of auxiliary each field evaluated at the current point
2403: . t - current time
2404: . u_tShift - the multiplier a for dF/dU_t
2405: . x - coordinates of the current point
2406: . n - normal at the current point
2407: . numConstants - number of constant parameters
2408: . constants - constant parameters
2409: - g0 - output values at the current point
2411: Level: intermediate
2413: Note:
2414: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2416: We are using a first order FEM model for the weak form\:
2417: \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
2419: .seealso: `PetscDS`, `PetscDSSetBdJacobian()`
2420: @*/
2421: 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[]))
2422: {
2423: PetscBdPointJac *tmp0, *tmp1, *tmp2, *tmp3;
2424: PetscInt n0, n1, n2, n3;
2426: PetscFunctionBegin;
2428: 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);
2429: 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);
2430: PetscCall(PetscWeakFormGetBdJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2431: *g0 = tmp0 ? tmp0[0] : NULL;
2432: *g1 = tmp1 ? tmp1[0] : NULL;
2433: *g2 = tmp2 ? tmp2[0] : NULL;
2434: *g3 = tmp3 ? tmp3[0] : NULL;
2435: PetscFunctionReturn(PETSC_SUCCESS);
2436: }
2438: /*@C
2439: PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field
2441: Not Collective
2443: Input Parameters:
2444: + ds - The PetscDS
2445: . f - The test field number
2446: . g - The field number
2447: . g0 - integrand for the test and basis function term
2448: . g1 - integrand for the test function and basis function gradient term
2449: . g2 - integrand for the test function gradient and basis function term
2450: - g3 - integrand for the test function gradient and basis function gradient term
2452: Calling sequence of `g0`:
2453: + dim - the spatial dimension
2454: . Nf - the number of fields
2455: . NfAux - the number of auxiliary fields
2456: . uOff - the offset into u[] and u_t[] for each field
2457: . uOff_x - the offset into u_x[] for each field
2458: . u - each field evaluated at the current point
2459: . u_t - the time derivative of each field evaluated at the current point
2460: . u_x - the gradient of each field evaluated at the current point
2461: . aOff - the offset into a[] and a_t[] for each auxiliary field
2462: . aOff_x - the offset into a_x[] for each auxiliary field
2463: . a - each auxiliary field evaluated at the current point
2464: . a_t - the time derivative of each auxiliary field evaluated at the current point
2465: . a_x - the gradient of auxiliary each field evaluated at the current point
2466: . t - current time
2467: . u_tShift - the multiplier a for dF/dU_t
2468: . x - coordinates of the current point
2469: . n - normal at the current point
2470: . numConstants - number of constant parameters
2471: . constants - constant parameters
2472: - g0 - output values at the current point
2474: Level: intermediate
2476: Note:
2477: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2479: We are using a first order FEM model for the weak form\:
2480: \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
2482: .seealso: `PetscDS`, `PetscDSGetBdJacobian()`
2483: @*/
2484: 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[]))
2485: {
2486: PetscFunctionBegin;
2492: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2493: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2494: PetscCall(PetscWeakFormSetIndexBdJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2495: PetscFunctionReturn(PETSC_SUCCESS);
2496: }
2498: /*@
2499: PetscDSHasBdJacobianPreconditioner - Signals that boundary Jacobian preconditioner functions have been set
2501: Not Collective
2503: Input Parameter:
2504: . ds - The `PetscDS`
2506: Output Parameter:
2507: . hasBdJacPre - flag that pointwise function for the boundary Jacobian preconditioner has been set
2509: Level: intermediate
2511: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2512: @*/
2513: PetscErrorCode PetscDSHasBdJacobianPreconditioner(PetscDS ds, PetscBool *hasBdJacPre)
2514: {
2515: PetscFunctionBegin;
2517: PetscAssertPointer(hasBdJacPre, 2);
2518: PetscCall(PetscWeakFormHasBdJacobianPreconditioner(ds->wf, hasBdJacPre));
2519: PetscFunctionReturn(PETSC_SUCCESS);
2520: }
2522: /*@C
2523: PetscDSGetBdJacobianPreconditioner - Get the pointwise boundary Jacobian preconditioner function for given test and basis field
2525: Not Collective; No Fortran Support
2527: Input Parameters:
2528: + ds - The `PetscDS`
2529: . f - The test field number
2530: - g - The field number
2532: Output Parameters:
2533: + g0 - integrand for the test and basis function term
2534: . g1 - integrand for the test function and basis function gradient term
2535: . g2 - integrand for the test function gradient and basis function term
2536: - g3 - integrand for the test function gradient and basis function gradient term
2538: Calling sequence of `g0`:
2539: + dim - the spatial dimension
2540: . Nf - the number of fields
2541: . NfAux - the number of auxiliary fields
2542: . uOff - the offset into u[] and u_t[] for each field
2543: . uOff_x - the offset into u_x[] for each field
2544: . u - each field evaluated at the current point
2545: . u_t - the time derivative of each field evaluated at the current point
2546: . u_x - the gradient of each field evaluated at the current point
2547: . aOff - the offset into a[] and a_t[] for each auxiliary field
2548: . aOff_x - the offset into a_x[] for each auxiliary field
2549: . a - each auxiliary field evaluated at the current point
2550: . a_t - the time derivative of each auxiliary field evaluated at the current point
2551: . a_x - the gradient of auxiliary each field evaluated at the current point
2552: . t - current time
2553: . u_tShift - the multiplier a for dF/dU_t
2554: . x - coordinates of the current point
2555: . n - normal at the current point
2556: . numConstants - number of constant parameters
2557: . constants - constant parameters
2558: - g0 - output values at the current point
2560: Level: intermediate
2562: Note:
2563: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2565: We are using a first order FEM model for the weak form\:
2566: \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
2568: .seealso: `PetscDS`, `PetscDSSetBdJacobianPreconditioner()`
2569: @*/
2570: 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[]))
2571: {
2572: PetscBdPointJac *tmp0, *tmp1, *tmp2, *tmp3;
2573: PetscInt n0, n1, n2, n3;
2575: PetscFunctionBegin;
2577: 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);
2578: 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);
2579: PetscCall(PetscWeakFormGetBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2580: *g0 = tmp0 ? tmp0[0] : NULL;
2581: *g1 = tmp1 ? tmp1[0] : NULL;
2582: *g2 = tmp2 ? tmp2[0] : NULL;
2583: *g3 = tmp3 ? tmp3[0] : NULL;
2584: PetscFunctionReturn(PETSC_SUCCESS);
2585: }
2587: /*@C
2588: PetscDSSetBdJacobianPreconditioner - Set the pointwise boundary Jacobian preconditioner function for given test and basis field
2590: Not Collective; No Fortran Support
2592: Input Parameters:
2593: + ds - The `PetscDS`
2594: . f - The test field number
2595: . g - The field number
2596: . g0 - integrand for the test and basis function term
2597: . g1 - integrand for the test function and basis function gradient term
2598: . g2 - integrand for the test function gradient and basis function term
2599: - g3 - integrand for the test function gradient and basis function gradient term
2601: Calling sequence of `g0':
2602: + dim - the spatial dimension
2603: . Nf - the number of fields
2604: . NfAux - the number of auxiliary fields
2605: . uOff - the offset into u[] and u_t[] for each field
2606: . uOff_x - the offset into u_x[] for each field
2607: . u - each field evaluated at the current point
2608: . u_t - the time derivative of each field evaluated at the current point
2609: . u_x - the gradient of each field evaluated at the current point
2610: . aOff - the offset into a[] and a_t[] for each auxiliary field
2611: . aOff_x - the offset into a_x[] for each auxiliary field
2612: . a - each auxiliary field evaluated at the current point
2613: . a_t - the time derivative of each auxiliary field evaluated at the current point
2614: . a_x - the gradient of auxiliary each field evaluated at the current point
2615: . t - current time
2616: . u_tShift - the multiplier a for dF/dU_t
2617: . x - coordinates of the current point
2618: . n - normal at the current point
2619: . numConstants - number of constant parameters
2620: . constants - constant parameters
2621: - g0 - output values at the current point
2623: Level: intermediate
2625: Note:
2626: `g1`, `g2`, and `g3` have identical calling sequences to `g0` and are omitted for brevity.
2628: We are using a first order FEM model for the weak form\:
2629: \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
2631: .seealso: `PetscDS`, `PetscDSGetBdJacobianPreconditioner()`
2632: @*/
2633: 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[]))
2634: {
2635: PetscFunctionBegin;
2641: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2642: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2643: PetscCall(PetscWeakFormSetIndexBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2644: PetscFunctionReturn(PETSC_SUCCESS);
2645: }
2647: /*@C
2648: PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field
2650: Not Collective
2652: Input Parameters:
2653: + prob - The PetscDS
2654: - f - The test field number
2656: Output Parameters:
2657: + sol - exact solution for the test field
2658: - ctx - exact solution context
2660: Calling sequence of `exactSol`:
2661: + dim - the spatial dimension
2662: . t - current time
2663: . x - coordinates of the current point
2664: . Nc - the number of field components
2665: . u - the solution field evaluated at the current point
2666: - ctx - a user context
2668: Level: intermediate
2670: .seealso: `PetscDS`, `PetscDSSetExactSolution()`, `PetscDSGetExactSolutionTimeDerivative()`
2671: @*/
2672: PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx)
2673: {
2674: PetscFunctionBegin;
2676: 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);
2677: if (sol) {
2678: PetscAssertPointer(sol, 3);
2679: *sol = prob->exactSol[f];
2680: }
2681: if (ctx) {
2682: PetscAssertPointer(ctx, 4);
2683: *ctx = prob->exactCtx[f];
2684: }
2685: PetscFunctionReturn(PETSC_SUCCESS);
2686: }
2688: /*@C
2689: PetscDSSetExactSolution - Set the pointwise exact solution function for a given test field
2691: Not Collective
2693: Input Parameters:
2694: + prob - The `PetscDS`
2695: . f - The test field number
2696: . sol - solution function for the test fields
2697: - ctx - solution context or `NULL`
2699: Calling sequence of `sol`:
2700: + dim - the spatial dimension
2701: . t - current time
2702: . x - coordinates of the current point
2703: . Nc - the number of field components
2704: . u - the solution field evaluated at the current point
2705: - ctx - a user context
2707: Level: intermediate
2709: .seealso: `PetscDS`, `PetscDSGetExactSolution()`
2710: @*/
2711: PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx)
2712: {
2713: PetscFunctionBegin;
2715: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2716: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
2717: if (sol) {
2719: prob->exactSol[f] = sol;
2720: }
2721: if (ctx) {
2723: prob->exactCtx[f] = ctx;
2724: }
2725: PetscFunctionReturn(PETSC_SUCCESS);
2726: }
2728: /*@C
2729: PetscDSGetExactSolutionTimeDerivative - Get the pointwise time derivative of the exact solution function for a given test field
2731: Not Collective
2733: Input Parameters:
2734: + prob - The `PetscDS`
2735: - f - The test field number
2737: Output Parameters:
2738: + sol - time derivative of the exact solution for the test field
2739: - ctx - time derivative of the exact solution context
2741: Calling sequence of `exactSol`:
2742: + dim - the spatial dimension
2743: . t - current time
2744: . x - coordinates of the current point
2745: . Nc - the number of field components
2746: . u - the solution field evaluated at the current point
2747: - ctx - a user context
2749: Level: intermediate
2751: .seealso: `PetscDS`, `PetscDSSetExactSolutionTimeDerivative()`, `PetscDSGetExactSolution()`
2752: @*/
2753: PetscErrorCode PetscDSGetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx)
2754: {
2755: PetscFunctionBegin;
2757: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
2758: if (sol) {
2759: PetscAssertPointer(sol, 3);
2760: *sol = prob->exactSol_t[f];
2761: }
2762: if (ctx) {
2763: PetscAssertPointer(ctx, 4);
2764: *ctx = prob->exactCtx_t[f];
2765: }
2766: PetscFunctionReturn(PETSC_SUCCESS);
2767: }
2769: /*@C
2770: PetscDSSetExactSolutionTimeDerivative - Set the pointwise time derivative of the exact solution function for a given test field
2772: Not Collective
2774: Input Parameters:
2775: + prob - The `PetscDS`
2776: . f - The test field number
2777: . sol - time derivative of the solution function for the test fields
2778: - ctx - time derivative of the solution context or `NULL`
2780: Calling sequence of `sol`:
2781: + dim - the spatial dimension
2782: . t - current time
2783: . x - coordinates of the current point
2784: . Nc - the number of field components
2785: . u - the solution field evaluated at the current point
2786: - ctx - a user context
2788: Level: intermediate
2790: .seealso: `PetscDS`, `PetscDSGetExactSolutionTimeDerivative()`, `PetscDSSetExactSolution()`
2791: @*/
2792: PetscErrorCode PetscDSSetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx)
2793: {
2794: PetscFunctionBegin;
2796: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2797: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
2798: if (sol) {
2800: prob->exactSol_t[f] = sol;
2801: }
2802: if (ctx) {
2804: prob->exactCtx_t[f] = ctx;
2805: }
2806: PetscFunctionReturn(PETSC_SUCCESS);
2807: }
2809: /*@C
2810: PetscDSGetConstants - Returns the array of constants passed to point functions
2812: Not Collective
2814: Input Parameter:
2815: . prob - The `PetscDS` object
2817: Output Parameters:
2818: + numConstants - The number of constants
2819: - constants - The array of constants, NULL if there are none
2821: Level: intermediate
2823: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSCreate()`
2824: @*/
2825: PetscErrorCode PetscDSGetConstants(PetscDS prob, PetscInt *numConstants, const PetscScalar *constants[])
2826: {
2827: PetscFunctionBegin;
2829: if (numConstants) {
2830: PetscAssertPointer(numConstants, 2);
2831: *numConstants = prob->numConstants;
2832: }
2833: if (constants) {
2834: PetscAssertPointer(constants, 3);
2835: *constants = prob->constants;
2836: }
2837: PetscFunctionReturn(PETSC_SUCCESS);
2838: }
2840: /*@C
2841: PetscDSSetConstants - Set the array of constants passed to point functions
2843: Not Collective
2845: Input Parameters:
2846: + prob - The `PetscDS` object
2847: . numConstants - The number of constants
2848: - constants - The array of constants, NULL if there are none
2850: Level: intermediate
2852: .seealso: `PetscDS`, `PetscDSGetConstants()`, `PetscDSCreate()`
2853: @*/
2854: PetscErrorCode PetscDSSetConstants(PetscDS prob, PetscInt numConstants, PetscScalar constants[])
2855: {
2856: PetscFunctionBegin;
2858: if (numConstants != prob->numConstants) {
2859: PetscCall(PetscFree(prob->constants));
2860: prob->numConstants = numConstants;
2861: if (prob->numConstants) {
2862: PetscCall(PetscMalloc1(prob->numConstants, &prob->constants));
2863: } else {
2864: prob->constants = NULL;
2865: }
2866: }
2867: if (prob->numConstants) {
2868: PetscAssertPointer(constants, 3);
2869: PetscCall(PetscArraycpy(prob->constants, constants, prob->numConstants));
2870: }
2871: PetscFunctionReturn(PETSC_SUCCESS);
2872: }
2874: /*@
2875: PetscDSGetFieldIndex - Returns the index of the given field
2877: Not Collective
2879: Input Parameters:
2880: + prob - The `PetscDS` object
2881: - disc - The discretization object
2883: Output Parameter:
2884: . f - The field number
2886: Level: beginner
2888: .seealso: `PetscDS`, `PetscGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2889: @*/
2890: PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f)
2891: {
2892: PetscInt g;
2894: PetscFunctionBegin;
2896: PetscAssertPointer(f, 3);
2897: *f = -1;
2898: for (g = 0; g < prob->Nf; ++g) {
2899: if (disc == prob->disc[g]) break;
2900: }
2901: PetscCheck(g != prob->Nf, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS.");
2902: *f = g;
2903: PetscFunctionReturn(PETSC_SUCCESS);
2904: }
2906: /*@
2907: PetscDSGetFieldSize - Returns the size of the given field in the full space basis
2909: Not Collective
2911: Input Parameters:
2912: + prob - The `PetscDS` object
2913: - f - The field number
2915: Output Parameter:
2916: . size - The size
2918: Level: beginner
2920: .seealso: `PetscDS`, `PetscDSGetFieldOffset()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2921: @*/
2922: PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size)
2923: {
2924: PetscFunctionBegin;
2926: PetscAssertPointer(size, 3);
2927: 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);
2928: PetscCall(PetscDSSetUp(prob));
2929: *size = prob->Nb[f];
2930: PetscFunctionReturn(PETSC_SUCCESS);
2931: }
2933: /*@
2934: PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis
2936: Not Collective
2938: Input Parameters:
2939: + prob - The `PetscDS` object
2940: - f - The field number
2942: Output Parameter:
2943: . off - The offset
2945: Level: beginner
2947: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2948: @*/
2949: PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off)
2950: {
2951: PetscInt size, g;
2953: PetscFunctionBegin;
2955: PetscAssertPointer(off, 3);
2956: 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);
2957: *off = 0;
2958: for (g = 0; g < f; ++g) {
2959: PetscCall(PetscDSGetFieldSize(prob, g, &size));
2960: *off += size;
2961: }
2962: PetscFunctionReturn(PETSC_SUCCESS);
2963: }
2965: /*@
2966: PetscDSGetFieldOffsetCohesive - Returns the offset of the given field in the full space basis on a cohesive cell
2968: Not Collective
2970: Input Parameters:
2971: + ds - The `PetscDS` object
2972: - f - The field number
2974: Output Parameter:
2975: . off - The offset
2977: Level: beginner
2979: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2980: @*/
2981: PetscErrorCode PetscDSGetFieldOffsetCohesive(PetscDS ds, PetscInt f, PetscInt *off)
2982: {
2983: PetscInt size, g;
2985: PetscFunctionBegin;
2987: PetscAssertPointer(off, 3);
2988: 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);
2989: *off = 0;
2990: for (g = 0; g < f; ++g) {
2991: PetscBool cohesive;
2993: PetscCall(PetscDSGetCohesive(ds, g, &cohesive));
2994: PetscCall(PetscDSGetFieldSize(ds, g, &size));
2995: *off += cohesive ? size : size * 2;
2996: }
2997: PetscFunctionReturn(PETSC_SUCCESS);
2998: }
3000: /*@
3001: PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point
3003: Not Collective
3005: Input Parameter:
3006: . prob - The `PetscDS` object
3008: Output Parameter:
3009: . dimensions - The number of dimensions
3011: Level: beginner
3013: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3014: @*/
3015: PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[])
3016: {
3017: PetscFunctionBegin;
3019: PetscCall(PetscDSSetUp(prob));
3020: PetscAssertPointer(dimensions, 2);
3021: *dimensions = prob->Nb;
3022: PetscFunctionReturn(PETSC_SUCCESS);
3023: }
3025: /*@
3026: PetscDSGetComponents - Returns the number of components for each field on an evaluation point
3028: Not Collective
3030: Input Parameter:
3031: . prob - The `PetscDS` object
3033: Output Parameter:
3034: . components - The number of components
3036: Level: beginner
3038: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3039: @*/
3040: PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[])
3041: {
3042: PetscFunctionBegin;
3044: PetscCall(PetscDSSetUp(prob));
3045: PetscAssertPointer(components, 2);
3046: *components = prob->Nc;
3047: PetscFunctionReturn(PETSC_SUCCESS);
3048: }
3050: /*@
3051: PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point
3053: Not Collective
3055: Input Parameters:
3056: + prob - The `PetscDS` object
3057: - f - The field number
3059: Output Parameter:
3060: . off - The offset
3062: Level: beginner
3064: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3065: @*/
3066: PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off)
3067: {
3068: PetscFunctionBegin;
3070: PetscAssertPointer(off, 3);
3071: 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);
3072: PetscCall(PetscDSSetUp(prob));
3073: *off = prob->off[f];
3074: PetscFunctionReturn(PETSC_SUCCESS);
3075: }
3077: /*@
3078: PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point
3080: Not Collective
3082: Input Parameter:
3083: . prob - The `PetscDS` object
3085: Output Parameter:
3086: . offsets - The offsets
3088: Level: beginner
3090: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3091: @*/
3092: PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[])
3093: {
3094: PetscFunctionBegin;
3096: PetscAssertPointer(offsets, 2);
3097: PetscCall(PetscDSSetUp(prob));
3098: *offsets = prob->off;
3099: PetscFunctionReturn(PETSC_SUCCESS);
3100: }
3102: /*@
3103: PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point
3105: Not Collective
3107: Input Parameter:
3108: . prob - The `PetscDS` object
3110: Output Parameter:
3111: . offsets - The offsets
3113: Level: beginner
3115: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3116: @*/
3117: PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[])
3118: {
3119: PetscFunctionBegin;
3121: PetscAssertPointer(offsets, 2);
3122: PetscCall(PetscDSSetUp(prob));
3123: *offsets = prob->offDer;
3124: PetscFunctionReturn(PETSC_SUCCESS);
3125: }
3127: /*@
3128: PetscDSGetComponentOffsetsCohesive - Returns the offset of each field on an evaluation point
3130: Not Collective
3132: Input Parameters:
3133: + ds - The `PetscDS` object
3134: - s - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive
3136: Output Parameter:
3137: . offsets - The offsets
3139: Level: beginner
3141: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3142: @*/
3143: PetscErrorCode PetscDSGetComponentOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
3144: {
3145: PetscFunctionBegin;
3147: PetscAssertPointer(offsets, 3);
3148: PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS");
3149: PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s);
3150: PetscCall(PetscDSSetUp(ds));
3151: *offsets = ds->offCohesive[s];
3152: PetscFunctionReturn(PETSC_SUCCESS);
3153: }
3155: /*@
3156: PetscDSGetComponentDerivativeOffsetsCohesive - Returns the offset of each field derivative on an evaluation point
3158: Not Collective
3160: Input Parameters:
3161: + ds - The `PetscDS` object
3162: - s - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive
3164: Output Parameter:
3165: . offsets - The offsets
3167: Level: beginner
3169: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
3170: @*/
3171: PetscErrorCode PetscDSGetComponentDerivativeOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
3172: {
3173: PetscFunctionBegin;
3175: PetscAssertPointer(offsets, 3);
3176: PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS");
3177: PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s);
3178: PetscCall(PetscDSSetUp(ds));
3179: *offsets = ds->offDerCohesive[s];
3180: PetscFunctionReturn(PETSC_SUCCESS);
3181: }
3183: /*@C
3184: PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization
3186: Not Collective
3188: Input Parameter:
3189: . prob - The `PetscDS` object
3191: Output Parameter:
3192: . T - The basis function and derivatives tabulation at quadrature points for each field
3194: Level: intermediate
3196: .seealso: `PetscDS`, `PetscTabulation`, `PetscDSCreate()`
3197: @*/
3198: PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscTabulation *T[])
3199: {
3200: PetscFunctionBegin;
3202: PetscAssertPointer(T, 2);
3203: PetscCall(PetscDSSetUp(prob));
3204: *T = prob->T;
3205: PetscFunctionReturn(PETSC_SUCCESS);
3206: }
3208: /*@C
3209: PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces
3211: Not Collective
3213: Input Parameter:
3214: . prob - The `PetscDS` object
3216: Output Parameter:
3217: . Tf - The basis function and derivative tabulation on each local face at quadrature points for each and field
3219: Level: intermediate
3221: .seealso: `PetscTabulation`, `PetscDS`, `PetscDSGetTabulation()`, `PetscDSCreate()`
3222: @*/
3223: PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscTabulation *Tf[])
3224: {
3225: PetscFunctionBegin;
3227: PetscAssertPointer(Tf, 2);
3228: PetscCall(PetscDSSetUp(prob));
3229: *Tf = prob->Tf;
3230: PetscFunctionReturn(PETSC_SUCCESS);
3231: }
3233: PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x)
3234: {
3235: PetscFunctionBegin;
3237: PetscCall(PetscDSSetUp(prob));
3238: if (u) {
3239: PetscAssertPointer(u, 2);
3240: *u = prob->u;
3241: }
3242: if (u_t) {
3243: PetscAssertPointer(u_t, 3);
3244: *u_t = prob->u_t;
3245: }
3246: if (u_x) {
3247: PetscAssertPointer(u_x, 4);
3248: *u_x = prob->u_x;
3249: }
3250: PetscFunctionReturn(PETSC_SUCCESS);
3251: }
3253: PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3)
3254: {
3255: PetscFunctionBegin;
3257: PetscCall(PetscDSSetUp(prob));
3258: if (f0) {
3259: PetscAssertPointer(f0, 2);
3260: *f0 = prob->f0;
3261: }
3262: if (f1) {
3263: PetscAssertPointer(f1, 3);
3264: *f1 = prob->f1;
3265: }
3266: if (g0) {
3267: PetscAssertPointer(g0, 4);
3268: *g0 = prob->g0;
3269: }
3270: if (g1) {
3271: PetscAssertPointer(g1, 5);
3272: *g1 = prob->g1;
3273: }
3274: if (g2) {
3275: PetscAssertPointer(g2, 6);
3276: *g2 = prob->g2;
3277: }
3278: if (g3) {
3279: PetscAssertPointer(g3, 7);
3280: *g3 = prob->g3;
3281: }
3282: PetscFunctionReturn(PETSC_SUCCESS);
3283: }
3285: PetscErrorCode PetscDSGetWorkspace(PetscDS prob, PetscReal **x, PetscScalar **basisReal, PetscScalar **basisDerReal, PetscScalar **testReal, PetscScalar **testDerReal)
3286: {
3287: PetscFunctionBegin;
3289: PetscCall(PetscDSSetUp(prob));
3290: if (x) {
3291: PetscAssertPointer(x, 2);
3292: *x = prob->x;
3293: }
3294: if (basisReal) {
3295: PetscAssertPointer(basisReal, 3);
3296: *basisReal = prob->basisReal;
3297: }
3298: if (basisDerReal) {
3299: PetscAssertPointer(basisDerReal, 4);
3300: *basisDerReal = prob->basisDerReal;
3301: }
3302: if (testReal) {
3303: PetscAssertPointer(testReal, 5);
3304: *testReal = prob->testReal;
3305: }
3306: if (testDerReal) {
3307: PetscAssertPointer(testDerReal, 6);
3308: *testDerReal = prob->testDerReal;
3309: }
3310: PetscFunctionReturn(PETSC_SUCCESS);
3311: }
3313: /*@C
3314: PetscDSAddBoundary - Add a boundary condition to the model.
3316: Collective
3318: Input Parameters:
3319: + ds - The PetscDS object
3320: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3321: . name - The BC name
3322: . label - The label defining constrained points
3323: . Nv - The number of `DMLabel` values for constrained points
3324: . values - An array of label values for constrained points
3325: . field - The field to constrain
3326: . Nc - The number of constrained field components (0 will constrain all fields)
3327: . comps - An array of constrained component numbers
3328: . bcFunc - A pointwise function giving boundary values
3329: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3330: - ctx - An optional user context for bcFunc
3332: Output Parameter:
3333: . bd - The boundary number
3335: Options Database Keys:
3336: + -bc_<boundary name> <num> - Overrides the boundary ids
3337: - -bc_<boundary name>_comp <num> - Overrides the boundary components
3339: Level: developer
3341: Note:
3342: Both `bcFunc` and `bcFunc_t` will depend on the boundary condition type. If the type if `DM_BC_ESSENTIAL`, then the calling sequence is\:
3344: $ void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])
3346: If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value, then the calling sequence is\:
3347: .vb
3348: void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3349: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3350: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3351: PetscReal time, const PetscReal x[], PetscScalar bcval[])
3352: .ve
3353: + dim - the spatial dimension
3354: . Nf - the number of fields
3355: . uOff - the offset into u[] and u_t[] for each field
3356: . uOff_x - the offset into u_x[] for each field
3357: . u - each field evaluated at the current point
3358: . u_t - the time derivative of each field evaluated at the current point
3359: . u_x - the gradient of each field evaluated at the current point
3360: . aOff - the offset into a[] and a_t[] for each auxiliary field
3361: . aOff_x - the offset into a_x[] for each auxiliary field
3362: . a - each auxiliary field evaluated at the current point
3363: . a_t - the time derivative of each auxiliary field evaluated at the current point
3364: . a_x - the gradient of auxiliary each field evaluated at the current point
3365: . t - current time
3366: . x - coordinates of the current point
3367: . numConstants - number of constant parameters
3368: . constants - constant parameters
3369: - bcval - output values at the current point
3371: Notes:
3372: The pointwise functions are used to provide boundary values for essential boundary
3373: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3374: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3375: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3377: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundaryByName()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3378: @*/
3379: 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)
3380: {
3381: DSBoundary head = ds->boundary, b;
3382: PetscInt n = 0;
3383: const char *lname;
3385: PetscFunctionBegin;
3388: PetscAssertPointer(name, 3);
3393: 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);
3394: if (Nc > 0) {
3395: PetscInt *fcomps;
3396: PetscInt c;
3398: PetscCall(PetscDSGetComponents(ds, &fcomps));
3399: 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);
3400: for (c = 0; c < Nc; ++c) {
3401: 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);
3402: }
3403: }
3404: PetscCall(PetscNew(&b));
3405: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3406: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf));
3407: PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf));
3408: PetscCall(PetscMalloc1(Nv, &b->values));
3409: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3410: PetscCall(PetscMalloc1(Nc, &b->comps));
3411: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3412: PetscCall(PetscObjectGetName((PetscObject)label, &lname));
3413: PetscCall(PetscStrallocpy(lname, (char **)&b->lname));
3414: b->type = type;
3415: b->label = label;
3416: b->Nv = Nv;
3417: b->field = field;
3418: b->Nc = Nc;
3419: b->func = bcFunc;
3420: b->func_t = bcFunc_t;
3421: b->ctx = ctx;
3422: b->next = NULL;
3423: /* Append to linked list so that we can preserve the order */
3424: if (!head) ds->boundary = b;
3425: while (head) {
3426: if (!head->next) {
3427: head->next = b;
3428: head = b;
3429: }
3430: head = head->next;
3431: ++n;
3432: }
3433: if (bd) {
3434: PetscAssertPointer(bd, 13);
3435: *bd = n;
3436: }
3437: PetscFunctionReturn(PETSC_SUCCESS);
3438: }
3440: // PetscClangLinter pragma ignore: -fdoc-section-header-unknown
3441: /*@C
3442: PetscDSAddBoundaryByName - Add a boundary condition to the model.
3444: Collective
3446: Input Parameters:
3447: + ds - The `PetscDS` object
3448: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3449: . name - The BC name
3450: . lname - The naem of the label defining constrained points
3451: . Nv - The number of `DMLabel` values for constrained points
3452: . values - An array of label values for constrained points
3453: . field - The field to constrain
3454: . Nc - The number of constrained field components (0 will constrain all fields)
3455: . comps - An array of constrained component numbers
3456: . bcFunc - A pointwise function giving boundary values
3457: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3458: - ctx - An optional user context for bcFunc
3460: Output Parameter:
3461: . bd - The boundary number
3463: Options Database Keys:
3464: + -bc_<boundary name> <num> - Overrides the boundary ids
3465: - -bc_<boundary name>_comp <num> - Overrides the boundary components
3467: Calling Sequence of `bcFunc` and `bcFunc_t`:
3468: If the type is `DM_BC_ESSENTIAL`
3469: .vb
3470: void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])
3471: .ve
3472: If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value,
3473: .vb
3474: void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3475: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3476: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3477: PetscReal time, const PetscReal x[], PetscScalar bcval[])
3478: .ve
3479: + dim - the spatial dimension
3480: . Nf - the number of fields
3481: . uOff - the offset into u[] and u_t[] for each field
3482: . uOff_x - the offset into u_x[] for each field
3483: . u - each field evaluated at the current point
3484: . u_t - the time derivative of each field evaluated at the current point
3485: . u_x - the gradient of each field evaluated at the current point
3486: . aOff - the offset into a[] and a_t[] for each auxiliary field
3487: . aOff_x - the offset into a_x[] for each auxiliary field
3488: . a - each auxiliary field evaluated at the current point
3489: . a_t - the time derivative of each auxiliary field evaluated at the current point
3490: . a_x - the gradient of auxiliary each field evaluated at the current point
3491: . t - current time
3492: . x - coordinates of the current point
3493: . numConstants - number of constant parameters
3494: . constants - constant parameters
3495: - bcval - output values at the current point
3497: Level: developer
3499: Notes:
3500: The pointwise functions are used to provide boundary values for essential boundary
3501: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3502: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3503: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3505: This function should only be used with `DMFOREST` currently, since labels cannot be defined before the underlying `DMPLEX` is built.
3507: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3508: @*/
3509: 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)
3510: {
3511: DSBoundary head = ds->boundary, b;
3512: PetscInt n = 0;
3514: PetscFunctionBegin;
3517: PetscAssertPointer(name, 3);
3518: PetscAssertPointer(lname, 4);
3522: PetscCall(PetscNew(&b));
3523: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3524: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf));
3525: PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf));
3526: PetscCall(PetscMalloc1(Nv, &b->values));
3527: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3528: PetscCall(PetscMalloc1(Nc, &b->comps));
3529: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3530: PetscCall(PetscStrallocpy(lname, (char **)&b->lname));
3531: b->type = type;
3532: b->label = NULL;
3533: b->Nv = Nv;
3534: b->field = field;
3535: b->Nc = Nc;
3536: b->func = bcFunc;
3537: b->func_t = bcFunc_t;
3538: b->ctx = ctx;
3539: b->next = NULL;
3540: /* Append to linked list so that we can preserve the order */
3541: if (!head) ds->boundary = b;
3542: while (head) {
3543: if (!head->next) {
3544: head->next = b;
3545: head = b;
3546: }
3547: head = head->next;
3548: ++n;
3549: }
3550: if (bd) {
3551: PetscAssertPointer(bd, 13);
3552: *bd = n;
3553: }
3554: PetscFunctionReturn(PETSC_SUCCESS);
3555: }
3557: /*@C
3558: PetscDSUpdateBoundary - Change a boundary condition for the model.
3560: Input Parameters:
3561: + ds - The `PetscDS` object
3562: . bd - The boundary condition number
3563: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3564: . name - The BC name
3565: . label - The label defining constrained points
3566: . Nv - The number of `DMLabel` ids for constrained points
3567: . values - An array of ids for constrained points
3568: . field - The field to constrain
3569: . Nc - The number of constrained field components
3570: . comps - An array of constrained component numbers
3571: . bcFunc - A pointwise function giving boundary values
3572: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or NULL
3573: - ctx - An optional user context for bcFunc
3575: Level: developer
3577: Notes:
3578: The pointwise functions are used to provide boundary values for essential boundary
3579: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3580: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3581: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3583: The boundary condition number is the order in which it was registered. The user can get the number of boundary conditions from `PetscDSGetNumBoundary()`.
3584: See `PetscDSAddBoundary()` for a description of the calling sequences for the callbacks.
3586: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSGetNumBoundary()`, `DMLabel`
3587: @*/
3588: 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)
3589: {
3590: DSBoundary b = ds->boundary;
3591: PetscInt n = 0;
3593: PetscFunctionBegin;
3595: while (b) {
3596: if (n == bd) break;
3597: b = b->next;
3598: ++n;
3599: }
3600: PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n);
3601: if (name) {
3602: PetscCall(PetscFree(b->name));
3603: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3604: }
3605: b->type = type;
3606: if (label) {
3607: const char *name;
3609: b->label = label;
3610: PetscCall(PetscFree(b->lname));
3611: PetscCall(PetscObjectGetName((PetscObject)label, &name));
3612: PetscCall(PetscStrallocpy(name, (char **)&b->lname));
3613: }
3614: if (Nv >= 0) {
3615: b->Nv = Nv;
3616: PetscCall(PetscFree(b->values));
3617: PetscCall(PetscMalloc1(Nv, &b->values));
3618: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3619: }
3620: if (field >= 0) b->field = field;
3621: if (Nc >= 0) {
3622: b->Nc = Nc;
3623: PetscCall(PetscFree(b->comps));
3624: PetscCall(PetscMalloc1(Nc, &b->comps));
3625: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3626: }
3627: if (bcFunc) b->func = bcFunc;
3628: if (bcFunc_t) b->func_t = bcFunc_t;
3629: if (ctx) b->ctx = ctx;
3630: PetscFunctionReturn(PETSC_SUCCESS);
3631: }
3633: /*@
3634: PetscDSGetNumBoundary - Get the number of registered BC
3636: Input Parameter:
3637: . ds - The `PetscDS` object
3639: Output Parameter:
3640: . numBd - The number of BC
3642: Level: intermediate
3644: .seealso: `PetscDS`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`
3645: @*/
3646: PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd)
3647: {
3648: DSBoundary b = ds->boundary;
3650: PetscFunctionBegin;
3652: PetscAssertPointer(numBd, 2);
3653: *numBd = 0;
3654: while (b) {
3655: ++(*numBd);
3656: b = b->next;
3657: }
3658: PetscFunctionReturn(PETSC_SUCCESS);
3659: }
3661: /*@C
3662: PetscDSGetBoundary - Gets a boundary condition to the model
3664: Input Parameters:
3665: + ds - The `PetscDS` object
3666: - bd - The BC number
3668: Output Parameters:
3669: + wf - The `PetscWeakForm` holding the pointwise functions
3670: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3671: . name - The BC name
3672: . label - The label defining constrained points
3673: . Nv - The number of `DMLabel` ids for constrained points
3674: . values - An array of ids for constrained points
3675: . field - The field to constrain
3676: . Nc - The number of constrained field components
3677: . comps - An array of constrained component numbers
3678: . func - A pointwise function giving boundary values
3679: . func_t - A pointwise function giving the time derivative of the boundary values
3680: - ctx - An optional user context for bcFunc
3682: Options Database Keys:
3683: + -bc_<boundary name> <num> - Overrides the boundary ids
3684: - -bc_<boundary name>_comp <num> - Overrides the boundary components
3686: Level: developer
3688: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `DMLabel`
3689: @*/
3690: 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)
3691: {
3692: DSBoundary b = ds->boundary;
3693: PetscInt n = 0;
3695: PetscFunctionBegin;
3697: while (b) {
3698: if (n == bd) break;
3699: b = b->next;
3700: ++n;
3701: }
3702: PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n);
3703: if (wf) {
3704: PetscAssertPointer(wf, 3);
3705: *wf = b->wf;
3706: }
3707: if (type) {
3708: PetscAssertPointer(type, 4);
3709: *type = b->type;
3710: }
3711: if (name) {
3712: PetscAssertPointer(name, 5);
3713: *name = b->name;
3714: }
3715: if (label) {
3716: PetscAssertPointer(label, 6);
3717: *label = b->label;
3718: }
3719: if (Nv) {
3720: PetscAssertPointer(Nv, 7);
3721: *Nv = b->Nv;
3722: }
3723: if (values) {
3724: PetscAssertPointer(values, 8);
3725: *values = b->values;
3726: }
3727: if (field) {
3728: PetscAssertPointer(field, 9);
3729: *field = b->field;
3730: }
3731: if (Nc) {
3732: PetscAssertPointer(Nc, 10);
3733: *Nc = b->Nc;
3734: }
3735: if (comps) {
3736: PetscAssertPointer(comps, 11);
3737: *comps = b->comps;
3738: }
3739: if (func) {
3740: PetscAssertPointer(func, 12);
3741: *func = b->func;
3742: }
3743: if (func_t) {
3744: PetscAssertPointer(func_t, 13);
3745: *func_t = b->func_t;
3746: }
3747: if (ctx) {
3748: PetscAssertPointer(ctx, 14);
3749: *ctx = b->ctx;
3750: }
3751: PetscFunctionReturn(PETSC_SUCCESS);
3752: }
3754: static PetscErrorCode DSBoundaryDuplicate_Internal(DSBoundary b, DSBoundary *bNew)
3755: {
3756: PetscFunctionBegin;
3757: PetscCall(PetscNew(bNew));
3758: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &(*bNew)->wf));
3759: PetscCall(PetscWeakFormCopy(b->wf, (*bNew)->wf));
3760: PetscCall(PetscStrallocpy(b->name, (char **)&((*bNew)->name)));
3761: PetscCall(PetscStrallocpy(b->lname, (char **)&((*bNew)->lname)));
3762: (*bNew)->type = b->type;
3763: (*bNew)->label = b->label;
3764: (*bNew)->Nv = b->Nv;
3765: PetscCall(PetscMalloc1(b->Nv, &(*bNew)->values));
3766: PetscCall(PetscArraycpy((*bNew)->values, b->values, b->Nv));
3767: (*bNew)->field = b->field;
3768: (*bNew)->Nc = b->Nc;
3769: PetscCall(PetscMalloc1(b->Nc, &(*bNew)->comps));
3770: PetscCall(PetscArraycpy((*bNew)->comps, b->comps, b->Nc));
3771: (*bNew)->func = b->func;
3772: (*bNew)->func_t = b->func_t;
3773: (*bNew)->ctx = b->ctx;
3774: PetscFunctionReturn(PETSC_SUCCESS);
3775: }
3777: /*@
3778: PetscDSCopyBoundary - Copy all boundary condition objects to the new problem
3780: Not Collective
3782: Input Parameters:
3783: + ds - The source `PetscDS` object
3784: . numFields - The number of selected fields, or `PETSC_DEFAULT` for all fields
3785: - fields - The selected fields, or NULL for all fields
3787: Output Parameter:
3788: . newds - The target `PetscDS`, now with a copy of the boundary conditions
3790: Level: intermediate
3792: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3793: @*/
3794: PetscErrorCode PetscDSCopyBoundary(PetscDS ds, PetscInt numFields, const PetscInt fields[], PetscDS newds)
3795: {
3796: DSBoundary b, *lastnext;
3798: PetscFunctionBegin;
3801: if (ds == newds) PetscFunctionReturn(PETSC_SUCCESS);
3802: PetscCall(PetscDSDestroyBoundary(newds));
3803: lastnext = &(newds->boundary);
3804: for (b = ds->boundary; b; b = b->next) {
3805: DSBoundary bNew;
3806: PetscInt fieldNew = -1;
3808: if (numFields > 0 && fields) {
3809: PetscInt f;
3811: for (f = 0; f < numFields; ++f)
3812: if (b->field == fields[f]) break;
3813: if (f == numFields) continue;
3814: fieldNew = f;
3815: }
3816: PetscCall(DSBoundaryDuplicate_Internal(b, &bNew));
3817: bNew->field = fieldNew < 0 ? b->field : fieldNew;
3818: *lastnext = bNew;
3819: lastnext = &(bNew->next);
3820: }
3821: PetscFunctionReturn(PETSC_SUCCESS);
3822: }
3824: /*@
3825: PetscDSDestroyBoundary - Remove all `DMBoundary` objects from the `PetscDS`
3827: Not Collective
3829: Input Parameter:
3830: . ds - The `PetscDS` object
3832: Level: intermediate
3834: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`
3835: @*/
3836: PetscErrorCode PetscDSDestroyBoundary(PetscDS ds)
3837: {
3838: DSBoundary next = ds->boundary;
3840: PetscFunctionBegin;
3841: while (next) {
3842: DSBoundary b = next;
3844: next = b->next;
3845: PetscCall(PetscWeakFormDestroy(&b->wf));
3846: PetscCall(PetscFree(b->name));
3847: PetscCall(PetscFree(b->lname));
3848: PetscCall(PetscFree(b->values));
3849: PetscCall(PetscFree(b->comps));
3850: PetscCall(PetscFree(b));
3851: }
3852: PetscFunctionReturn(PETSC_SUCCESS);
3853: }
3855: /*@
3856: PetscDSSelectDiscretizations - Copy discretizations to the new problem with different field layout
3858: Not Collective
3860: Input Parameters:
3861: + prob - The `PetscDS` object
3862: . numFields - Number of new fields
3863: - fields - Old field number for each new field
3865: Output Parameter:
3866: . newprob - The `PetscDS` copy
3868: Level: intermediate
3870: .seealso: `PetscDS`, `PetscDSSelectEquations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3871: @*/
3872: PetscErrorCode PetscDSSelectDiscretizations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
3873: {
3874: PetscInt Nf, Nfn, fn;
3876: PetscFunctionBegin;
3878: if (fields) PetscAssertPointer(fields, 3);
3880: PetscCall(PetscDSGetNumFields(prob, &Nf));
3881: PetscCall(PetscDSGetNumFields(newprob, &Nfn));
3882: numFields = numFields < 0 ? Nf : numFields;
3883: for (fn = 0; fn < numFields; ++fn) {
3884: const PetscInt f = fields ? fields[fn] : fn;
3885: PetscObject disc;
3887: if (f >= Nf) continue;
3888: PetscCall(PetscDSGetDiscretization(prob, f, &disc));
3889: PetscCall(PetscDSSetDiscretization(newprob, fn, disc));
3890: }
3891: PetscFunctionReturn(PETSC_SUCCESS);
3892: }
3894: /*@
3895: PetscDSSelectEquations - Copy pointwise function pointers to the new problem with different field layout
3897: Not Collective
3899: Input Parameters:
3900: + prob - The `PetscDS` object
3901: . numFields - Number of new fields
3902: - fields - Old field number for each new field
3904: Output Parameter:
3905: . newprob - The `PetscDS` copy
3907: Level: intermediate
3909: .seealso: `PetscDS`, `PetscDSSelectDiscretizations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3910: @*/
3911: PetscErrorCode PetscDSSelectEquations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
3912: {
3913: PetscInt Nf, Nfn, fn, gn;
3915: PetscFunctionBegin;
3917: if (fields) PetscAssertPointer(fields, 3);
3919: PetscCall(PetscDSGetNumFields(prob, &Nf));
3920: PetscCall(PetscDSGetNumFields(newprob, &Nfn));
3921: 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);
3922: for (fn = 0; fn < numFields; ++fn) {
3923: const PetscInt f = fields ? fields[fn] : fn;
3924: PetscPointFunc obj;
3925: PetscPointFunc f0, f1;
3926: PetscBdPointFunc f0Bd, f1Bd;
3927: PetscRiemannFunc r;
3929: if (f >= Nf) continue;
3930: PetscCall(PetscDSGetObjective(prob, f, &obj));
3931: PetscCall(PetscDSGetResidual(prob, f, &f0, &f1));
3932: PetscCall(PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd));
3933: PetscCall(PetscDSGetRiemannSolver(prob, f, &r));
3934: PetscCall(PetscDSSetObjective(newprob, fn, obj));
3935: PetscCall(PetscDSSetResidual(newprob, fn, f0, f1));
3936: PetscCall(PetscDSSetBdResidual(newprob, fn, f0Bd, f1Bd));
3937: PetscCall(PetscDSSetRiemannSolver(newprob, fn, r));
3938: for (gn = 0; gn < numFields; ++gn) {
3939: const PetscInt g = fields ? fields[gn] : gn;
3940: PetscPointJac g0, g1, g2, g3;
3941: PetscPointJac g0p, g1p, g2p, g3p;
3942: PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd;
3944: if (g >= Nf) continue;
3945: PetscCall(PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3));
3946: PetscCall(PetscDSGetJacobianPreconditioner(prob, f, g, &g0p, &g1p, &g2p, &g3p));
3947: PetscCall(PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd));
3948: PetscCall(PetscDSSetJacobian(newprob, fn, gn, g0, g1, g2, g3));
3949: PetscCall(PetscDSSetJacobianPreconditioner(newprob, fn, gn, g0p, g1p, g2p, g3p));
3950: PetscCall(PetscDSSetBdJacobian(newprob, fn, gn, g0Bd, g1Bd, g2Bd, g3Bd));
3951: }
3952: }
3953: PetscFunctionReturn(PETSC_SUCCESS);
3954: }
3956: /*@
3957: PetscDSCopyEquations - Copy all pointwise function pointers to another `PetscDS`
3959: Not Collective
3961: Input Parameter:
3962: . prob - The `PetscDS` object
3964: Output Parameter:
3965: . newprob - The `PetscDS` copy
3967: Level: intermediate
3969: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3970: @*/
3971: PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob)
3972: {
3973: PetscWeakForm wf, newwf;
3974: PetscInt Nf, Ng;
3976: PetscFunctionBegin;
3979: PetscCall(PetscDSGetNumFields(prob, &Nf));
3980: PetscCall(PetscDSGetNumFields(newprob, &Ng));
3981: PetscCheck(Nf == Ng, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %" PetscInt_FMT " != %" PetscInt_FMT, Nf, Ng);
3982: PetscCall(PetscDSGetWeakForm(prob, &wf));
3983: PetscCall(PetscDSGetWeakForm(newprob, &newwf));
3984: PetscCall(PetscWeakFormCopy(wf, newwf));
3985: PetscFunctionReturn(PETSC_SUCCESS);
3986: }
3988: /*@
3989: PetscDSCopyConstants - Copy all constants to another `PetscDS`
3991: Not Collective
3993: Input Parameter:
3994: . prob - The `PetscDS` object
3996: Output Parameter:
3997: . newprob - The `PetscDS` copy
3999: Level: intermediate
4001: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
4002: @*/
4003: PetscErrorCode PetscDSCopyConstants(PetscDS prob, PetscDS newprob)
4004: {
4005: PetscInt Nc;
4006: const PetscScalar *constants;
4008: PetscFunctionBegin;
4011: PetscCall(PetscDSGetConstants(prob, &Nc, &constants));
4012: PetscCall(PetscDSSetConstants(newprob, Nc, (PetscScalar *)constants));
4013: PetscFunctionReturn(PETSC_SUCCESS);
4014: }
4016: /*@
4017: PetscDSCopyExactSolutions - Copy all exact solutions to another `PetscDS`
4019: Not Collective
4021: Input Parameter:
4022: . ds - The `PetscDS` object
4024: Output Parameter:
4025: . newds - The `PetscDS` copy
4027: Level: intermediate
4029: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
4030: @*/
4031: PetscErrorCode PetscDSCopyExactSolutions(PetscDS ds, PetscDS newds)
4032: {
4033: PetscSimplePointFunc sol;
4034: void *ctx;
4035: PetscInt Nf, f;
4037: PetscFunctionBegin;
4040: PetscCall(PetscDSGetNumFields(ds, &Nf));
4041: for (f = 0; f < Nf; ++f) {
4042: PetscCall(PetscDSGetExactSolution(ds, f, &sol, &ctx));
4043: PetscCall(PetscDSSetExactSolution(newds, f, sol, ctx));
4044: PetscCall(PetscDSGetExactSolutionTimeDerivative(ds, f, &sol, &ctx));
4045: PetscCall(PetscDSSetExactSolutionTimeDerivative(newds, f, sol, ctx));
4046: }
4047: PetscFunctionReturn(PETSC_SUCCESS);
4048: }
4050: PetscErrorCode PetscDSCopy(PetscDS ds, DM dmNew, PetscDS dsNew)
4051: {
4052: DSBoundary b;
4053: PetscInt cdim, Nf, f, d;
4054: PetscBool isCohesive;
4055: void *ctx;
4057: PetscFunctionBegin;
4058: PetscCall(PetscDSCopyConstants(ds, dsNew));
4059: PetscCall(PetscDSCopyExactSolutions(ds, dsNew));
4060: PetscCall(PetscDSSelectDiscretizations(ds, PETSC_DETERMINE, NULL, dsNew));
4061: PetscCall(PetscDSCopyEquations(ds, dsNew));
4062: PetscCall(PetscDSGetNumFields(ds, &Nf));
4063: for (f = 0; f < Nf; ++f) {
4064: PetscCall(PetscDSGetContext(ds, f, &ctx));
4065: PetscCall(PetscDSSetContext(dsNew, f, ctx));
4066: PetscCall(PetscDSGetCohesive(ds, f, &isCohesive));
4067: PetscCall(PetscDSSetCohesive(dsNew, f, isCohesive));
4068: PetscCall(PetscDSGetJetDegree(ds, f, &d));
4069: PetscCall(PetscDSSetJetDegree(dsNew, f, d));
4070: }
4071: if (Nf) {
4072: PetscCall(PetscDSGetCoordinateDimension(ds, &cdim));
4073: PetscCall(PetscDSSetCoordinateDimension(dsNew, cdim));
4074: }
4075: PetscCall(PetscDSCopyBoundary(ds, PETSC_DETERMINE, NULL, dsNew));
4076: for (b = dsNew->boundary; b; b = b->next) {
4077: PetscCall(DMGetLabel(dmNew, b->lname, &b->label));
4078: /* Do not check if label exists here, since p4est calls this for the reference tree which does not have the labels */
4079: //PetscCheck(b->label,PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Label %s missing in new DM", name);
4080: }
4081: PetscFunctionReturn(PETSC_SUCCESS);
4082: }
4084: PetscErrorCode PetscDSGetHeightSubspace(PetscDS prob, PetscInt height, PetscDS *subprob)
4085: {
4086: PetscInt dim, Nf, f;
4088: PetscFunctionBegin;
4090: PetscAssertPointer(subprob, 3);
4091: if (height == 0) {
4092: *subprob = prob;
4093: PetscFunctionReturn(PETSC_SUCCESS);
4094: }
4095: PetscCall(PetscDSGetNumFields(prob, &Nf));
4096: PetscCall(PetscDSGetSpatialDimension(prob, &dim));
4097: PetscCheck(height <= dim, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_OUTOFRANGE, "DS can only handle height in [0, %" PetscInt_FMT "], not %" PetscInt_FMT, dim, height);
4098: if (!prob->subprobs) PetscCall(PetscCalloc1(dim, &prob->subprobs));
4099: if (!prob->subprobs[height - 1]) {
4100: PetscInt cdim;
4102: PetscCall(PetscDSCreate(PetscObjectComm((PetscObject)prob), &prob->subprobs[height - 1]));
4103: PetscCall(PetscDSGetCoordinateDimension(prob, &cdim));
4104: PetscCall(PetscDSSetCoordinateDimension(prob->subprobs[height - 1], cdim));
4105: for (f = 0; f < Nf; ++f) {
4106: PetscFE subfe;
4107: PetscObject obj;
4108: PetscClassId id;
4110: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
4111: PetscCall(PetscObjectGetClassId(obj, &id));
4112: if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetHeightSubspace((PetscFE)obj, height, &subfe));
4113: else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unsupported discretization type for field %" PetscInt_FMT, f);
4114: PetscCall(PetscDSSetDiscretization(prob->subprobs[height - 1], f, (PetscObject)subfe));
4115: }
4116: }
4117: *subprob = prob->subprobs[height - 1];
4118: PetscFunctionReturn(PETSC_SUCCESS);
4119: }
4121: PetscErrorCode PetscDSPermuteQuadPoint(PetscDS ds, PetscInt ornt, PetscInt field, PetscInt q, PetscInt *qperm)
4122: {
4123: IS permIS;
4124: PetscQuadrature quad;
4125: DMPolytopeType ct;
4126: const PetscInt *perm;
4127: PetscInt Na, Nq;
4129: PetscFunctionBeginHot;
4130: PetscCall(PetscFEGetQuadrature((PetscFE)ds->disc[field], &quad));
4131: PetscCall(PetscQuadratureGetData(quad, NULL, NULL, &Nq, NULL, NULL));
4132: PetscCall(PetscQuadratureGetCellType(quad, &ct));
4133: PetscCheck(q >= 0 && q < Nq, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Quadrature point %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", q, Nq);
4134: Na = DMPolytopeTypeGetNumArrangments(ct) / 2;
4135: 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);
4136: if (!ds->quadPerm[(PetscInt)ct]) PetscCall(PetscQuadratureComputePermutations(quad, NULL, &ds->quadPerm[(PetscInt)ct]));
4137: permIS = ds->quadPerm[(PetscInt)ct][ornt + Na];
4138: PetscCall(ISGetIndices(permIS, &perm));
4139: *qperm = perm[q];
4140: PetscCall(ISRestoreIndices(permIS, &perm));
4141: PetscFunctionReturn(PETSC_SUCCESS);
4142: }
4144: PetscErrorCode PetscDSGetDiscType_Internal(PetscDS ds, PetscInt f, PetscDiscType *disctype)
4145: {
4146: PetscObject obj;
4147: PetscClassId id;
4148: PetscInt Nf;
4150: PetscFunctionBegin;
4152: PetscAssertPointer(disctype, 3);
4153: *disctype = PETSC_DISC_NONE;
4154: PetscCall(PetscDSGetNumFields(ds, &Nf));
4155: PetscCheck(f < Nf, PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_SIZ, "Field %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, Nf);
4156: PetscCall(PetscDSGetDiscretization(ds, f, &obj));
4157: if (obj) {
4158: PetscCall(PetscObjectGetClassId(obj, &id));
4159: if (id == PETSCFE_CLASSID) *disctype = PETSC_DISC_FE;
4160: else *disctype = PETSC_DISC_FV;
4161: }
4162: PetscFunctionReturn(PETSC_SUCCESS);
4163: }
4165: static PetscErrorCode PetscDSDestroy_Basic(PetscDS ds)
4166: {
4167: PetscFunctionBegin;
4168: PetscCall(PetscFree(ds->data));
4169: PetscFunctionReturn(PETSC_SUCCESS);
4170: }
4172: static PetscErrorCode PetscDSInitialize_Basic(PetscDS ds)
4173: {
4174: PetscFunctionBegin;
4175: ds->ops->setfromoptions = NULL;
4176: ds->ops->setup = NULL;
4177: ds->ops->view = NULL;
4178: ds->ops->destroy = PetscDSDestroy_Basic;
4179: PetscFunctionReturn(PETSC_SUCCESS);
4180: }
4182: /*MC
4183: PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions
4185: Level: intermediate
4187: .seealso: `PetscDSType`, `PetscDSCreate()`, `PetscDSSetType()`
4188: M*/
4190: PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS ds)
4191: {
4192: PetscDS_Basic *b;
4194: PetscFunctionBegin;
4196: PetscCall(PetscNew(&b));
4197: ds->data = b;
4199: PetscCall(PetscDSInitialize_Basic(ds));
4200: PetscFunctionReturn(PETSC_SUCCESS);
4201: }