Actual source code: dtds.c
1: #include <petsc/private/petscdsimpl.h>
3: PetscClassId PETSCDS_CLASSID = 0;
5: PetscFunctionList PetscDSList = NULL;
6: PetscBool PetscDSRegisterAllCalled = PETSC_FALSE;
8: /*@
9: PetscDSRegister - Adds a new `PetscDS` implementation
11: Not Collective; No Fortran Support
13: Input Parameters:
14: + sname - The name of a new user-defined creation routine
15: - function - The creation routine itself
17: Example Usage:
18: .vb
19: PetscDSRegister("my_ds", MyPetscDSCreate);
20: .ve
22: Then, your PetscDS type can be chosen with the procedural interface via
23: .vb
24: PetscDSCreate(MPI_Comm, PetscDS *);
25: PetscDSSetType(PetscDS, "my_ds");
26: .ve
27: or at runtime via the option
28: .vb
29: -petscds_type my_ds
30: .ve
32: Level: advanced
34: Note:
35: `PetscDSRegister()` may be called multiple times to add several user-defined `PetscDSs`
37: .seealso: `PetscDSType`, `PetscDS`, `PetscDSRegisterAll()`
38: @*/
39: PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS))
40: {
41: PetscFunctionBegin;
42: PetscCall(PetscFunctionListAdd(&PetscDSList, sname, function));
43: PetscFunctionReturn(PETSC_SUCCESS);
44: }
46: /*@
47: PetscDSSetType - Builds a particular `PetscDS`
49: Collective; No Fortran Support
51: Input Parameters:
52: + prob - The `PetscDS` object
53: - name - The `PetscDSType`
55: Options Database Key:
56: . -petscds_type type - Sets the PetscDS type; use -help for a list of available types
58: Level: intermediate
60: .seealso: `PetscDSType`, `PetscDS`, `PetscDSGetType()`, `PetscDSCreate()`
61: @*/
62: PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name)
63: {
64: PetscErrorCode (*r)(PetscDS);
65: PetscBool match;
67: PetscFunctionBegin;
69: PetscCall(PetscObjectTypeCompare((PetscObject)prob, name, &match));
70: if (match) PetscFunctionReturn(PETSC_SUCCESS);
72: PetscCall(PetscDSRegisterAll());
73: PetscCall(PetscFunctionListFind(PetscDSList, name, &r));
74: PetscCheck(r, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscDS type: %s", name);
76: PetscTryTypeMethod(prob, destroy);
77: prob->ops->destroy = NULL;
79: PetscCall((*r)(prob));
80: PetscCall(PetscObjectChangeTypeName((PetscObject)prob, name));
81: PetscFunctionReturn(PETSC_SUCCESS);
82: }
84: /*@
85: PetscDSGetType - Gets the `PetscDSType` name (as a string) from the `PetscDS`
87: Not Collective; No Fortran Support
89: Input Parameter:
90: . prob - The `PetscDS`
92: Output Parameter:
93: . name - The `PetscDSType` name
95: Level: intermediate
97: .seealso: `PetscDSType`, `PetscDS`, `PetscDSSetType()`, `PetscDSCreate()`
98: @*/
99: PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name)
100: {
101: PetscFunctionBegin;
103: PetscAssertPointer(name, 2);
104: PetscCall(PetscDSRegisterAll());
105: *name = ((PetscObject)prob)->type_name;
106: PetscFunctionReturn(PETSC_SUCCESS);
107: }
109: static PetscErrorCode PetscDSView_Ascii(PetscDS ds, PetscViewer viewer)
110: {
111: PetscViewerFormat format;
112: const PetscScalar *constants;
113: PetscInt Nf, numConstants, f;
115: PetscFunctionBegin;
116: PetscCall(PetscDSGetNumFields(ds, &Nf));
117: PetscCall(PetscViewerGetFormat(viewer, &format));
118: PetscCall(PetscViewerASCIIPrintf(viewer, "Discrete System with %" PetscInt_FMT " fields\n", Nf));
119: PetscCall(PetscViewerASCIIPushTab(viewer));
120: PetscCall(PetscViewerASCIIPrintf(viewer, " cell total dim %" PetscInt_FMT " total comp %" PetscInt_FMT "\n", ds->totDim, ds->totComp));
121: if (ds->isCohesive) PetscCall(PetscViewerASCIIPrintf(viewer, " cohesive cell\n"));
122: for (f = 0; f < Nf; ++f) {
123: DSBoundary b;
124: PetscObject obj;
125: PetscClassId id;
126: PetscQuadrature q;
127: const char *name;
128: PetscInt Nc, Nq, Nqc;
130: PetscCall(PetscDSGetDiscretization(ds, f, &obj));
131: PetscCall(PetscObjectGetClassId(obj, &id));
132: PetscCall(PetscObjectGetName(obj, &name));
133: PetscCall(PetscViewerASCIIPrintf(viewer, "Field %s", name ? name : "<unknown>"));
134: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
135: if (id == PETSCFE_CLASSID) {
136: PetscCall(PetscFEGetNumComponents((PetscFE)obj, &Nc));
137: PetscCall(PetscFEGetQuadrature((PetscFE)obj, &q));
138: PetscCall(PetscViewerASCIIPrintf(viewer, " FEM"));
139: } else if (id == PETSCFV_CLASSID) {
140: PetscCall(PetscFVGetNumComponents((PetscFV)obj, &Nc));
141: PetscCall(PetscFVGetQuadrature((PetscFV)obj, &q));
142: PetscCall(PetscViewerASCIIPrintf(viewer, " FVM"));
143: } else SETERRQ(PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f);
144: if (Nc > 1) PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " components", Nc));
145: else PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT " component ", Nc));
146: if (ds->implicit[f]) PetscCall(PetscViewerASCIIPrintf(viewer, " (implicit)"));
147: else PetscCall(PetscViewerASCIIPrintf(viewer, " (explicit)"));
148: if (q) {
149: PetscCall(PetscQuadratureGetData(q, NULL, &Nqc, &Nq, NULL, NULL));
150: PetscCall(PetscViewerASCIIPrintf(viewer, " (Nq %" PetscInt_FMT " Nqc %" PetscInt_FMT ")", Nq, Nqc));
151: }
152: PetscCall(PetscViewerASCIIPrintf(viewer, " %" PetscInt_FMT "-jet", ds->jetDegree[f]));
153: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
154: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
155: PetscCall(PetscViewerASCIIPushTab(viewer));
156: if (id == PETSCFE_CLASSID) PetscCall(PetscFEView((PetscFE)obj, viewer));
157: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVView((PetscFV)obj, viewer));
158: PetscCall(PetscViewerASCIIPopTab(viewer));
160: for (b = ds->boundary; b; b = b->next) {
161: char *name;
163: if (b->field != f) continue;
164: PetscCall(PetscViewerASCIIPushTab(viewer));
165: PetscCall(PetscViewerASCIIPrintf(viewer, "Boundary %s (%s) %s\n", b->name, b->lname, DMBoundaryConditionTypes[b->type]));
166: if (!b->Nc) {
167: PetscCall(PetscViewerASCIIPrintf(viewer, " all components\n"));
168: } else {
169: PetscCall(PetscViewerASCIIPrintf(viewer, " components: "));
170: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
171: for (PetscInt c = 0; c < b->Nc; ++c) {
172: if (c > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", "));
173: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->comps[c]));
174: }
175: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
176: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
177: }
178: PetscCall(PetscViewerASCIIPrintf(viewer, " values: "));
179: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
180: for (PetscInt i = 0; i < b->Nv; ++i) {
181: if (i > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", "));
182: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->values[i]));
183: }
184: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
185: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
186: #if defined(__clang__)
187: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat-pedantic")
188: #elif defined(__GNUC__) || defined(__GNUG__)
189: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat")
190: #endif
191: if (b->func) {
192: PetscCall(PetscDLAddr(b->func, &name));
193: if (name) PetscCall(PetscViewerASCIIPrintf(viewer, " func: %s\n", name));
194: else PetscCall(PetscViewerASCIIPrintf(viewer, " func: %p\n", b->func));
195: PetscCall(PetscFree(name));
196: }
197: if (b->func_t) {
198: PetscCall(PetscDLAddr(b->func_t, &name));
199: if (name) PetscCall(PetscViewerASCIIPrintf(viewer, " func_t: %s\n", name));
200: else PetscCall(PetscViewerASCIIPrintf(viewer, " func_t: %p\n", b->func_t));
201: PetscCall(PetscFree(name));
202: }
203: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_END()
204: PetscCall(PetscWeakFormView(b->wf, viewer));
205: PetscCall(PetscViewerASCIIPopTab(viewer));
206: }
207: }
208: PetscCall(PetscDSGetConstants(ds, &numConstants, &constants));
209: if (numConstants) {
210: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " constants\n", numConstants));
211: PetscCall(PetscViewerASCIIPushTab(viewer));
212: for (f = 0; f < numConstants; ++f) PetscCall(PetscViewerASCIIPrintf(viewer, "%g\n", (double)PetscRealPart(constants[f])));
213: PetscCall(PetscViewerASCIIPopTab(viewer));
214: }
215: PetscCall(PetscWeakFormView(ds->wf, viewer));
216: PetscCall(PetscViewerASCIIPopTab(viewer));
217: PetscFunctionReturn(PETSC_SUCCESS);
218: }
220: /*@
221: PetscDSViewFromOptions - View a `PetscDS` based on values in the options database
223: Collective
225: Input Parameters:
226: + A - the `PetscDS` object
227: . obj - Optional object that provides the options prefix used in the search of the options database
228: - name - command line option
230: Options Database Key:
231: . -name [viewertype][:...] - option name and values. See `PetscObjectViewFromOptions()` for the possible arguments
233: Level: intermediate
235: .seealso: `PetscDSType`, `PetscDS`, `PetscDSView()`, `PetscObjectViewFromOptions()`, `PetscDSCreate()`
236: @*/
237: PetscErrorCode PetscDSViewFromOptions(PetscDS A, PetscObject obj, const char name[])
238: {
239: PetscFunctionBegin;
241: PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
242: PetscFunctionReturn(PETSC_SUCCESS);
243: }
245: /*@
246: PetscDSView - Views a `PetscDS`
248: Collective
250: Input Parameters:
251: + prob - the `PetscDS` object to view
252: - v - the viewer
254: Level: developer
256: .seealso: `PetscDSType`, `PetscDS`, `PetscViewer`, `PetscDSDestroy()`, `PetscDSViewFromOptions()`
257: @*/
258: PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v)
259: {
260: PetscBool isascii;
262: PetscFunctionBegin;
264: if (!v) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)prob), &v));
266: PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &isascii));
267: if (isascii) PetscCall(PetscDSView_Ascii(prob, v));
268: PetscTryTypeMethod(prob, view, v);
269: PetscFunctionReturn(PETSC_SUCCESS);
270: }
272: /*@
273: PetscDSSetFromOptions - sets parameters in a `PetscDS` from the options database
275: Collective
277: Input Parameter:
278: . prob - the `PetscDS` object to set options for
280: Options Database Keys:
281: + -petscds_type type - Set the `PetscDS` type
282: . -petscds_view - View the `PetscDS`
283: . -petscds_jac_pre (true|false) - Turn formation of a separate Jacobian preconditioner on or off
284: . -bc_NAME ids - comma separated list of label ids for the boundary condition NAME
285: - -bc_NAME_comp comps - comma separated list of field components to constrain for the boundary condition NAME
287: Level: intermediate
289: .seealso: `PetscDS`, `PetscDSView()`
290: @*/
291: PetscErrorCode PetscDSSetFromOptions(PetscDS prob)
292: {
293: DSBoundary b;
294: const char *defaultType;
295: char name[256];
296: PetscBool flg;
298: PetscFunctionBegin;
300: if (!((PetscObject)prob)->type_name) {
301: defaultType = PETSCDSBASIC;
302: } else {
303: defaultType = ((PetscObject)prob)->type_name;
304: }
305: PetscCall(PetscDSRegisterAll());
307: PetscObjectOptionsBegin((PetscObject)prob);
308: for (b = prob->boundary; b; b = b->next) {
309: char optname[1024];
310: PetscInt ids[1024], len = 1024;
311: PetscBool flg;
313: PetscCall(PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name));
314: PetscCall(PetscMemzero(ids, sizeof(ids)));
315: PetscCall(PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg));
316: if (flg) {
317: b->Nv = len;
318: PetscCall(PetscFree(b->values));
319: PetscCall(PetscMalloc1(len, &b->values));
320: PetscCall(PetscArraycpy(b->values, ids, len));
321: PetscCall(PetscWeakFormRewriteKeys(b->wf, b->label, len, b->values));
322: }
323: len = 1024;
324: PetscCall(PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name));
325: PetscCall(PetscMemzero(ids, sizeof(ids)));
326: PetscCall(PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg));
327: if (flg) {
328: b->Nc = len;
329: PetscCall(PetscFree(b->comps));
330: PetscCall(PetscMalloc1(len, &b->comps));
331: PetscCall(PetscArraycpy(b->comps, ids, len));
332: }
333: }
334: PetscCall(PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg));
335: if (flg) {
336: PetscCall(PetscDSSetType(prob, name));
337: } else if (!((PetscObject)prob)->type_name) {
338: PetscCall(PetscDSSetType(prob, defaultType));
339: }
340: PetscCall(PetscOptionsBool("-petscds_jac_pre", "Discrete System", "PetscDSUseJacobianPreconditioner", prob->useJacPre, &prob->useJacPre, &flg));
341: PetscCall(PetscOptionsBool("-petscds_force_quad", "Discrete System", "PetscDSSetForceQuad", prob->forceQuad, &prob->forceQuad, &flg));
342: PetscCall(PetscOptionsInt("-petscds_print_integrate", "Discrete System", "", prob->printIntegrate, &prob->printIntegrate, NULL));
343: PetscTryTypeMethod(prob, setfromoptions);
344: /* process any options handlers added with PetscObjectAddOptionsHandler() */
345: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)prob, PetscOptionsObject));
346: PetscOptionsEnd();
347: if (prob->Nf) PetscCall(PetscDSViewFromOptions(prob, NULL, "-petscds_view"));
348: PetscFunctionReturn(PETSC_SUCCESS);
349: }
351: /*@
352: PetscDSSetUp - Construct data structures for the `PetscDS`
354: Collective
356: Input Parameter:
357: . prob - the `PetscDS` object to setup
359: Level: developer
361: .seealso: `PetscDS`, `PetscDSView()`, `PetscDSDestroy()`
362: @*/
363: PetscErrorCode PetscDSSetUp(PetscDS prob)
364: {
365: const PetscInt Nf = prob->Nf;
366: PetscBool hasH = PETSC_FALSE;
367: PetscInt maxOrder[4] = {-2, -2, -2, -2};
368: PetscInt dim, dimEmbed, NbMax = 0, NcMax = 0, NqMax = 0, NsMax = 1, f;
370: PetscFunctionBegin;
372: if (prob->setup) PetscFunctionReturn(PETSC_SUCCESS);
373: /* Calculate sizes */
374: PetscCall(PetscDSGetSpatialDimension(prob, &dim));
375: PetscCall(PetscDSGetCoordinateDimension(prob, &dimEmbed));
376: prob->totDim = prob->totComp = 0;
377: PetscCall(PetscMalloc2(Nf, &prob->Nc, Nf, &prob->Nb));
378: PetscCall(PetscCalloc2(Nf + 1, &prob->off, Nf + 1, &prob->offDer));
379: 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]));
380: PetscCall(PetscMalloc2(Nf, &prob->T, Nf, &prob->Tf));
381: if (prob->forceQuad) {
382: // Note: This assumes we have one kind of cell at each dimension.
383: // We can fix this by having quadrature hold the celltype
384: PetscQuadrature maxQuad[4] = {NULL, NULL, NULL, NULL};
386: for (f = 0; f < Nf; ++f) {
387: PetscObject obj;
388: PetscClassId id;
389: PetscQuadrature q = NULL, fq = NULL;
390: PetscInt dim = -1, order = -1, forder = -1;
392: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
393: if (!obj) continue;
394: PetscCall(PetscObjectGetClassId(obj, &id));
395: if (id == PETSCFE_CLASSID) {
396: PetscFE fe = (PetscFE)obj;
398: PetscCall(PetscFEGetQuadrature(fe, &q));
399: PetscCall(PetscFEGetFaceQuadrature(fe, &fq));
400: } else if (id == PETSCFV_CLASSID) {
401: PetscFV fv = (PetscFV)obj;
403: PetscCall(PetscFVGetQuadrature(fv, &q));
404: }
405: if (q) {
406: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
407: PetscCall(PetscQuadratureGetOrder(q, &order));
408: if (order > maxOrder[dim]) {
409: maxOrder[dim] = order;
410: maxQuad[dim] = q;
411: }
412: }
413: if (fq) {
414: PetscCall(PetscQuadratureGetData(fq, &dim, NULL, NULL, NULL, NULL));
415: PetscCall(PetscQuadratureGetOrder(fq, &forder));
416: if (forder > maxOrder[dim]) {
417: maxOrder[dim] = forder;
418: maxQuad[dim] = fq;
419: }
420: }
421: }
422: for (f = 0; f < Nf; ++f) {
423: PetscObject obj;
424: PetscClassId id;
425: PetscQuadrature q;
426: PetscInt dim;
428: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
429: if (!obj) continue;
430: PetscCall(PetscObjectGetClassId(obj, &id));
431: if (id == PETSCFE_CLASSID) {
432: PetscFE fe = (PetscFE)obj;
434: PetscCall(PetscFEGetQuadrature(fe, &q));
435: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
436: PetscCall(PetscFESetQuadrature(fe, maxQuad[dim]));
437: PetscCall(PetscFESetFaceQuadrature(fe, dim ? maxQuad[dim - 1] : NULL));
438: } else if (id == PETSCFV_CLASSID) {
439: PetscFV fv = (PetscFV)obj;
441: PetscCall(PetscFVGetQuadrature(fv, &q));
442: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
443: PetscCall(PetscFVSetQuadrature(fv, maxQuad[dim]));
444: }
445: }
446: }
447: for (f = 0; f < Nf; ++f) {
448: PetscObject obj;
449: PetscClassId id;
450: PetscQuadrature q = NULL;
451: PetscInt Nq = 0, Nb, Nc;
453: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
454: if (prob->jetDegree[f] > 1) hasH = PETSC_TRUE;
455: if (!obj) {
456: /* Empty mesh */
457: Nb = Nc = 0;
458: prob->T[f] = prob->Tf[f] = NULL;
459: } else {
460: PetscCall(PetscObjectGetClassId(obj, &id));
461: if (id == PETSCFE_CLASSID) {
462: PetscFE fe = (PetscFE)obj;
464: PetscCall(PetscFEGetQuadrature(fe, &q));
465: {
466: PetscQuadrature fq;
467: PetscInt dim, order;
469: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
470: PetscCall(PetscQuadratureGetOrder(q, &order));
471: if (maxOrder[dim] < 0) maxOrder[dim] = order;
472: 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]);
473: PetscCall(PetscFEGetFaceQuadrature(fe, &fq));
474: if (fq) {
475: PetscCall(PetscQuadratureGetData(fq, &dim, NULL, NULL, NULL, NULL));
476: PetscCall(PetscQuadratureGetOrder(fq, &order));
477: if (maxOrder[dim] < 0) maxOrder[dim] = order;
478: 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]);
479: }
480: }
481: PetscCall(PetscFEGetDimension(fe, &Nb));
482: PetscCall(PetscFEGetNumComponents(fe, &Nc));
483: PetscCall(PetscFEGetCellTabulation(fe, prob->jetDegree[f], &prob->T[f]));
484: PetscCall(PetscFEGetFaceTabulation(fe, prob->jetDegree[f], &prob->Tf[f]));
485: } else if (id == PETSCFV_CLASSID) {
486: PetscFV fv = (PetscFV)obj;
488: PetscCall(PetscFVGetQuadrature(fv, &q));
489: PetscCall(PetscFVGetNumComponents(fv, &Nc));
490: Nb = Nc;
491: PetscCall(PetscFVGetCellTabulation(fv, &prob->T[f]));
492: /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */
493: } else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f);
494: }
495: prob->Nc[f] = Nc;
496: prob->Nb[f] = Nb;
497: prob->off[f + 1] = Nc + prob->off[f];
498: prob->offDer[f + 1] = Nc * dim + prob->offDer[f];
499: prob->offCohesive[0][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[0][f];
500: prob->offDerCohesive[0][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[0][f];
501: prob->offCohesive[1][f] = (prob->cohesive[f] ? 0 : Nc) + prob->offCohesive[0][f];
502: prob->offDerCohesive[1][f] = (prob->cohesive[f] ? 0 : Nc) * dimEmbed + prob->offDerCohesive[0][f];
503: prob->offCohesive[2][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[2][f];
504: prob->offDerCohesive[2][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[2][f];
505: if (q) PetscCall(PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL));
506: NqMax = PetscMax(NqMax, Nq);
507: NbMax = PetscMax(NbMax, Nb);
508: NcMax = PetscMax(NcMax, Nc);
509: prob->totDim += Nb;
510: prob->totComp += Nc;
511: /* There are two faces for all fields on a cohesive cell, except for cohesive fields */
512: if (prob->isCohesive && !prob->cohesive[f]) prob->totDim += Nb;
513: }
514: prob->offCohesive[1][Nf] = prob->offCohesive[0][Nf];
515: prob->offDerCohesive[1][Nf] = prob->offDerCohesive[0][Nf];
516: /* Allocate works space */
517: NsMax = 2; /* A non-cohesive discretizations can be used on a cohesive cell, so we need this extra workspace for all DS */
518: 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));
519: PetscCall(PetscMalloc5(dimEmbed, &prob->x, NbMax * NcMax, &prob->basisReal, NbMax * NcMax * dimEmbed, &prob->basisDerReal, NbMax * NcMax, &prob->testReal, NbMax * NcMax * dimEmbed, &prob->testDerReal));
520: 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,
521: &prob->g2, NsMax * NsMax * NqMax * NcMax * NcMax * dimEmbed * dimEmbed, &prob->g3));
522: PetscTryTypeMethod(prob, setup);
523: prob->setup = PETSC_TRUE;
524: PetscFunctionReturn(PETSC_SUCCESS);
525: }
527: static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob)
528: {
529: PetscFunctionBegin;
530: PetscCall(PetscFree2(prob->Nc, prob->Nb));
531: PetscCall(PetscFree2(prob->off, prob->offDer));
532: PetscCall(PetscFree6(prob->offCohesive[0], prob->offCohesive[1], prob->offCohesive[2], prob->offDerCohesive[0], prob->offDerCohesive[1], prob->offDerCohesive[2]));
533: PetscCall(PetscFree2(prob->T, prob->Tf));
534: PetscCall(PetscFree3(prob->u, prob->u_t, prob->u_x));
535: PetscCall(PetscFree5(prob->x, prob->basisReal, prob->basisDerReal, prob->testReal, prob->testDerReal));
536: PetscCall(PetscFree6(prob->f0, prob->f1, prob->g0, prob->g1, prob->g2, prob->g3));
537: PetscFunctionReturn(PETSC_SUCCESS);
538: }
540: static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew)
541: {
542: PetscObject *tmpd;
543: PetscBool *tmpi;
544: PetscInt *tmpk;
545: PetscBool *tmpc;
546: PetscPointFn **tmpup;
547: PetscSimplePointFn **tmpexactSol, **tmpexactSol_t, **tmplowerBound, **tmpupperBound;
548: void **tmpexactCtx, **tmpexactCtx_t, **tmplowerCtx, **tmpupperCtx;
549: void **tmpctx;
550: PetscInt Nf = prob->Nf, f;
552: PetscFunctionBegin;
553: if (Nf >= NfNew) PetscFunctionReturn(PETSC_SUCCESS);
554: prob->setup = PETSC_FALSE;
555: PetscCall(PetscDSDestroyStructs_Static(prob));
556: PetscCall(PetscMalloc4(NfNew, &tmpd, NfNew, &tmpi, NfNew, &tmpc, NfNew, &tmpk));
557: for (f = 0; f < Nf; ++f) {
558: tmpd[f] = prob->disc[f];
559: tmpi[f] = prob->implicit[f];
560: tmpc[f] = prob->cohesive[f];
561: tmpk[f] = prob->jetDegree[f];
562: }
563: for (f = Nf; f < NfNew; ++f) {
564: tmpd[f] = NULL;
565: tmpi[f] = PETSC_TRUE, tmpc[f] = PETSC_FALSE;
566: tmpk[f] = 1;
567: }
568: PetscCall(PetscFree4(prob->disc, prob->implicit, prob->cohesive, prob->jetDegree));
569: PetscCall(PetscWeakFormSetNumFields(prob->wf, NfNew));
570: prob->Nf = NfNew;
571: prob->disc = tmpd;
572: prob->implicit = tmpi;
573: prob->cohesive = tmpc;
574: prob->jetDegree = tmpk;
575: PetscCall(PetscCalloc2(NfNew, &tmpup, NfNew, &tmpctx));
576: for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f];
577: for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f];
578: for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL;
579: for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL;
580: PetscCall(PetscFree2(prob->update, prob->ctx));
581: prob->update = tmpup;
582: prob->ctx = tmpctx;
583: PetscCall(PetscCalloc4(NfNew, &tmpexactSol, NfNew, &tmpexactCtx, NfNew, &tmpexactSol_t, NfNew, &tmpexactCtx_t));
584: PetscCall(PetscCalloc4(NfNew, &tmplowerBound, NfNew, &tmplowerCtx, NfNew, &tmpupperBound, NfNew, &tmpupperCtx));
585: for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f];
586: for (f = 0; f < Nf; ++f) tmpexactCtx[f] = prob->exactCtx[f];
587: for (f = 0; f < Nf; ++f) tmpexactSol_t[f] = prob->exactSol_t[f];
588: for (f = 0; f < Nf; ++f) tmpexactCtx_t[f] = prob->exactCtx_t[f];
589: for (f = 0; f < Nf; ++f) tmplowerBound[f] = prob->lowerBound[f];
590: for (f = 0; f < Nf; ++f) tmplowerCtx[f] = prob->lowerCtx[f];
591: for (f = 0; f < Nf; ++f) tmpupperBound[f] = prob->upperBound[f];
592: for (f = 0; f < Nf; ++f) tmpupperCtx[f] = prob->upperCtx[f];
593: for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL;
594: for (f = Nf; f < NfNew; ++f) tmpexactCtx[f] = NULL;
595: for (f = Nf; f < NfNew; ++f) tmpexactSol_t[f] = NULL;
596: for (f = Nf; f < NfNew; ++f) tmpexactCtx_t[f] = NULL;
597: for (f = Nf; f < NfNew; ++f) tmplowerBound[f] = NULL;
598: for (f = Nf; f < NfNew; ++f) tmplowerCtx[f] = NULL;
599: for (f = Nf; f < NfNew; ++f) tmpupperBound[f] = NULL;
600: for (f = Nf; f < NfNew; ++f) tmpupperCtx[f] = NULL;
601: PetscCall(PetscFree4(prob->exactSol, prob->exactCtx, prob->exactSol_t, prob->exactCtx_t));
602: PetscCall(PetscFree4(prob->lowerBound, prob->lowerCtx, prob->upperBound, prob->upperCtx));
603: prob->exactSol = tmpexactSol;
604: prob->exactCtx = tmpexactCtx;
605: prob->exactSol_t = tmpexactSol_t;
606: prob->exactCtx_t = tmpexactCtx_t;
607: prob->lowerBound = tmplowerBound;
608: prob->lowerCtx = tmplowerCtx;
609: prob->upperBound = tmpupperBound;
610: prob->upperCtx = tmpupperCtx;
611: PetscFunctionReturn(PETSC_SUCCESS);
612: }
614: /*@
615: PetscDSDestroy - Destroys a `PetscDS` object
617: Collective
619: Input Parameter:
620: . ds - the `PetscDS` object to destroy
622: Level: developer
624: .seealso: `PetscDSView()`
625: @*/
626: PetscErrorCode PetscDSDestroy(PetscDS *ds)
627: {
628: PetscFunctionBegin;
629: if (!*ds) PetscFunctionReturn(PETSC_SUCCESS);
632: if (--((PetscObject)*ds)->refct > 0) {
633: *ds = NULL;
634: PetscFunctionReturn(PETSC_SUCCESS);
635: }
636: ((PetscObject)*ds)->refct = 0;
637: if ((*ds)->subprobs) {
638: PetscInt dim;
640: PetscCall(PetscDSGetSpatialDimension(*ds, &dim));
641: for (PetscInt d = 0; d < dim; ++d) PetscCall(PetscDSDestroy(&(*ds)->subprobs[d]));
642: }
643: PetscCall(PetscFree((*ds)->subprobs));
644: PetscCall(PetscDSDestroyStructs_Static(*ds));
645: for (PetscInt f = 0; f < (*ds)->Nf; ++f) PetscCall(PetscObjectDereference((*ds)->disc[f]));
646: PetscCall(PetscFree4((*ds)->disc, (*ds)->implicit, (*ds)->cohesive, (*ds)->jetDegree));
647: PetscCall(PetscWeakFormDestroy(&(*ds)->wf));
648: PetscCall(PetscFree2((*ds)->update, (*ds)->ctx));
649: PetscCall(PetscFree4((*ds)->exactSol, (*ds)->exactCtx, (*ds)->exactSol_t, (*ds)->exactCtx_t));
650: PetscCall(PetscFree4((*ds)->lowerBound, (*ds)->lowerCtx, (*ds)->upperBound, (*ds)->upperCtx));
651: PetscTryTypeMethod(*ds, destroy);
652: PetscCall(PetscDSDestroyBoundary(*ds));
653: PetscCall(PetscFree((*ds)->constants));
654: for (PetscInt c = 0; c < DM_NUM_POLYTOPES; ++c) {
655: const PetscInt Na = DMPolytopeTypeGetNumArrangements((DMPolytopeType)c);
656: if ((*ds)->quadPerm[c])
657: for (PetscInt o = 0; o < Na; ++o) PetscCall(ISDestroy(&(*ds)->quadPerm[c][o]));
658: PetscCall(PetscFree((*ds)->quadPerm[c]));
659: (*ds)->quadPerm[c] = NULL;
660: }
661: PetscCall(PetscHeaderDestroy(ds));
662: PetscFunctionReturn(PETSC_SUCCESS);
663: }
665: /*@
666: PetscDSCreate - Creates an empty `PetscDS` object. The type can then be set with `PetscDSSetType()`.
668: Collective
670: Input Parameter:
671: . comm - The communicator for the `PetscDS` object
673: Output Parameter:
674: . ds - The `PetscDS` object
676: Level: beginner
678: .seealso: `PetscDS`, `PetscDSSetType()`, `PETSCDSBASIC`, `PetscDSType`, `PetscDSDestroy()`
679: @*/
680: PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *ds)
681: {
682: PetscDS p;
684: PetscFunctionBegin;
685: PetscAssertPointer(ds, 2);
686: PetscCall(PetscDSInitializePackage());
688: PetscCall(PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView));
689: p->Nf = 0;
690: p->setup = PETSC_FALSE;
691: p->numConstants = 0;
692: p->numFuncConstants = 3; // Row and col fields, cell size
693: p->dimEmbed = -1;
694: p->useJacPre = PETSC_TRUE;
695: p->forceQuad = PETSC_TRUE;
696: PetscCall(PetscMalloc1(p->numConstants + p->numFuncConstants, &p->constants));
697: PetscCall(PetscWeakFormCreate(comm, &p->wf));
698: PetscCall(PetscArrayzero(p->quadPerm, DM_NUM_POLYTOPES));
699: *ds = p;
700: PetscFunctionReturn(PETSC_SUCCESS);
701: }
703: /*@
704: PetscDSGetNumFields - Returns the number of fields in the `PetscDS`
706: Not Collective
708: Input Parameter:
709: . prob - The `PetscDS` object
711: Output Parameter:
712: . Nf - The number of fields
714: Level: beginner
716: .seealso: `PetscDS`, `PetscDSGetSpatialDimension()`, `PetscDSCreate()`
717: @*/
718: PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf)
719: {
720: PetscFunctionBegin;
722: PetscAssertPointer(Nf, 2);
723: *Nf = prob->Nf;
724: PetscFunctionReturn(PETSC_SUCCESS);
725: }
727: /*@
728: PetscDSGetSpatialDimension - Returns the spatial dimension of the `PetscDS`, meaning the topological dimension of the discretizations
730: Not Collective
732: Input Parameter:
733: . prob - The `PetscDS` object
735: Output Parameter:
736: . dim - The spatial dimension
738: Level: beginner
740: .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
741: @*/
742: PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim)
743: {
744: PetscFunctionBegin;
746: PetscAssertPointer(dim, 2);
747: *dim = 0;
748: if (prob->Nf) {
749: PetscObject obj;
750: PetscClassId id;
752: PetscCall(PetscDSGetDiscretization(prob, 0, &obj));
753: if (obj) {
754: PetscCall(PetscObjectGetClassId(obj, &id));
755: if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetSpatialDimension((PetscFE)obj, dim));
756: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetSpatialDimension((PetscFV)obj, dim));
757: else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
758: }
759: }
760: PetscFunctionReturn(PETSC_SUCCESS);
761: }
763: /*@
764: PetscDSGetCoordinateDimension - Returns the coordinate dimension of the `PetscDS`, meaning the dimension of the space into which the discretiaztions are embedded
766: Not Collective
768: Input Parameter:
769: . prob - The `PetscDS` object
771: Output Parameter:
772: . dimEmbed - The coordinate dimension
774: Level: beginner
776: .seealso: `PetscDS`, `PetscDSSetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
777: @*/
778: PetscErrorCode PetscDSGetCoordinateDimension(PetscDS prob, PetscInt *dimEmbed)
779: {
780: PetscFunctionBegin;
782: PetscAssertPointer(dimEmbed, 2);
783: PetscCheck(prob->dimEmbed >= 0, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONGSTATE, "No coordinate dimension set for this DS");
784: *dimEmbed = prob->dimEmbed;
785: PetscFunctionReturn(PETSC_SUCCESS);
786: }
788: /*@
789: PetscDSSetCoordinateDimension - Set the coordinate dimension of the `PetscDS`, meaning the dimension of the space into which the discretiaztions are embedded
791: Logically Collective
793: Input Parameters:
794: + prob - The `PetscDS` object
795: - dimEmbed - The coordinate dimension
797: Level: beginner
799: .seealso: `PetscDS`, `PetscDSGetCoordinateDimension()`, `PetscDSGetSpatialDimension()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
800: @*/
801: PetscErrorCode PetscDSSetCoordinateDimension(PetscDS prob, PetscInt dimEmbed)
802: {
803: PetscFunctionBegin;
805: PetscCheck(dimEmbed >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Coordinate dimension must be non-negative, not %" PetscInt_FMT, dimEmbed);
806: prob->dimEmbed = dimEmbed;
807: PetscFunctionReturn(PETSC_SUCCESS);
808: }
810: /*@
811: PetscDSGetForceQuad - Returns the flag to force matching quadratures among the field discretizations
813: Not collective
815: Input Parameter:
816: . ds - The `PetscDS` object
818: Output Parameter:
819: . forceQuad - The flag
821: Level: intermediate
823: .seealso: `PetscDS`, `PetscDSSetForceQuad()`, `PetscDSGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
824: @*/
825: PetscErrorCode PetscDSGetForceQuad(PetscDS ds, PetscBool *forceQuad)
826: {
827: PetscFunctionBegin;
829: PetscAssertPointer(forceQuad, 2);
830: *forceQuad = ds->forceQuad;
831: PetscFunctionReturn(PETSC_SUCCESS);
832: }
834: /*@
835: PetscDSSetForceQuad - Set the flag to force matching quadratures among the field discretizations
837: Logically collective on ds
839: Input Parameters:
840: + ds - The `PetscDS` object
841: - forceQuad - The flag
843: Level: intermediate
845: .seealso: `PetscDS`, `PetscDSGetForceQuad()`, `PetscDSGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
846: @*/
847: PetscErrorCode PetscDSSetForceQuad(PetscDS ds, PetscBool forceQuad)
848: {
849: PetscFunctionBegin;
851: ds->forceQuad = forceQuad;
852: PetscFunctionReturn(PETSC_SUCCESS);
853: }
855: /*@
856: PetscDSIsCohesive - Returns the flag indicating that this `PetscDS` is for a cohesive cell
858: Not Collective
860: Input Parameter:
861: . ds - The `PetscDS` object
863: Output Parameter:
864: . isCohesive - The flag
866: Level: developer
868: .seealso: `PetscDS`, `PetscDSGetNumCohesive()`, `PetscDSGetCohesive()`, `PetscDSSetCohesive()`, `PetscDSCreate()`
869: @*/
870: PetscErrorCode PetscDSIsCohesive(PetscDS ds, PetscBool *isCohesive)
871: {
872: PetscFunctionBegin;
874: PetscAssertPointer(isCohesive, 2);
875: *isCohesive = ds->isCohesive;
876: PetscFunctionReturn(PETSC_SUCCESS);
877: }
879: /*@
880: PetscDSGetNumCohesive - Returns the number of cohesive fields, meaning those defined on the interior of a cohesive cell
882: Not Collective
884: Input Parameter:
885: . ds - The `PetscDS` object
887: Output Parameter:
888: . numCohesive - The number of cohesive fields
890: Level: developer
892: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSCreate()`
893: @*/
894: PetscErrorCode PetscDSGetNumCohesive(PetscDS ds, PetscInt *numCohesive)
895: {
896: PetscFunctionBegin;
898: PetscAssertPointer(numCohesive, 2);
899: *numCohesive = 0;
900: for (PetscInt f = 0; f < ds->Nf; ++f) *numCohesive += ds->cohesive[f] ? 1 : 0;
901: PetscFunctionReturn(PETSC_SUCCESS);
902: }
904: /*@
905: PetscDSGetCohesive - Returns the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell
907: Not Collective
909: Input Parameters:
910: + ds - The `PetscDS` object
911: - f - The field index
913: Output Parameter:
914: . isCohesive - The flag
916: Level: developer
918: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
919: @*/
920: PetscErrorCode PetscDSGetCohesive(PetscDS ds, PetscInt f, PetscBool *isCohesive)
921: {
922: PetscFunctionBegin;
924: PetscAssertPointer(isCohesive, 3);
925: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
926: *isCohesive = ds->cohesive[f];
927: PetscFunctionReturn(PETSC_SUCCESS);
928: }
930: /*@
931: PetscDSSetCohesive - Set the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell
933: Not Collective
935: Input Parameters:
936: + ds - The `PetscDS` object
937: . f - The field index
938: - isCohesive - The flag for a cohesive field
940: Level: developer
942: .seealso: `PetscDS`, `PetscDSGetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
943: @*/
944: PetscErrorCode PetscDSSetCohesive(PetscDS ds, PetscInt f, PetscBool isCohesive)
945: {
946: PetscFunctionBegin;
948: 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);
949: ds->cohesive[f] = isCohesive;
950: ds->isCohesive = PETSC_FALSE;
951: for (PetscInt i = 0; i < ds->Nf; ++i) ds->isCohesive = ds->isCohesive || ds->cohesive[f] ? PETSC_TRUE : PETSC_FALSE;
952: PetscFunctionReturn(PETSC_SUCCESS);
953: }
955: /*@
956: PetscDSGetTotalDimension - Returns the total size of the approximation space for this system
958: Not Collective
960: Input Parameter:
961: . prob - The `PetscDS` object
963: Output Parameter:
964: . dim - The total problem dimension
966: Level: beginner
968: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
969: @*/
970: PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim)
971: {
972: PetscFunctionBegin;
974: PetscCall(PetscDSSetUp(prob));
975: PetscAssertPointer(dim, 2);
976: *dim = prob->totDim;
977: PetscFunctionReturn(PETSC_SUCCESS);
978: }
980: /*@
981: PetscDSGetTotalComponents - Returns the total number of components in this system
983: Not Collective
985: Input Parameter:
986: . prob - The `PetscDS` object
988: Output Parameter:
989: . Nc - The total number of components
991: Level: beginner
993: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
994: @*/
995: PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc)
996: {
997: PetscFunctionBegin;
999: PetscCall(PetscDSSetUp(prob));
1000: PetscAssertPointer(Nc, 2);
1001: *Nc = prob->totComp;
1002: PetscFunctionReturn(PETSC_SUCCESS);
1003: }
1005: /*@
1006: PetscDSGetDiscretization - Returns the discretization object for the given field
1008: Not Collective
1010: Input Parameters:
1011: + prob - The `PetscDS` object
1012: - f - The field number
1014: Output Parameter:
1015: . disc - The discretization object, this can be a `PetscFE` or a `PetscFV`
1017: Level: beginner
1019: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1020: @*/
1021: PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc)
1022: {
1023: PetscFunctionBeginHot;
1025: PetscAssertPointer(disc, 3);
1026: 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);
1027: *disc = prob->disc[f];
1028: PetscFunctionReturn(PETSC_SUCCESS);
1029: }
1031: /*@
1032: PetscDSSetDiscretization - Sets the discretization object for the given field
1034: Not Collective
1036: Input Parameters:
1037: + prob - The `PetscDS` object
1038: . f - The field number
1039: - disc - The discretization object, this can be a `PetscFE` or a `PetscFV`
1041: Level: beginner
1043: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSGetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1044: @*/
1045: PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc)
1046: {
1047: PetscFunctionBegin;
1049: if (disc) PetscAssertPointer(disc, 3);
1050: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1051: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
1052: PetscCall(PetscObjectDereference(prob->disc[f]));
1053: prob->disc[f] = disc;
1054: PetscCall(PetscObjectReference(disc));
1055: if (disc) {
1056: PetscClassId id;
1058: PetscCall(PetscObjectGetClassId(disc, &id));
1059: if (id == PETSCFE_CLASSID) {
1060: PetscCall(PetscDSSetImplicit(prob, f, PETSC_TRUE));
1061: } else if (id == PETSCFV_CLASSID) {
1062: PetscCall(PetscDSSetImplicit(prob, f, PETSC_FALSE));
1063: }
1064: PetscCall(PetscDSSetJetDegree(prob, f, 1));
1065: }
1066: PetscFunctionReturn(PETSC_SUCCESS);
1067: }
1069: /*@
1070: PetscDSGetWeakForm - Returns the weak form object from within the `PetscDS`
1072: Not Collective
1074: Input Parameter:
1075: . ds - The `PetscDS` object
1077: Output Parameter:
1078: . wf - The weak form object
1080: Level: beginner
1082: .seealso: `PetscWeakForm`, `PetscDSSetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1083: @*/
1084: PetscErrorCode PetscDSGetWeakForm(PetscDS ds, PetscWeakForm *wf)
1085: {
1086: PetscFunctionBegin;
1088: PetscAssertPointer(wf, 2);
1089: *wf = ds->wf;
1090: PetscFunctionReturn(PETSC_SUCCESS);
1091: }
1093: /*@
1094: PetscDSSetWeakForm - Sets the weak form object to be used by the `PetscDS`
1096: Not Collective
1098: Input Parameters:
1099: + ds - The `PetscDS` object
1100: - wf - The weak form object
1102: Level: beginner
1104: .seealso: `PetscWeakForm`, `PetscDSGetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1105: @*/
1106: PetscErrorCode PetscDSSetWeakForm(PetscDS ds, PetscWeakForm wf)
1107: {
1108: PetscFunctionBegin;
1111: PetscCall(PetscObjectDereference((PetscObject)ds->wf));
1112: ds->wf = wf;
1113: PetscCall(PetscObjectReference((PetscObject)wf));
1114: PetscCall(PetscWeakFormSetNumFields(wf, ds->Nf));
1115: PetscFunctionReturn(PETSC_SUCCESS);
1116: }
1118: /*@
1119: PetscDSAddDiscretization - Adds a discretization object
1121: Not Collective
1123: Input Parameters:
1124: + prob - The `PetscDS` object
1125: - disc - The discretization object, this can be a `PetscFE` or `PetscFV`
1127: Level: beginner
1129: .seealso: `PetscWeakForm`, `PetscFE`, `PetscFV`, `PetscDSGetDiscretization()`, `PetscDSSetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1130: @*/
1131: PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc)
1132: {
1133: PetscFunctionBegin;
1134: PetscCall(PetscDSSetDiscretization(prob, prob->Nf, disc));
1135: PetscFunctionReturn(PETSC_SUCCESS);
1136: }
1138: /*@
1139: PetscDSGetQuadrature - Returns the quadrature, which must agree for all fields in the `PetscDS`
1141: Not Collective
1143: Input Parameter:
1144: . prob - The `PetscDS` object
1146: Output Parameter:
1147: . q - The quadrature object
1149: Level: intermediate
1151: .seealso: `PetscDS`, `PetscQuadrature`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1152: @*/
1153: PetscErrorCode PetscDSGetQuadrature(PetscDS prob, PetscQuadrature *q)
1154: {
1155: PetscObject obj;
1156: PetscClassId id;
1158: PetscFunctionBegin;
1159: *q = NULL;
1160: if (!prob->Nf) PetscFunctionReturn(PETSC_SUCCESS);
1161: PetscCall(PetscDSGetDiscretization(prob, 0, &obj));
1162: PetscCall(PetscObjectGetClassId(obj, &id));
1163: if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetQuadrature((PetscFE)obj, q));
1164: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetQuadrature((PetscFV)obj, q));
1165: else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
1166: PetscFunctionReturn(PETSC_SUCCESS);
1167: }
1169: /*@
1170: PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for `TSARKIMEX`
1172: Not Collective
1174: Input Parameters:
1175: + prob - The `PetscDS` object
1176: - f - The field number
1178: Output Parameter:
1179: . implicit - The flag indicating what kind of solve to use for this field
1181: Level: developer
1183: .seealso: `TSARKIMEX`, `PetscDS`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1184: @*/
1185: PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit)
1186: {
1187: PetscFunctionBegin;
1189: PetscAssertPointer(implicit, 3);
1190: 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);
1191: *implicit = prob->implicit[f];
1192: PetscFunctionReturn(PETSC_SUCCESS);
1193: }
1195: /*@
1196: PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for `TSARKIMEX`
1198: Not Collective
1200: Input Parameters:
1201: + prob - The `PetscDS` object
1202: . f - The field number
1203: - implicit - The flag indicating what kind of solve to use for this field
1205: Level: developer
1207: .seealso: `TSARKIMEX`, `PetscDSGetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1208: @*/
1209: PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit)
1210: {
1211: PetscFunctionBegin;
1213: 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);
1214: prob->implicit[f] = implicit;
1215: PetscFunctionReturn(PETSC_SUCCESS);
1216: }
1218: /*@
1219: PetscDSGetJetDegree - Returns the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate.
1221: Not Collective
1223: Input Parameters:
1224: + ds - The `PetscDS` object
1225: - f - The field number
1227: Output Parameter:
1228: . k - The highest derivative we need to tabulate
1230: Level: developer
1232: .seealso: `PetscDS`, `PetscDSSetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1233: @*/
1234: PetscErrorCode PetscDSGetJetDegree(PetscDS ds, PetscInt f, PetscInt *k)
1235: {
1236: PetscFunctionBegin;
1238: PetscAssertPointer(k, 3);
1239: 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);
1240: *k = ds->jetDegree[f];
1241: PetscFunctionReturn(PETSC_SUCCESS);
1242: }
1244: /*@
1245: PetscDSSetJetDegree - Set the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate.
1247: Not Collective
1249: Input Parameters:
1250: + ds - The `PetscDS` object
1251: . f - The field number
1252: - k - The highest derivative we need to tabulate
1254: Level: developer
1256: .seealso: `PetscDS`, `PetscDSGetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1257: @*/
1258: PetscErrorCode PetscDSSetJetDegree(PetscDS ds, PetscInt f, PetscInt k)
1259: {
1260: PetscFunctionBegin;
1262: 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);
1263: ds->jetDegree[f] = k;
1264: PetscFunctionReturn(PETSC_SUCCESS);
1265: }
1267: /*@
1268: PetscDSGetObjective - Get the pointwise objective function for a given test field that was provided with `PetscDSSetObjective()`
1270: Not Collective
1272: Input Parameters:
1273: + ds - The `PetscDS`
1274: - f - The test field number
1276: Output Parameter:
1277: . obj - integrand for the test function term, see `PetscPointFn`
1279: Level: intermediate
1281: Note:
1282: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi\,\mathrm{obj}(u, u_t, \nabla u, x, t)$
1284: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSSetObjective()`, `PetscDSGetResidual()`
1285: @*/
1286: PetscErrorCode PetscDSGetObjective(PetscDS ds, PetscInt f, PetscPointFn **obj)
1287: {
1288: PetscPointFn **tmp;
1289: PetscInt n;
1291: PetscFunctionBegin;
1293: PetscAssertPointer(obj, 3);
1294: 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);
1295: PetscCall(PetscWeakFormGetObjective(ds->wf, NULL, 0, f, 0, &n, &tmp));
1296: *obj = tmp ? tmp[0] : NULL;
1297: PetscFunctionReturn(PETSC_SUCCESS);
1298: }
1300: /*@
1301: PetscDSSetObjective - Set the pointwise objective function for a given test field
1303: Not Collective
1305: Input Parameters:
1306: + ds - The `PetscDS`
1307: . f - The test field number
1308: - obj - integrand for the test function term, see `PetscPointFn`
1310: Level: intermediate
1312: Note:
1313: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi\,\mathrm{obj}(u, u_t, \nabla u, x, t)$
1315: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSGetObjective()`, `PetscDSSetResidual()`
1316: @*/
1317: PetscErrorCode PetscDSSetObjective(PetscDS ds, PetscInt f, PetscPointFn *obj)
1318: {
1319: PetscFunctionBegin;
1322: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1323: PetscCall(PetscWeakFormSetIndexObjective(ds->wf, NULL, 0, f, 0, 0, obj));
1324: PetscFunctionReturn(PETSC_SUCCESS);
1325: }
1327: /*@
1328: PetscDSGetResidual - Get the pointwise residual function for a given test field
1330: Not Collective
1332: Input Parameters:
1333: + ds - The `PetscDS`
1334: - f - The test field number
1336: Output Parameters:
1337: + f0 - integrand for the test function term, see `PetscPointFn`
1338: - f1 - integrand for the test function gradient term, see `PetscPointFn`
1340: Level: intermediate
1342: Note:
1343: 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)$
1345: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSSetResidual()`
1346: @*/
1347: PetscErrorCode PetscDSGetResidual(PetscDS ds, PetscInt f, PetscPointFn **f0, PetscPointFn **f1)
1348: {
1349: PetscPointFn **tmp0, **tmp1;
1350: PetscInt n0, n1;
1352: PetscFunctionBegin;
1354: 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);
1355: PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1));
1356: *f0 = tmp0 ? tmp0[0] : NULL;
1357: *f1 = tmp1 ? tmp1[0] : NULL;
1358: PetscFunctionReturn(PETSC_SUCCESS);
1359: }
1361: /*@
1362: PetscDSSetResidual - Set the pointwise residual function for a given test field
1364: Not Collective
1366: Input Parameters:
1367: + ds - The `PetscDS`
1368: . f - The test field number
1369: . f0 - integrand for the test function term, see `PetscPointFn`
1370: - f1 - integrand for the test function gradient term, see `PetscPointFn`
1372: Level: intermediate
1374: Note:
1375: 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)$
1377: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSGetResidual()`
1378: @*/
1379: PetscErrorCode PetscDSSetResidual(PetscDS ds, PetscInt f, PetscPointFn *f0, PetscPointFn *f1)
1380: {
1381: PetscFunctionBegin;
1385: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1386: PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1));
1387: PetscFunctionReturn(PETSC_SUCCESS);
1388: }
1390: /*@
1391: PetscDSGetRHSResidual - Get the pointwise RHS residual function for explicit timestepping for a given test field
1393: Not Collective
1395: Input Parameters:
1396: + ds - The `PetscDS`
1397: - f - The test field number
1399: Output Parameters:
1400: + f0 - integrand for the test function term, see `PetscPointFn`
1401: - f1 - integrand for the test function gradient term, see `PetscPointFn`
1403: Level: intermediate
1405: Note:
1406: 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)$
1408: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSSetRHSResidual()`
1409: @*/
1410: PetscErrorCode PetscDSGetRHSResidual(PetscDS ds, PetscInt f, PetscPointFn **f0, PetscPointFn **f1)
1411: {
1412: PetscPointFn **tmp0, **tmp1;
1413: PetscInt n0, n1;
1415: PetscFunctionBegin;
1417: 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);
1418: PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 100, &n0, &tmp0, &n1, &tmp1));
1419: *f0 = tmp0 ? tmp0[0] : NULL;
1420: *f1 = tmp1 ? tmp1[0] : NULL;
1421: PetscFunctionReturn(PETSC_SUCCESS);
1422: }
1424: /*@
1425: PetscDSSetRHSResidual - Set the pointwise residual function for explicit timestepping for a given test field
1427: Not Collective
1429: Input Parameters:
1430: + ds - The `PetscDS`
1431: . f - The test field number
1432: . f0 - integrand for the test function term, see `PetscPointFn`
1433: - f1 - integrand for the test function gradient term, see `PetscPointFn`
1435: Level: intermediate
1437: Note:
1438: 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)$
1440: .seealso: `PetscDS`, `PetscDSGetResidual()`
1441: @*/
1442: PetscErrorCode PetscDSSetRHSResidual(PetscDS ds, PetscInt f, PetscPointFn *f0, PetscPointFn *f1)
1443: {
1444: PetscFunctionBegin;
1448: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1449: PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 100, 0, f0, 0, f1));
1450: PetscFunctionReturn(PETSC_SUCCESS);
1451: }
1453: /*@
1454: PetscDSHasJacobian - Checks that the Jacobian functions have been set
1456: Not Collective
1458: Input Parameter:
1459: . ds - The `PetscDS`
1461: Output Parameter:
1462: . hasJac - flag that indicates the pointwise function for the Jacobian has been set
1464: Level: intermediate
1466: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1467: @*/
1468: PetscErrorCode PetscDSHasJacobian(PetscDS ds, PetscBool *hasJac)
1469: {
1470: PetscFunctionBegin;
1472: PetscCall(PetscWeakFormHasJacobian(ds->wf, hasJac));
1473: PetscFunctionReturn(PETSC_SUCCESS);
1474: }
1476: /*@
1477: PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field
1479: Not Collective
1481: Input Parameters:
1482: + ds - The `PetscDS`
1483: . f - The test field number
1484: - g - The field number
1486: Output Parameters:
1487: + g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1488: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1489: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1490: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1492: Level: intermediate
1494: Note:
1495: We are using a first order FEM model for the weak form\:
1497: $$
1498: \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
1499: + \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
1500: $$
1502: .seealso: `PetscDS`, `PetscDSSetJacobian()`, `PetscPointJacFn`
1503: @*/
1504: PetscErrorCode PetscDSGetJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn **g0, PetscPointJacFn **g1, PetscPointJacFn **g2, PetscPointJacFn **g3)
1505: {
1506: PetscPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
1507: PetscInt n0, n1, n2, n3;
1509: PetscFunctionBegin;
1511: 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);
1512: 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);
1513: PetscCall(PetscWeakFormGetJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1514: *g0 = tmp0 ? tmp0[0] : NULL;
1515: *g1 = tmp1 ? tmp1[0] : NULL;
1516: *g2 = tmp2 ? tmp2[0] : NULL;
1517: *g3 = tmp3 ? tmp3[0] : NULL;
1518: PetscFunctionReturn(PETSC_SUCCESS);
1519: }
1521: /*@
1522: PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields
1524: Not Collective
1526: Input Parameters:
1527: + ds - The `PetscDS`
1528: . f - The test field number
1529: . g - The field number
1530: . g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1531: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1532: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1533: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1535: Level: intermediate
1537: Note:
1538: We are using a first order FEM model for the weak form\:
1540: $$
1541: \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
1542: + \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
1543: $$
1545: .seealso: `PetscDS`, `PetscDSGetJacobian()`, `PetscPointJacFn`
1546: @*/
1547: PetscErrorCode PetscDSSetJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn *g0, PetscPointJacFn *g1, PetscPointJacFn *g2, PetscPointJacFn *g3)
1548: {
1549: PetscFunctionBegin;
1555: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1556: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1557: PetscCall(PetscWeakFormSetIndexJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1558: PetscFunctionReturn(PETSC_SUCCESS);
1559: }
1561: /*@
1562: PetscDSUseJacobianPreconditioner - Set whether to construct a Jacobian preconditioner
1564: Not Collective
1566: Input Parameters:
1567: + prob - The `PetscDS`
1568: - useJacPre - flag that enables construction of a Jacobian preconditioner
1570: Level: intermediate
1572: Developer Note:
1573: Should be called `PetscDSSetUseJacobianPreconditioner()`
1575: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1576: @*/
1577: PetscErrorCode PetscDSUseJacobianPreconditioner(PetscDS prob, PetscBool useJacPre)
1578: {
1579: PetscFunctionBegin;
1581: prob->useJacPre = useJacPre;
1582: PetscFunctionReturn(PETSC_SUCCESS);
1583: }
1585: /*@
1586: PetscDSHasJacobianPreconditioner - Checks if a Jacobian matrix for constructing a preconditioner has been set
1588: Not Collective
1590: Input Parameter:
1591: . ds - The `PetscDS`
1593: Output Parameter:
1594: . hasJacPre - the flag
1596: Level: intermediate
1598: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1599: @*/
1600: PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS ds, PetscBool *hasJacPre)
1601: {
1602: PetscFunctionBegin;
1604: *hasJacPre = PETSC_FALSE;
1605: if (!ds->useJacPre) PetscFunctionReturn(PETSC_SUCCESS);
1606: PetscCall(PetscWeakFormHasJacobianPreconditioner(ds->wf, hasJacPre));
1607: PetscFunctionReturn(PETSC_SUCCESS);
1608: }
1610: /*@
1611: PetscDSGetJacobianPreconditioner - Get the pointwise Jacobian function for given test and basis field that constructs the matrix used
1612: to compute the preconditioner. If this is missing, the system matrix is used to build the preconditioner.
1614: Not Collective
1616: Input Parameters:
1617: + ds - The `PetscDS`
1618: . f - The test field number
1619: - g - The field number
1621: Output Parameters:
1622: + g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1623: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1624: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1625: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1627: Level: intermediate
1629: Note:
1630: We are using a first order FEM model for the weak form\:
1632: $$
1633: \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
1634: + \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
1635: $$
1637: Developer Note:
1638: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
1640: .seealso: `PetscDS`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`, `PetscPointJacFn`
1641: @*/
1642: PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn **g0, PetscPointJacFn **g1, PetscPointJacFn **g2, PetscPointJacFn **g3)
1643: {
1644: PetscPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
1645: PetscInt n0, n1, n2, n3;
1647: PetscFunctionBegin;
1649: 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);
1650: 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);
1651: PetscCall(PetscWeakFormGetJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1652: *g0 = tmp0 ? tmp0[0] : NULL;
1653: *g1 = tmp1 ? tmp1[0] : NULL;
1654: *g2 = tmp2 ? tmp2[0] : NULL;
1655: *g3 = tmp3 ? tmp3[0] : NULL;
1656: PetscFunctionReturn(PETSC_SUCCESS);
1657: }
1659: /*@
1660: PetscDSSetJacobianPreconditioner - Set the pointwise Jacobian function for given test and basis fields that constructs the matrix used
1661: to compute the preconditioner. If this is missing, the system matrix is used to build the preconditioner.
1663: Not Collective
1665: Input Parameters:
1666: + ds - The `PetscDS`
1667: . f - The test field number
1668: . g - The field number
1669: . g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1670: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1671: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1672: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1674: Level: intermediate
1676: Note:
1677: We are using a first order FEM model for the weak form\:
1679: $$
1680: \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
1681: + \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
1682: $$
1684: Developer Note:
1685: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
1687: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobian()`, `PetscPointJacFn`
1688: @*/
1689: PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn *g0, PetscPointJacFn *g1, PetscPointJacFn *g2, PetscPointJacFn *g3)
1690: {
1691: PetscFunctionBegin;
1697: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1698: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1699: PetscCall(PetscWeakFormSetIndexJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1700: PetscFunctionReturn(PETSC_SUCCESS);
1701: }
1703: /*@
1704: PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, $dF/du_t$, has been set
1706: Not Collective
1708: Input Parameter:
1709: . ds - The `PetscDS`
1711: Output Parameter:
1712: . hasDynJac - flag that pointwise function for dynamic Jacobian has been set
1714: Level: intermediate
1716: .seealso: `PetscDS`, `PetscDSGetDynamicJacobian()`, `PetscDSSetDynamicJacobian()`, `PetscDSGetJacobian()`
1717: @*/
1718: PetscErrorCode PetscDSHasDynamicJacobian(PetscDS ds, PetscBool *hasDynJac)
1719: {
1720: PetscFunctionBegin;
1722: PetscCall(PetscWeakFormHasDynamicJacobian(ds->wf, hasDynJac));
1723: PetscFunctionReturn(PETSC_SUCCESS);
1724: }
1726: /*@
1727: PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, $dF/du_t$, function for given test and basis field
1729: Not Collective
1731: Input Parameters:
1732: + ds - The `PetscDS`
1733: . f - The test field number
1734: - g - The field number
1736: Output Parameters:
1737: + g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1738: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1739: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1740: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1742: Level: intermediate
1744: Note:
1745: We are using a first order FEM model for the weak form\:
1747: $$
1748: \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
1749: + \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
1750: $$
1752: .seealso: `PetscDS`, `PetscDSSetJacobian()`, `PetscDSSetDynamicJacobian()`, `PetscPointJacFn`
1753: @*/
1754: PetscErrorCode PetscDSGetDynamicJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn **g0, PetscPointJacFn **g1, PetscPointJacFn **g2, PetscPointJacFn **g3)
1755: {
1756: PetscPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
1757: PetscInt n0, n1, n2, n3;
1759: PetscFunctionBegin;
1761: 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);
1762: 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);
1763: PetscCall(PetscWeakFormGetDynamicJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1764: *g0 = tmp0 ? tmp0[0] : NULL;
1765: *g1 = tmp1 ? tmp1[0] : NULL;
1766: *g2 = tmp2 ? tmp2[0] : NULL;
1767: *g3 = tmp3 ? tmp3[0] : NULL;
1768: PetscFunctionReturn(PETSC_SUCCESS);
1769: }
1771: /*@
1772: PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, $dF/du_t$, function for given test and basis fields
1774: Not Collective
1776: Input Parameters:
1777: + ds - The `PetscDS`
1778: . f - The test field number
1779: . g - The field number
1780: . g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1781: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1782: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1783: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1785: Level: intermediate
1787: Note:
1788: We are using a first order FEM model for the weak form\:
1790: $$
1791: \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
1792: + \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
1793: $$
1795: .seealso: `PetscDS`, `PetscDSGetDynamicJacobian()`, `PetscDSGetJacobian()`, `PetscPointJacFn`
1796: @*/
1797: PetscErrorCode PetscDSSetDynamicJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn *g0, PetscPointJacFn *g1, PetscPointJacFn *g2, PetscPointJacFn *g3)
1798: {
1799: PetscFunctionBegin;
1805: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1806: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1807: PetscCall(PetscWeakFormSetIndexDynamicJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1808: PetscFunctionReturn(PETSC_SUCCESS);
1809: }
1811: /*@
1812: PetscDSGetRiemannSolver - Returns the Riemann solver for the given field
1814: Not Collective
1816: Input Parameters:
1817: + ds - The `PetscDS` object
1818: - f - The field number
1820: Output Parameter:
1821: . r - Riemann solver, see `PetscRiemannFn`
1823: Level: intermediate
1825: .seealso: `PetscDS`, `PetscRiemannFn`, `PetscDSSetRiemannSolver()`
1826: @*/
1827: PetscErrorCode PetscDSGetRiemannSolver(PetscDS ds, PetscInt f, PetscRiemannFn **r)
1828: {
1829: PetscRiemannFn **tmp;
1830: PetscInt n;
1832: PetscFunctionBegin;
1834: PetscAssertPointer(r, 3);
1835: 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);
1836: PetscCall(PetscWeakFormGetRiemannSolver(ds->wf, NULL, 0, f, 0, &n, &tmp));
1837: *r = tmp ? tmp[0] : NULL;
1838: PetscFunctionReturn(PETSC_SUCCESS);
1839: }
1841: /*@
1842: PetscDSSetRiemannSolver - Sets the Riemann solver for the given field
1844: Not Collective
1846: Input Parameters:
1847: + ds - The `PetscDS` object
1848: . f - The field number
1849: - r - Riemann solver, see `PetscRiemannFn`
1851: Level: intermediate
1853: .seealso: `PetscDS`, `PetscRiemannFn`, `PetscDSGetRiemannSolver()`
1854: @*/
1855: PetscErrorCode PetscDSSetRiemannSolver(PetscDS ds, PetscInt f, PetscRiemannFn *r)
1856: {
1857: PetscFunctionBegin;
1860: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1861: PetscCall(PetscWeakFormSetIndexRiemannSolver(ds->wf, NULL, 0, f, 0, 0, r));
1862: PetscFunctionReturn(PETSC_SUCCESS);
1863: }
1865: /*@
1866: PetscDSGetUpdate - Get the pointwise update function for a given field
1868: Not Collective
1870: Input Parameters:
1871: + ds - The `PetscDS`
1872: - f - The field number
1874: Output Parameter:
1875: . update - update function, see `PetscPointFn`
1877: Level: intermediate
1879: .seealso: `PetscDS`, `PetscPointFn`, `PetscDSSetUpdate()`, `PetscDSSetResidual()`
1880: @*/
1881: PetscErrorCode PetscDSGetUpdate(PetscDS ds, PetscInt f, PetscPointFn **update)
1882: {
1883: PetscFunctionBegin;
1885: 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);
1886: if (update) {
1887: PetscAssertPointer(update, 3);
1888: *update = ds->update[f];
1889: }
1890: PetscFunctionReturn(PETSC_SUCCESS);
1891: }
1893: /*@
1894: PetscDSSetUpdate - Set the pointwise update function for a given field
1896: Not Collective
1898: Input Parameters:
1899: + ds - The `PetscDS`
1900: . f - The field number
1901: - update - update function, see `PetscPointFn`
1903: Level: intermediate
1905: .seealso: `PetscDS`, `PetscPointFn`, `PetscDSGetResidual()`
1906: @*/
1907: PetscErrorCode PetscDSSetUpdate(PetscDS ds, PetscInt f, PetscPointFn *update)
1908: {
1909: PetscFunctionBegin;
1912: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1913: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
1914: ds->update[f] = update;
1915: PetscFunctionReturn(PETSC_SUCCESS);
1916: }
1918: /*@
1919: PetscDSGetContext - Returns the context that was passed by `PetscDSSetContext()`
1921: Not Collective
1923: Input Parameters:
1924: + ds - The `PetscDS`
1925: . f - The field number
1926: - ctx - the context
1928: Level: intermediate
1930: Fortran Notes:
1931: This only works when the context is a Fortran derived type or a `PetscObject`. Define `ctx` with
1932: .vb
1933: type(tUsertype), pointer :: ctx
1934: .ve
1936: .seealso: `PetscDS`, `PetscPointFn`, `PetscDSSetContext()`
1937: @*/
1938: PetscErrorCode PetscDSGetContext(PetscDS ds, PetscInt f, PetscCtxRt ctx)
1939: {
1940: PetscFunctionBegin;
1942: 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);
1943: PetscAssertPointer(ctx, 3);
1944: *(void **)ctx = ds->ctx[f];
1945: PetscFunctionReturn(PETSC_SUCCESS);
1946: }
1948: /*@
1949: PetscDSSetContext - Sets the context that is passed back to some of the pointwise function callbacks used by this `PetscDS`
1951: Not Collective
1953: Input Parameters:
1954: + ds - The `PetscDS`
1955: . f - The field number
1956: - ctx - the context
1958: Level: intermediate
1960: .seealso: `PetscDS`, `PetscPointFn`, `PetscDSGetContext()`
1961: @*/
1962: PetscErrorCode PetscDSSetContext(PetscDS ds, PetscInt f, PetscCtx ctx)
1963: {
1964: PetscFunctionBegin;
1966: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1967: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
1968: ds->ctx[f] = ctx;
1969: PetscFunctionReturn(PETSC_SUCCESS);
1970: }
1972: /*@
1973: PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field
1975: Not Collective
1977: Input Parameters:
1978: + ds - The PetscDS
1979: - f - The test field number
1981: Output Parameters:
1982: + f0 - boundary integrand for the test function term, see `PetscBdPointFn`
1983: - f1 - boundary integrand for the test function gradient term, see `PetscBdPointFn`
1985: Level: intermediate
1987: Note:
1988: We are using a first order FEM model for the weak form\:
1990: $$
1991: \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
1992: $$
1994: .seealso: `PetscDS`, `PetscBdPointFn`, `PetscDSSetBdResidual()`
1995: @*/
1996: PetscErrorCode PetscDSGetBdResidual(PetscDS ds, PetscInt f, PetscBdPointFn **f0, PetscBdPointFn **f1)
1997: {
1998: PetscBdPointFn **tmp0, **tmp1;
1999: PetscInt n0, n1;
2001: PetscFunctionBegin;
2003: 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);
2004: PetscCall(PetscWeakFormGetBdResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1));
2005: *f0 = tmp0 ? tmp0[0] : NULL;
2006: *f1 = tmp1 ? tmp1[0] : NULL;
2007: PetscFunctionReturn(PETSC_SUCCESS);
2008: }
2010: /*@
2011: PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field
2013: Not Collective
2015: Input Parameters:
2016: + ds - The `PetscDS`
2017: . f - The test field number
2018: . f0 - boundary integrand for the test function term, see `PetscBdPointFn`
2019: - f1 - boundary integrand for the test function gradient term, see `PetscBdPointFn`
2021: Level: intermediate
2023: Note:
2024: We are using a first order FEM model for the weak form\:
2026: $$
2027: \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
2028: $$
2030: .seealso: `PetscDS`, `PetscBdPointFn`, `PetscDSGetBdResidual()`
2031: @*/
2032: PetscErrorCode PetscDSSetBdResidual(PetscDS ds, PetscInt f, PetscBdPointFn *f0, PetscBdPointFn *f1)
2033: {
2034: PetscFunctionBegin;
2036: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2037: PetscCall(PetscWeakFormSetIndexBdResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1));
2038: PetscFunctionReturn(PETSC_SUCCESS);
2039: }
2041: /*@
2042: PetscDSHasBdJacobian - Indicates that boundary Jacobian functions have been set
2044: Not Collective
2046: Input Parameter:
2047: . ds - The `PetscDS`
2049: Output Parameter:
2050: . hasBdJac - flag that pointwise function for the boundary Jacobian has been set
2052: Level: intermediate
2054: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2055: @*/
2056: PetscErrorCode PetscDSHasBdJacobian(PetscDS ds, PetscBool *hasBdJac)
2057: {
2058: PetscFunctionBegin;
2060: PetscAssertPointer(hasBdJac, 2);
2061: PetscCall(PetscWeakFormHasBdJacobian(ds->wf, hasBdJac));
2062: PetscFunctionReturn(PETSC_SUCCESS);
2063: }
2065: /*@
2066: PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field
2068: Not Collective
2070: Input Parameters:
2071: + ds - The `PetscDS`
2072: . f - The test field number
2073: - g - The field number
2075: Output Parameters:
2076: + g0 - integrand for the test and basis function term, see `PetscBdPointJacFn`
2077: . g1 - integrand for the test function and basis function gradient term, see `PetscBdPointJacFn`
2078: . g2 - integrand for the test function gradient and basis function term, see `PetscBdPointJacFn`
2079: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscBdPointJacFn`
2081: Level: intermediate
2083: Note:
2084: We are using a first order FEM model for the weak form\:
2086: $$
2087: \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
2088: + \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
2089: $$
2091: .seealso: `PetscDS`, `PetscBdPointJacFn`, `PetscDSSetBdJacobian()`
2092: @*/
2093: PetscErrorCode PetscDSGetBdJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscBdPointJacFn **g0, PetscBdPointJacFn **g1, PetscBdPointJacFn **g2, PetscBdPointJacFn **g3)
2094: {
2095: PetscBdPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
2096: PetscInt n0, n1, n2, n3;
2098: PetscFunctionBegin;
2100: 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);
2101: 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);
2102: PetscCall(PetscWeakFormGetBdJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2103: *g0 = tmp0 ? tmp0[0] : NULL;
2104: *g1 = tmp1 ? tmp1[0] : NULL;
2105: *g2 = tmp2 ? tmp2[0] : NULL;
2106: *g3 = tmp3 ? tmp3[0] : NULL;
2107: PetscFunctionReturn(PETSC_SUCCESS);
2108: }
2110: /*@
2111: PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field
2113: Not Collective
2115: Input Parameters:
2116: + ds - The PetscDS
2117: . f - The test field number
2118: . g - The field number
2119: . g0 - integrand for the test and basis function term, see `PetscBdPointJacFn`
2120: . g1 - integrand for the test function and basis function gradient term, see `PetscBdPointJacFn`
2121: . g2 - integrand for the test function gradient and basis function term, see `PetscBdPointJacFn`
2122: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscBdPointJacFn`
2124: Level: intermediate
2126: Note:
2127: We are using a first order FEM model for the weak form\:
2129: $$
2130: \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
2131: + \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
2132: $$
2134: .seealso: `PetscDS`, `PetscBdPointJacFn`, `PetscDSGetBdJacobian()`
2135: @*/
2136: PetscErrorCode PetscDSSetBdJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscBdPointJacFn *g0, PetscBdPointJacFn *g1, PetscBdPointJacFn *g2, PetscBdPointJacFn *g3)
2137: {
2138: PetscFunctionBegin;
2144: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2145: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2146: PetscCall(PetscWeakFormSetIndexBdJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2147: PetscFunctionReturn(PETSC_SUCCESS);
2148: }
2150: /*@
2151: PetscDSHasBdJacobianPreconditioner - Signals that boundary Jacobian preconditioner functions have been set with `PetscDSSetBdJacobianPreconditioner()`
2153: Not Collective
2155: Input Parameter:
2156: . ds - The `PetscDS`
2158: Output Parameter:
2159: . hasBdJacPre - flag that pointwise function for the boundary Jacobian matrix to construct the preconditioner has been set
2161: Level: intermediate
2163: Developer Note:
2164: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
2166: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2167: @*/
2168: PetscErrorCode PetscDSHasBdJacobianPreconditioner(PetscDS ds, PetscBool *hasBdJacPre)
2169: {
2170: PetscFunctionBegin;
2172: PetscAssertPointer(hasBdJacPre, 2);
2173: PetscCall(PetscWeakFormHasBdJacobianPreconditioner(ds->wf, hasBdJacPre));
2174: PetscFunctionReturn(PETSC_SUCCESS);
2175: }
2177: /*@
2178: PetscDSGetBdJacobianPreconditioner - Get the pointwise boundary Jacobian function for given test and basis field that constructs the
2179: matrix used to construct the preconditioner
2181: Not Collective; No Fortran Support
2183: Input Parameters:
2184: + ds - The `PetscDS`
2185: . f - The test field number
2186: - g - The field number
2188: Output Parameters:
2189: + g0 - integrand for the test and basis function term, see `PetscBdPointJacFn`
2190: . g1 - integrand for the test function and basis function gradient term, see `PetscBdPointJacFn`
2191: . g2 - integrand for the test function gradient and basis function term, see `PetscBdPointJacFn`
2192: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscBdPointJacFn`
2194: Level: intermediate
2196: Note:
2197: We are using a first order FEM model for the weak form\:
2199: $$
2200: \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
2201: + \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
2202: $$
2204: Developer Note:
2205: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
2207: .seealso: `PetscDS`, `PetscBdPointJacFn`, `PetscDSSetBdJacobianPreconditioner()`
2208: @*/
2209: PetscErrorCode PetscDSGetBdJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, PetscBdPointJacFn **g0, PetscBdPointJacFn **g1, PetscBdPointJacFn **g2, PetscBdPointJacFn **g3)
2210: {
2211: PetscBdPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
2212: PetscInt n0, n1, n2, n3;
2214: PetscFunctionBegin;
2216: 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);
2217: 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);
2218: PetscCall(PetscWeakFormGetBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2219: *g0 = tmp0 ? tmp0[0] : NULL;
2220: *g1 = tmp1 ? tmp1[0] : NULL;
2221: *g2 = tmp2 ? tmp2[0] : NULL;
2222: *g3 = tmp3 ? tmp3[0] : NULL;
2223: PetscFunctionReturn(PETSC_SUCCESS);
2224: }
2226: /*@
2227: PetscDSSetBdJacobianPreconditioner - Set the pointwise boundary Jacobian preconditioner function for given test and basis field that constructs the
2228: matrix used to construct the preconditioner
2230: Not Collective; No Fortran Support
2232: Input Parameters:
2233: + ds - The `PetscDS`
2234: . f - The test field number
2235: . g - The field number
2236: . g0 - integrand for the test and basis function term, see `PetscBdPointJacFn`
2237: . g1 - integrand for the test function and basis function gradient term, see `PetscBdPointJacFn`
2238: . g2 - integrand for the test function gradient and basis function term, see `PetscBdPointJacFn`
2239: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscBdPointJacFn`
2241: Level: intermediate
2243: Note:
2244: We are using a first order FEM model for the weak form\:
2246: $$
2247: \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
2248: + \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
2249: $$
2251: Developer Note:
2252: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
2254: .seealso: `PetscDS`, `PetscBdPointJacFn`, `PetscDSGetBdJacobianPreconditioner()`
2255: @*/
2256: PetscErrorCode PetscDSSetBdJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, PetscBdPointJacFn *g0, PetscBdPointJacFn *g1, PetscBdPointJacFn *g2, PetscBdPointJacFn *g3)
2257: {
2258: PetscFunctionBegin;
2264: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2265: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2266: PetscCall(PetscWeakFormSetIndexBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2267: PetscFunctionReturn(PETSC_SUCCESS);
2268: }
2270: /*@
2271: PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field
2273: Not Collective
2275: Input Parameters:
2276: + prob - The `PetscDS`
2277: - f - The test field number
2279: Output Parameters:
2280: + sol - exact solution function for the test field, see `PetscPointExactSolutionFn`
2281: - ctx - exact solution context
2283: Level: intermediate
2285: .seealso: `PetscDS`, `PetscPointExactSolutionFn`, `PetscDSSetExactSolution()`, `PetscDSGetExactSolutionTimeDerivative()`
2286: @*/
2287: PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscPointExactSolutionFn **sol, void **ctx)
2288: {
2289: PetscFunctionBegin;
2291: 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);
2292: if (sol) {
2293: PetscAssertPointer(sol, 3);
2294: *sol = prob->exactSol[f];
2295: }
2296: if (ctx) {
2297: PetscAssertPointer(ctx, 4);
2298: *ctx = prob->exactCtx[f];
2299: }
2300: PetscFunctionReturn(PETSC_SUCCESS);
2301: }
2303: /*@
2304: PetscDSSetExactSolution - Set the pointwise exact solution function for a given test field
2306: Not Collective
2308: Input Parameters:
2309: + prob - The `PetscDS`
2310: . f - The test field number
2311: . sol - solution function for the test fields, see `PetscPointExactSolutionFn`
2312: - ctx - solution context or `NULL`
2314: Level: intermediate
2316: .seealso: `PetscDS`, `PetscPointExactSolutionFn`, `PetscDSGetExactSolution()`
2317: @*/
2318: PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscPointExactSolutionFn *sol, PetscCtx ctx)
2319: {
2320: PetscFunctionBegin;
2322: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2323: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
2324: if (sol) {
2326: prob->exactSol[f] = sol;
2327: }
2328: if (ctx) {
2330: prob->exactCtx[f] = ctx;
2331: }
2332: PetscFunctionReturn(PETSC_SUCCESS);
2333: }
2335: /*@
2336: PetscDSGetExactSolutionTimeDerivative - Get the pointwise time derivative of the exact solution function for a given test field
2338: Not Collective
2340: Input Parameters:
2341: + prob - The `PetscDS`
2342: - f - The test field number
2344: Output Parameters:
2345: + sol - time derivative of the exact solution for the test field, see `PetscPointExactSolutionFn`
2346: - ctx - the exact solution context
2348: Level: intermediate
2350: .seealso: `PetscDS`, `PetscPointExactSolutionFn`, `PetscDSSetExactSolutionTimeDerivative()`, `PetscDSGetExactSolution()`
2351: @*/
2352: PetscErrorCode PetscDSGetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscPointExactSolutionFn **sol, void **ctx)
2353: {
2354: PetscFunctionBegin;
2356: 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);
2357: if (sol) {
2358: PetscAssertPointer(sol, 3);
2359: *sol = prob->exactSol_t[f];
2360: }
2361: if (ctx) {
2362: PetscAssertPointer(ctx, 4);
2363: *ctx = prob->exactCtx_t[f];
2364: }
2365: PetscFunctionReturn(PETSC_SUCCESS);
2366: }
2368: /*@
2369: PetscDSSetExactSolutionTimeDerivative - Set the pointwise time derivative of the exact solution function for a given test field
2371: Not Collective
2373: Input Parameters:
2374: + prob - The `PetscDS`
2375: . f - The test field number
2376: . sol - time derivative of the solution function for the test fields, see `PetscPointExactSolutionFn`
2377: - ctx - the solution context or `NULL`
2379: Level: intermediate
2381: .seealso: `PetscDS`, `PetscPointExactSolutionFn`, `PetscDSGetExactSolutionTimeDerivative()`, `PetscDSSetExactSolution()`
2382: @*/
2383: PetscErrorCode PetscDSSetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscPointExactSolutionFn *sol, PetscCtx ctx)
2384: {
2385: PetscFunctionBegin;
2387: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2388: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
2389: if (sol) {
2391: prob->exactSol_t[f] = sol;
2392: }
2393: if (ctx) {
2395: prob->exactCtx_t[f] = ctx;
2396: }
2397: PetscFunctionReturn(PETSC_SUCCESS);
2398: }
2400: /*@
2401: PetscDSGetLowerBound - Get the pointwise lower bound function for a given field
2403: Not Collective
2405: Input Parameters:
2406: + ds - The PetscDS
2407: - f - The field number
2409: Output Parameters:
2410: + lb - lower bound function for the field, see `PetscPointBoundFn`
2411: - ctx - lower bound context that was set with `PetscDSSetLowerBound()`
2413: Level: intermediate
2415: .seealso: `PetscDS`, `PetscPointBoundFn`, `PetscDSSetLowerBound()`, `PetscDSGetUpperBound()`, `PetscDSGetExactSolution()`
2416: @*/
2417: PetscErrorCode PetscDSGetLowerBound(PetscDS ds, PetscInt f, PetscPointBoundFn **lb, void **ctx)
2418: {
2419: PetscFunctionBegin;
2421: 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);
2422: if (lb) {
2423: PetscAssertPointer(lb, 3);
2424: *lb = ds->lowerBound[f];
2425: }
2426: if (ctx) {
2427: PetscAssertPointer(ctx, 4);
2428: *ctx = ds->lowerCtx[f];
2429: }
2430: PetscFunctionReturn(PETSC_SUCCESS);
2431: }
2433: /*@
2434: PetscDSSetLowerBound - Set the pointwise lower bound function for a given field
2436: Not Collective
2438: Input Parameters:
2439: + ds - The `PetscDS`
2440: . f - The field number
2441: . lb - lower bound function for the test fields, see `PetscPointBoundFn`
2442: - ctx - lower bound context or `NULL` which will be passed to `lb`
2444: Level: intermediate
2446: .seealso: `PetscDS`, `PetscPointBoundFn`, `PetscDSGetLowerBound()`, `PetscDSGetUpperBound()`, `PetscDSGetExactSolution()`
2447: @*/
2448: PetscErrorCode PetscDSSetLowerBound(PetscDS ds, PetscInt f, PetscPointBoundFn *lb, PetscCtx ctx)
2449: {
2450: PetscFunctionBegin;
2452: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2453: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
2454: if (lb) {
2456: ds->lowerBound[f] = lb;
2457: }
2458: if (ctx) {
2460: ds->lowerCtx[f] = ctx;
2461: }
2462: PetscFunctionReturn(PETSC_SUCCESS);
2463: }
2465: /*@
2466: PetscDSGetUpperBound - Get the pointwise upper bound function for a given field
2468: Not Collective
2470: Input Parameters:
2471: + ds - The `PetscDS`
2472: - f - The field number
2474: Output Parameters:
2475: + ub - upper bound function for the field, see `PetscPointBoundFn`
2476: - ctx - upper bound context that was set with `PetscDSSetUpperBound()`
2478: Level: intermediate
2480: .seealso: `PetscDS`, `PetscPointBoundFn`, `PetscDSSetUpperBound()`, `PetscDSGetLowerBound()`, `PetscDSGetExactSolution()`
2481: @*/
2482: PetscErrorCode PetscDSGetUpperBound(PetscDS ds, PetscInt f, PetscPointBoundFn **ub, void **ctx)
2483: {
2484: PetscFunctionBegin;
2486: 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);
2487: if (ub) {
2488: PetscAssertPointer(ub, 3);
2489: *ub = ds->upperBound[f];
2490: }
2491: if (ctx) {
2492: PetscAssertPointer(ctx, 4);
2493: *ctx = ds->upperCtx[f];
2494: }
2495: PetscFunctionReturn(PETSC_SUCCESS);
2496: }
2498: /*@
2499: PetscDSSetUpperBound - Set the pointwise upper bound function for a given field
2501: Not Collective
2503: Input Parameters:
2504: + ds - The `PetscDS`
2505: . f - The field number
2506: . ub - upper bound function for the test fields, see `PetscPointBoundFn`
2507: - ctx - context or `NULL` that will be passed to `ub`
2509: Level: intermediate
2511: .seealso: `PetscDS`, `PetscPointBoundFn`, `PetscDSGetUpperBound()`, `PetscDSGetLowerBound()`, `PetscDSGetExactSolution()`
2512: @*/
2513: PetscErrorCode PetscDSSetUpperBound(PetscDS ds, PetscInt f, PetscPointBoundFn *ub, PetscCtx ctx)
2514: {
2515: PetscFunctionBegin;
2517: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2518: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
2519: if (ub) {
2521: ds->upperBound[f] = ub;
2522: }
2523: if (ctx) {
2525: ds->upperCtx[f] = ctx;
2526: }
2527: PetscFunctionReturn(PETSC_SUCCESS);
2528: }
2530: /*@
2531: PetscDSGetConstants - Returns the array of constants passed to point functions from a `PetscDS` object
2533: Not Collective
2535: Input Parameter:
2536: . ds - The `PetscDS` object
2538: Output Parameters:
2539: + numConstants - The number of constants, or pass in `NULL` if not required
2540: - constants - The array of constants, `NULL` if there are none
2542: Level: intermediate
2544: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSCreate()`
2545: @*/
2546: PetscErrorCode PetscDSGetConstants(PetscDS ds, PeOp PetscInt *numConstants, PeOp const PetscScalar *constants[])
2547: {
2548: PetscFunctionBegin;
2550: if (numConstants) {
2551: PetscAssertPointer(numConstants, 2);
2552: *numConstants = ds->numConstants;
2553: }
2554: if (constants) {
2555: PetscAssertPointer(constants, 3);
2556: *constants = ds->constants;
2557: }
2558: PetscFunctionReturn(PETSC_SUCCESS);
2559: }
2561: /*@
2562: PetscDSSetConstants - Set the array of constants passed to point functions from a `PetscDS`
2564: Not Collective
2566: Input Parameters:
2567: + ds - The `PetscDS` object
2568: . numConstants - The number of constants
2569: - constants - The array of constants, `NULL` if there are none
2571: Level: intermediate
2573: .seealso: `PetscDS`, `PetscDSGetConstants()`, `PetscDSCreate()`
2574: @*/
2575: PetscErrorCode PetscDSSetConstants(PetscDS ds, PetscInt numConstants, PetscScalar constants[])
2576: {
2577: PetscFunctionBegin;
2579: if (numConstants != ds->numConstants) {
2580: PetscCall(PetscFree(ds->constants));
2581: ds->numConstants = numConstants;
2582: PetscCall(PetscMalloc1(ds->numConstants + ds->numFuncConstants, &ds->constants));
2583: }
2584: if (ds->numConstants) {
2585: PetscAssertPointer(constants, 3);
2586: PetscCall(PetscArraycpy(ds->constants, constants, ds->numConstants));
2587: }
2588: PetscFunctionReturn(PETSC_SUCCESS);
2589: }
2591: /*@
2592: PetscDSSetIntegrationParameters - Set the parameters for a particular integration
2594: Not Collective
2596: Input Parameters:
2597: + ds - The `PetscDS` object
2598: . fieldI - The test field for a given point function, or `PETSC_DETERMINE`
2599: - fieldJ - The basis field for a given point function, or `PETSC_DETERMINE`
2601: Level: intermediate
2603: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSGetConstants()`, `PetscDSCreate()`
2604: @*/
2605: PetscErrorCode PetscDSSetIntegrationParameters(PetscDS ds, PetscInt fieldI, PetscInt fieldJ)
2606: {
2607: PetscFunctionBegin;
2609: ds->constants[ds->numConstants] = fieldI;
2610: ds->constants[ds->numConstants + 1] = fieldJ;
2611: PetscFunctionReturn(PETSC_SUCCESS);
2612: }
2614: /*@
2615: PetscDSSetCellParameters - Set the parameters for a particular cell
2617: Not Collective
2619: Input Parameters:
2620: + ds - The `PetscDS` object
2621: - volume - The cell volume
2623: Level: intermediate
2625: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSGetConstants()`, `PetscDSCreate()`
2626: @*/
2627: PetscErrorCode PetscDSSetCellParameters(PetscDS ds, PetscReal volume)
2628: {
2629: PetscFunctionBegin;
2631: ds->constants[ds->numConstants + 2] = volume;
2632: PetscFunctionReturn(PETSC_SUCCESS);
2633: }
2635: /*@
2636: PetscDSGetFieldIndex - Returns the index of the given field
2638: Not Collective
2640: Input Parameters:
2641: + prob - The `PetscDS` object
2642: - disc - The discretization object
2644: Output Parameter:
2645: . f - The field number
2647: Level: beginner
2649: .seealso: `PetscDS`, `PetscGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2650: @*/
2651: PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f)
2652: {
2653: PetscInt g;
2655: PetscFunctionBegin;
2657: PetscAssertPointer(f, 3);
2658: *f = -1;
2659: for (g = 0; g < prob->Nf; ++g) {
2660: if (disc == prob->disc[g]) break;
2661: }
2662: PetscCheck(g != prob->Nf, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS.");
2663: *f = g;
2664: PetscFunctionReturn(PETSC_SUCCESS);
2665: }
2667: /*@
2668: PetscDSGetFieldSize - Returns the size of the given field in the full space basis
2670: Not Collective
2672: Input Parameters:
2673: + prob - The `PetscDS` object
2674: - f - The field number
2676: Output Parameter:
2677: . size - The size
2679: Level: beginner
2681: .seealso: `PetscDS`, `PetscDSGetFieldOffset()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2682: @*/
2683: PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size)
2684: {
2685: PetscFunctionBegin;
2687: PetscAssertPointer(size, 3);
2688: 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);
2689: PetscCall(PetscDSSetUp(prob));
2690: *size = prob->Nb[f];
2691: PetscFunctionReturn(PETSC_SUCCESS);
2692: }
2694: /*@
2695: PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis
2697: Not Collective
2699: Input Parameters:
2700: + prob - The `PetscDS` object
2701: - f - The field number
2703: Output Parameter:
2704: . off - The offset
2706: Level: beginner
2708: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2709: @*/
2710: PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off)
2711: {
2712: PetscInt size;
2714: PetscFunctionBegin;
2716: PetscAssertPointer(off, 3);
2717: 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);
2718: *off = 0;
2719: for (PetscInt g = 0; g < f; ++g) {
2720: PetscCall(PetscDSGetFieldSize(prob, g, &size));
2721: *off += size;
2722: }
2723: PetscFunctionReturn(PETSC_SUCCESS);
2724: }
2726: /*@
2727: PetscDSGetFieldOffsetCohesive - Returns the offset of the given field in the full space basis on a cohesive cell
2729: Not Collective
2731: Input Parameters:
2732: + ds - The `PetscDS` object
2733: - f - The field number
2735: Output Parameter:
2736: . off - The offset
2738: Level: beginner
2740: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2741: @*/
2742: PetscErrorCode PetscDSGetFieldOffsetCohesive(PetscDS ds, PetscInt f, PetscInt *off)
2743: {
2744: PetscInt size;
2746: PetscFunctionBegin;
2748: PetscAssertPointer(off, 3);
2749: 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);
2750: *off = 0;
2751: for (PetscInt g = 0; g < f; ++g) {
2752: PetscBool cohesive;
2754: PetscCall(PetscDSGetCohesive(ds, g, &cohesive));
2755: PetscCall(PetscDSGetFieldSize(ds, g, &size));
2756: *off += cohesive ? size : size * 2;
2757: }
2758: PetscFunctionReturn(PETSC_SUCCESS);
2759: }
2761: /*@
2762: PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point
2764: Not Collective
2766: Input Parameter:
2767: . prob - The `PetscDS` object
2769: Output Parameter:
2770: . dimensions - The number of dimensions
2772: Level: beginner
2774: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2775: @*/
2776: PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[])
2777: {
2778: PetscFunctionBegin;
2780: PetscCall(PetscDSSetUp(prob));
2781: PetscAssertPointer(dimensions, 2);
2782: *dimensions = prob->Nb;
2783: PetscFunctionReturn(PETSC_SUCCESS);
2784: }
2786: /*@
2787: PetscDSGetComponents - Returns the number of components for each field on an evaluation point
2789: Not Collective
2791: Input Parameter:
2792: . prob - The `PetscDS` object
2794: Output Parameter:
2795: . components - The number of components
2797: Level: beginner
2799: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2800: @*/
2801: PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[])
2802: {
2803: PetscFunctionBegin;
2805: PetscCall(PetscDSSetUp(prob));
2806: PetscAssertPointer(components, 2);
2807: *components = prob->Nc;
2808: PetscFunctionReturn(PETSC_SUCCESS);
2809: }
2811: /*@
2812: PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point
2814: Not Collective
2816: Input Parameters:
2817: + prob - The `PetscDS` object
2818: - f - The field number
2820: Output Parameter:
2821: . off - The offset
2823: Level: beginner
2825: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2826: @*/
2827: PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off)
2828: {
2829: PetscFunctionBegin;
2831: PetscAssertPointer(off, 3);
2832: 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);
2833: PetscCall(PetscDSSetUp(prob));
2834: *off = prob->off[f];
2835: PetscFunctionReturn(PETSC_SUCCESS);
2836: }
2838: /*@
2839: PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point
2841: Not Collective
2843: Input Parameter:
2844: . prob - The `PetscDS` object
2846: Output Parameter:
2847: . offsets - The offsets
2849: Level: beginner
2851: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2852: @*/
2853: PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[])
2854: {
2855: PetscFunctionBegin;
2857: PetscAssertPointer(offsets, 2);
2858: PetscCall(PetscDSSetUp(prob));
2859: *offsets = prob->off;
2860: PetscFunctionReturn(PETSC_SUCCESS);
2861: }
2863: /*@
2864: PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point
2866: Not Collective
2868: Input Parameter:
2869: . prob - The `PetscDS` object
2871: Output Parameter:
2872: . offsets - The offsets
2874: Level: beginner
2876: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2877: @*/
2878: PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[])
2879: {
2880: PetscFunctionBegin;
2882: PetscAssertPointer(offsets, 2);
2883: PetscCall(PetscDSSetUp(prob));
2884: *offsets = prob->offDer;
2885: PetscFunctionReturn(PETSC_SUCCESS);
2886: }
2888: /*@
2889: PetscDSGetComponentOffsetsCohesive - Returns the offset of each field on an evaluation point
2891: Not Collective
2893: Input Parameters:
2894: + ds - The `PetscDS` object
2895: - s - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive
2897: Output Parameter:
2898: . offsets - The offsets
2900: Level: beginner
2902: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2903: @*/
2904: PetscErrorCode PetscDSGetComponentOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
2905: {
2906: PetscFunctionBegin;
2908: PetscAssertPointer(offsets, 3);
2909: PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS");
2910: PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s);
2911: PetscCall(PetscDSSetUp(ds));
2912: *offsets = ds->offCohesive[s];
2913: PetscFunctionReturn(PETSC_SUCCESS);
2914: }
2916: /*@
2917: PetscDSGetComponentDerivativeOffsetsCohesive - Returns the offset of each field derivative on an evaluation point
2919: Not Collective
2921: Input Parameters:
2922: + ds - The `PetscDS` object
2923: - s - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive
2925: Output Parameter:
2926: . offsets - The offsets
2928: Level: beginner
2930: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2931: @*/
2932: PetscErrorCode PetscDSGetComponentDerivativeOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
2933: {
2934: PetscFunctionBegin;
2936: PetscAssertPointer(offsets, 3);
2937: PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS");
2938: PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s);
2939: PetscCall(PetscDSSetUp(ds));
2940: *offsets = ds->offDerCohesive[s];
2941: PetscFunctionReturn(PETSC_SUCCESS);
2942: }
2944: /*@
2945: PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization
2947: Not Collective
2949: Input Parameter:
2950: . prob - The `PetscDS` object
2952: Output Parameter:
2953: . T - The basis function and derivatives tabulation at quadrature points for each field, see `PetscTabulation` for its details
2955: Level: intermediate
2957: Note:
2958: The tabulation is only valid so long as the `PetscDS` has not be destroyed. There is no `PetscDSRestoreTabulation()` in C.
2960: Fortran Note:
2961: Use the declaration
2962: .vb
2963: PetscTabulation, pointer :: tab(:)
2964: .ve
2965: and access the values using, for example,
2966: .vb
2967: tab(i)%ptr%K
2968: tab(i)%ptr%T(j)%ptr
2969: .ve
2970: where $ i = 1, 2, ..., Nf $ and $ j = 1, 2, ..., tab(i)%ptr%K+1 $.
2972: Use `PetscDSRestoreTabulation()` to restore the array
2974: Developer Note:
2975: The Fortran language syntax does not directly support arrays of pointers, the '%ptr' notation allows mimicking their use in Fortran.
2977: .seealso: `PetscDS`, `PetscTabulation`, `PetscDSCreate()`
2978: @*/
2979: PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscTabulation *T[]) PeNS
2980: {
2981: PetscFunctionBegin;
2983: PetscAssertPointer(T, 2);
2984: PetscCall(PetscDSSetUp(prob));
2985: *T = prob->T;
2986: PetscFunctionReturn(PETSC_SUCCESS);
2987: }
2989: /*@
2990: PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces
2992: Not Collective
2994: Input Parameter:
2995: . prob - The `PetscDS` object
2997: Output Parameter:
2998: . Tf - The basis function and derivative tabulation on each local face at quadrature points for each field
3000: Level: intermediate
3002: Note:
3003: The tabulation is only valid so long as the `PetscDS` has not be destroyed. There is no `PetscDSRestoreFaceTabulation()` in C.
3005: .seealso: `PetscTabulation`, `PetscDS`, `PetscDSGetTabulation()`, `PetscDSCreate()`
3006: @*/
3007: PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscTabulation *Tf[])
3008: {
3009: PetscFunctionBegin;
3011: PetscAssertPointer(Tf, 2);
3012: PetscCall(PetscDSSetUp(prob));
3013: *Tf = prob->Tf;
3014: PetscFunctionReturn(PETSC_SUCCESS);
3015: }
3017: /*@
3018: PetscDSGetEvaluationArrays - Get scratch arrays used to evaluate fields, time derivatives, and field gradients at quadrature points.
3020: Not Collective
3022: Input Parameter:
3023: . prob - the `PetscDS`
3025: Output Parameters:
3026: + u - array for the field values, or `NULL` if not needed
3027: . u_t - array for the field time derivatives, or `NULL` if not needed
3028: - u_x - array for the field gradients, or `NULL` if not needed
3030: Level: developer
3032: Note:
3033: The returned arrays are owned by the `PetscDS` and must not be freed by the caller.
3035: .seealso: `PetscDS`, `PetscDSGetWeakFormArrays()`, `PetscDSGetWorkspace()`
3036: @*/
3037: PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar *u[], PetscScalar *u_t[], PetscScalar *u_x[])
3038: {
3039: PetscFunctionBegin;
3041: PetscCall(PetscDSSetUp(prob));
3042: if (u) {
3043: PetscAssertPointer(u, 2);
3044: *u = prob->u;
3045: }
3046: if (u_t) {
3047: PetscAssertPointer(u_t, 3);
3048: *u_t = prob->u_t;
3049: }
3050: if (u_x) {
3051: PetscAssertPointer(u_x, 4);
3052: *u_x = prob->u_x;
3053: }
3054: PetscFunctionReturn(PETSC_SUCCESS);
3055: }
3057: PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar *f0[], PetscScalar *f1[], PetscScalar *g0[], PetscScalar *g1[], PetscScalar *g2[], PetscScalar *g3[])
3058: {
3059: PetscFunctionBegin;
3061: PetscCall(PetscDSSetUp(prob));
3062: if (f0) {
3063: PetscAssertPointer(f0, 2);
3064: *f0 = prob->f0;
3065: }
3066: if (f1) {
3067: PetscAssertPointer(f1, 3);
3068: *f1 = prob->f1;
3069: }
3070: if (g0) {
3071: PetscAssertPointer(g0, 4);
3072: *g0 = prob->g0;
3073: }
3074: if (g1) {
3075: PetscAssertPointer(g1, 5);
3076: *g1 = prob->g1;
3077: }
3078: if (g2) {
3079: PetscAssertPointer(g2, 6);
3080: *g2 = prob->g2;
3081: }
3082: if (g3) {
3083: PetscAssertPointer(g3, 7);
3084: *g3 = prob->g3;
3085: }
3086: PetscFunctionReturn(PETSC_SUCCESS);
3087: }
3089: /*@
3090: PetscDSGetWorkspace - Get scratch storage used during discretization computations.
3092: Not Collective
3094: Input Parameter:
3095: . prob - the `PetscDS`
3097: Output Parameters:
3098: + x - array for real-valued quadrature point coordinates, or `NULL` if not needed
3099: . basisReal - array for the real-valued basis function values, or `NULL` if not needed
3100: . basisDerReal - array for the real-valued basis function derivatives, or `NULL` if not needed
3101: . testReal - array for the real-valued test function values, or `NULL` if not needed
3102: - testDerReal - array for the real-valued test function derivatives, or `NULL` if not needed
3104: Level: developer
3106: Note:
3107: The returned arrays are owned by the `PetscDS` and must not be freed by the caller.
3109: .seealso: `PetscDS`, `PetscDSGetEvaluationArrays()`, `PetscDSGetWeakFormArrays()`
3110: @*/
3111: PetscErrorCode PetscDSGetWorkspace(PetscDS prob, PetscReal **x, PetscScalar **basisReal, PetscScalar **basisDerReal, PetscScalar **testReal, PetscScalar **testDerReal)
3112: {
3113: PetscFunctionBegin;
3115: PetscCall(PetscDSSetUp(prob));
3116: if (x) {
3117: PetscAssertPointer(x, 2);
3118: *x = prob->x;
3119: }
3120: if (basisReal) {
3121: PetscAssertPointer(basisReal, 3);
3122: *basisReal = prob->basisReal;
3123: }
3124: if (basisDerReal) {
3125: PetscAssertPointer(basisDerReal, 4);
3126: *basisDerReal = prob->basisDerReal;
3127: }
3128: if (testReal) {
3129: PetscAssertPointer(testReal, 5);
3130: *testReal = prob->testReal;
3131: }
3132: if (testDerReal) {
3133: PetscAssertPointer(testDerReal, 6);
3134: *testDerReal = prob->testDerReal;
3135: }
3136: PetscFunctionReturn(PETSC_SUCCESS);
3137: }
3139: /*@
3140: PetscDSAddBoundary - Add a boundary condition to the model.
3142: Collective
3144: Input Parameters:
3145: + ds - The `PetscDS` object
3146: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3147: . name - The name for the boundary condition
3148: . label - The label defining constrained points
3149: . Nv - The number of `DMLabel` values for constrained points
3150: . values - An array of label values for constrained points
3151: . field - The field to constrain
3152: . Nc - The number of constrained field components (0 will constrain all fields)
3153: . comps - An array of constrained component numbers
3154: . bcFunc - A pointwise function giving boundary values
3155: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or `NULL`
3156: - ctx - An optional application context for `bcFunc`
3158: Output Parameter:
3159: . bd - The boundary number
3161: Options Database Keys:
3162: + -bc_NAME values - comma separated list of values for the boundary condition NAME
3163: - -bc_NAME_comp comps - comma separated list of components for the boundary condition NAME
3165: Level: developer
3167: Note:
3168: Both `bcFunc` and `bcFunc_t` will depend on the boundary condition type. If the type if `DM_BC_ESSENTIAL`, then the calling sequence is\:
3169: .vb
3170: void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])
3171: .ve
3173: If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value, then the calling sequence is\:
3174: .vb
3175: void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3176: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3177: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3178: PetscReal time, const PetscReal x[], PetscScalar bcval[])
3179: .ve
3180: + dim - the coordinate dimension
3181: . Nf - the number of fields
3182: . uOff - the offset into u[] and u_t[] for each field
3183: . uOff_x - the offset into u_x[] for each field
3184: . u - each field evaluated at the current point
3185: . u_t - the time derivative of each field evaluated at the current point
3186: . u_x - the gradient of each field evaluated at the current point
3187: . aOff - the offset into a[] and a_t[] for each auxiliary field
3188: . aOff_x - the offset into a_x[] for each auxiliary field
3189: . a - each auxiliary field evaluated at the current point
3190: . a_t - the time derivative of each auxiliary field evaluated at the current point
3191: . a_x - the gradient of auxiliary each field evaluated at the current point
3192: . t - current time
3193: . x - coordinates of the current point
3194: . numConstants - number of constant parameters
3195: . constants - constant parameters
3196: - bcval - output values at the current point
3198: Notes:
3199: The pointwise functions are used to provide boundary values for essential boundary
3200: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3201: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3202: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3204: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundaryByName()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3205: @*/
3206: PetscErrorCode PetscDSAddBoundary(PetscDS ds, DMBoundaryConditionType type, const char name[], DMLabel label, PetscInt Nv, const PetscInt values[], PetscInt field, PetscInt Nc, const PetscInt comps[], PetscVoidFn *bcFunc, PetscVoidFn *bcFunc_t, PetscCtx ctx, PetscInt *bd)
3207: {
3208: DSBoundary head = ds->boundary, b;
3209: PetscInt n = 0;
3210: const char *lname;
3212: PetscFunctionBegin;
3215: PetscAssertPointer(name, 3);
3220: 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);
3221: if (Nc > 0) {
3222: PetscInt *fcomps;
3224: PetscCall(PetscDSGetComponents(ds, &fcomps));
3225: 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);
3226: for (PetscInt c = 0; c < Nc; ++c) {
3227: 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);
3228: }
3229: }
3230: PetscCall(PetscNew(&b));
3231: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3232: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf));
3233: PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf));
3234: PetscCall(PetscMalloc1(Nv, &b->values));
3235: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3236: PetscCall(PetscMalloc1(Nc, &b->comps));
3237: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3238: PetscCall(PetscObjectGetName((PetscObject)label, &lname));
3239: PetscCall(PetscStrallocpy(lname, (char **)&b->lname));
3240: b->type = type;
3241: b->label = label;
3242: b->Nv = Nv;
3243: b->field = field;
3244: b->Nc = Nc;
3245: b->func = bcFunc;
3246: b->func_t = bcFunc_t;
3247: b->ctx = ctx;
3248: b->next = NULL;
3249: /* Append to linked list so that we can preserve the order */
3250: if (!head) ds->boundary = b;
3251: while (head) {
3252: if (!head->next) {
3253: head->next = b;
3254: head = b;
3255: }
3256: head = head->next;
3257: ++n;
3258: }
3259: if (bd) {
3260: PetscAssertPointer(bd, 13);
3261: *bd = n;
3262: }
3263: PetscFunctionReturn(PETSC_SUCCESS);
3264: }
3266: // PetscClangLinter pragma ignore: -fdoc-section-header-unknown
3267: /*@
3268: PetscDSAddBoundaryByName - Add a boundary condition to the model.
3270: Collective
3272: Input Parameters:
3273: + ds - The `PetscDS` object
3274: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3275: . name - The boundary condition name
3276: . lname - The name of the label defining constrained points
3277: . Nv - The number of `DMLabel` values for constrained points
3278: . values - An array of label values for constrained points
3279: . field - The field to constrain
3280: . Nc - The number of constrained field components (0 will constrain all fields)
3281: . comps - An array of constrained component numbers
3282: . bcFunc - A pointwise function giving boundary values
3283: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or `NULL`
3284: - ctx - An optional application context for `bcFunc`
3286: Output Parameter:
3287: . bd - The boundary number
3289: Options Database Keys:
3290: + -bc_NAME values - comma separated list of values for the boundary condition NAME
3291: - -bc_NAME_comp comps - comma separated list of components for the boundary condition NAME
3293: Calling Sequence of `bcFunc` and `bcFunc_t`:
3294: If the type is `DM_BC_ESSENTIAL`
3295: .vb
3296: void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])
3297: .ve
3298: If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value,
3299: .vb
3300: void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3301: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3302: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3303: PetscReal time, const PetscReal x[], PetscScalar bcval[])
3304: .ve
3305: + dim - the coordinate dimension
3306: . Nf - the number of fields
3307: . uOff - the offset into `u`[] and `u_t`[] for each field
3308: . uOff_x - the offset into `u_x`[] for each field
3309: . u - each field evaluated at the current point
3310: . u_t - the time derivative of each field evaluated at the current point
3311: . u_x - the gradient of each field evaluated at the current point
3312: . aOff - the offset into `a`[] and `a_t`[] for each auxiliary field
3313: . aOff_x - the offset into `a_x`[] for each auxiliary field
3314: . a - each auxiliary field evaluated at the current point
3315: . a_t - the time derivative of each auxiliary field evaluated at the current point
3316: . a_x - the gradient of auxiliary each field evaluated at the current point
3317: . t - current time
3318: . x - coordinates of the current point
3319: . numConstants - number of constant parameters
3320: . constants - constant parameters
3321: - bcval - output values at the current point
3323: Level: developer
3325: Notes:
3326: The pointwise functions are used to provide boundary values for essential boundary
3327: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3328: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3329: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3331: This function should only be used with `DMFOREST` currently, since labels cannot be defined before the underlying `DMPLEX` is built.
3333: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3334: @*/
3335: PetscErrorCode PetscDSAddBoundaryByName(PetscDS ds, DMBoundaryConditionType type, const char name[], const char lname[], PetscInt Nv, const PetscInt values[], PetscInt field, PetscInt Nc, const PetscInt comps[], PetscVoidFn *bcFunc, PetscVoidFn *bcFunc_t, PetscCtx ctx, PetscInt *bd)
3336: {
3337: DSBoundary head = ds->boundary, b;
3338: PetscInt n = 0;
3340: PetscFunctionBegin;
3343: PetscAssertPointer(name, 3);
3344: PetscAssertPointer(lname, 4);
3348: PetscCall(PetscNew(&b));
3349: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3350: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf));
3351: PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf));
3352: PetscCall(PetscMalloc1(Nv, &b->values));
3353: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3354: PetscCall(PetscMalloc1(Nc, &b->comps));
3355: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3356: PetscCall(PetscStrallocpy(lname, (char **)&b->lname));
3357: b->type = type;
3358: b->label = NULL;
3359: b->Nv = Nv;
3360: b->field = field;
3361: b->Nc = Nc;
3362: b->func = bcFunc;
3363: b->func_t = bcFunc_t;
3364: b->ctx = ctx;
3365: b->next = NULL;
3366: /* Append to linked list so that we can preserve the order */
3367: if (!head) ds->boundary = b;
3368: while (head) {
3369: if (!head->next) {
3370: head->next = b;
3371: head = b;
3372: }
3373: head = head->next;
3374: ++n;
3375: }
3376: if (bd) {
3377: PetscAssertPointer(bd, 13);
3378: *bd = n;
3379: }
3380: PetscFunctionReturn(PETSC_SUCCESS);
3381: }
3383: /*@
3384: PetscDSUpdateBoundary - Change a boundary condition for the model.
3386: Input Parameters:
3387: + ds - The `PetscDS` object
3388: . bd - The boundary condition number
3389: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3390: . name - The boundary condition name
3391: . label - The label defining constrained points
3392: . Nv - The number of `DMLabel` ids for constrained points
3393: . values - An array of ids for constrained points
3394: . field - The field to constrain
3395: . Nc - The number of constrained field components
3396: . comps - An array of constrained component numbers
3397: . bcFunc - A pointwise function giving boundary values
3398: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or `NULL`
3399: - ctx - An optional application context for `bcFunc`
3401: Level: developer
3403: Notes:
3404: The pointwise functions are used to provide boundary values for essential boundary
3405: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3406: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3407: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3409: The boundary condition number is the order in which it was registered. The user can get the number of boundary conditions from `PetscDSGetNumBoundary()`.
3410: See `PetscDSAddBoundary()` for a description of the calling sequences for the callbacks.
3412: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSGetNumBoundary()`, `DMLabel`
3413: @*/
3414: PetscErrorCode PetscDSUpdateBoundary(PetscDS ds, PetscInt bd, DMBoundaryConditionType type, const char name[], DMLabel label, PetscInt Nv, const PetscInt values[], PetscInt field, PetscInt Nc, const PetscInt comps[], PetscVoidFn *bcFunc, PetscVoidFn *bcFunc_t, PetscCtx ctx)
3415: {
3416: DSBoundary b = ds->boundary;
3417: PetscInt n = 0;
3419: PetscFunctionBegin;
3421: while (b) {
3422: if (n == bd) break;
3423: b = b->next;
3424: ++n;
3425: }
3426: PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n);
3427: if (name) {
3428: PetscCall(PetscFree(b->name));
3429: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3430: }
3431: b->type = type;
3432: if (label) {
3433: const char *name;
3435: b->label = label;
3436: PetscCall(PetscFree(b->lname));
3437: PetscCall(PetscObjectGetName((PetscObject)label, &name));
3438: PetscCall(PetscStrallocpy(name, (char **)&b->lname));
3439: }
3440: if (Nv >= 0) {
3441: b->Nv = Nv;
3442: PetscCall(PetscFree(b->values));
3443: PetscCall(PetscMalloc1(Nv, &b->values));
3444: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3445: }
3446: if (field >= 0) b->field = field;
3447: if (Nc >= 0) {
3448: b->Nc = Nc;
3449: PetscCall(PetscFree(b->comps));
3450: PetscCall(PetscMalloc1(Nc, &b->comps));
3451: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3452: }
3453: if (bcFunc) b->func = bcFunc;
3454: if (bcFunc_t) b->func_t = bcFunc_t;
3455: if (ctx) b->ctx = ctx;
3456: PetscFunctionReturn(PETSC_SUCCESS);
3457: }
3459: /*@
3460: PetscDSGetNumBoundary - Get the number of registered boundary conditions
3462: Input Parameter:
3463: . ds - The `PetscDS` object
3465: Output Parameter:
3466: . numBd - The number of boundary conditions
3468: Level: intermediate
3470: .seealso: `PetscDS`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`
3471: @*/
3472: PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd)
3473: {
3474: DSBoundary b = ds->boundary;
3476: PetscFunctionBegin;
3478: PetscAssertPointer(numBd, 2);
3479: *numBd = 0;
3480: while (b) {
3481: ++(*numBd);
3482: b = b->next;
3483: }
3484: PetscFunctionReturn(PETSC_SUCCESS);
3485: }
3487: /*@
3488: PetscDSGetBoundary - Gets a boundary condition from the model
3490: Input Parameters:
3491: + ds - The `PetscDS` object
3492: - bd - The boundary condition number
3494: Output Parameters:
3495: + wf - The `PetscWeakForm` holding the pointwise functions
3496: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3497: . name - The boundary condition name
3498: . label - The label defining constrained points
3499: . Nv - The number of `DMLabel` ids for constrained points
3500: . values - An array of ids for constrained points
3501: . field - The field to constrain
3502: . Nc - The number of constrained field components
3503: . comps - An array of constrained component numbers
3504: . func - A pointwise function giving boundary values
3505: . func_t - A pointwise function giving the time derivative of the boundary values
3506: - ctx - An optional application context for `func`
3508: Level: developer
3510: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `DMLabel`
3511: @*/
3512: PetscErrorCode PetscDSGetBoundary(PetscDS ds, PetscInt bd, PetscWeakForm *wf, DMBoundaryConditionType *type, const char *name[], DMLabel *label, PetscInt *Nv, const PetscInt *values[], PetscInt *field, PetscInt *Nc, const PetscInt *comps[], PetscVoidFn **func, PetscVoidFn **func_t, PetscCtxRt ctx)
3513: {
3514: DSBoundary b = ds->boundary;
3515: PetscInt n = 0;
3517: PetscFunctionBegin;
3519: while (b) {
3520: if (n == bd) break;
3521: b = b->next;
3522: ++n;
3523: }
3524: PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n);
3525: if (wf) {
3526: PetscAssertPointer(wf, 3);
3527: *wf = b->wf;
3528: }
3529: if (type) {
3530: PetscAssertPointer(type, 4);
3531: *type = b->type;
3532: }
3533: if (name) {
3534: PetscAssertPointer(name, 5);
3535: *name = b->name;
3536: }
3537: if (label) {
3538: PetscAssertPointer(label, 6);
3539: *label = b->label;
3540: }
3541: if (Nv) {
3542: PetscAssertPointer(Nv, 7);
3543: *Nv = b->Nv;
3544: }
3545: if (values) {
3546: PetscAssertPointer(values, 8);
3547: *values = b->values;
3548: }
3549: if (field) {
3550: PetscAssertPointer(field, 9);
3551: *field = b->field;
3552: }
3553: if (Nc) {
3554: PetscAssertPointer(Nc, 10);
3555: *Nc = b->Nc;
3556: }
3557: if (comps) {
3558: PetscAssertPointer(comps, 11);
3559: *comps = b->comps;
3560: }
3561: if (func) {
3562: PetscAssertPointer(func, 12);
3563: *func = b->func;
3564: }
3565: if (func_t) {
3566: PetscAssertPointer(func_t, 13);
3567: *func_t = b->func_t;
3568: }
3569: if (ctx) {
3570: PetscAssertPointer(ctx, 14);
3571: *(void **)ctx = b->ctx;
3572: }
3573: PetscFunctionReturn(PETSC_SUCCESS);
3574: }
3576: /*@
3577: PetscDSUpdateBoundaryLabels - Update `DMLabel` in each boundary condition using the label name and the input `DM`
3579: Not Collective
3581: Input Parameters:
3582: + ds - The source `PetscDS` object
3583: - dm - The `DM` holding labels
3585: Level: intermediate
3587: .seealso: `PetscDS`, `DMBoundary`, `DM`, `PetscDSCopyBoundary()`, `PetscDSCreate()`, `DMGetLabel()`
3588: @*/
3589: PetscErrorCode PetscDSUpdateBoundaryLabels(PetscDS ds, DM dm)
3590: {
3591: DSBoundary b;
3593: PetscFunctionBegin;
3596: for (b = ds->boundary; b; b = b->next) {
3597: if (b->lname) PetscCall(DMGetLabel(dm, b->lname, &b->label));
3598: }
3599: PetscFunctionReturn(PETSC_SUCCESS);
3600: }
3602: static PetscErrorCode DSBoundaryDuplicate_Internal(DSBoundary b, DSBoundary *bNew)
3603: {
3604: PetscFunctionBegin;
3605: PetscCall(PetscNew(bNew));
3606: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &(*bNew)->wf));
3607: PetscCall(PetscWeakFormCopy(b->wf, (*bNew)->wf));
3608: PetscCall(PetscStrallocpy(b->name, (char **)&((*bNew)->name)));
3609: PetscCall(PetscStrallocpy(b->lname, (char **)&((*bNew)->lname)));
3610: (*bNew)->type = b->type;
3611: (*bNew)->label = b->label;
3612: (*bNew)->Nv = b->Nv;
3613: PetscCall(PetscMalloc1(b->Nv, &(*bNew)->values));
3614: PetscCall(PetscArraycpy((*bNew)->values, b->values, b->Nv));
3615: (*bNew)->field = b->field;
3616: (*bNew)->Nc = b->Nc;
3617: PetscCall(PetscMalloc1(b->Nc, &(*bNew)->comps));
3618: PetscCall(PetscArraycpy((*bNew)->comps, b->comps, b->Nc));
3619: (*bNew)->func = b->func;
3620: (*bNew)->func_t = b->func_t;
3621: (*bNew)->ctx = b->ctx;
3622: PetscFunctionReturn(PETSC_SUCCESS);
3623: }
3625: /*@
3626: PetscDSCopyBoundary - Copy all boundary condition objects to the new `PetscDS`
3628: Not Collective
3630: Input Parameters:
3631: + ds - The source `PetscDS` object
3632: . numFields - The number of selected fields, or `PETSC_DEFAULT` for all fields
3633: - fields - The selected fields, or `NULL` for all fields
3635: Output Parameter:
3636: . newds - The target `PetscDS`, now with a copy of the boundary conditions
3638: Level: intermediate
3640: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3641: @*/
3642: PetscErrorCode PetscDSCopyBoundary(PetscDS ds, PetscInt numFields, const PetscInt fields[], PetscDS newds)
3643: {
3644: DSBoundary b, *lastnext;
3646: PetscFunctionBegin;
3649: if (ds == newds) PetscFunctionReturn(PETSC_SUCCESS);
3650: PetscCall(PetscDSDestroyBoundary(newds));
3651: lastnext = &newds->boundary;
3652: for (b = ds->boundary; b; b = b->next) {
3653: DSBoundary bNew;
3654: PetscInt fieldNew = -1;
3656: if (numFields > 0 && fields) {
3657: PetscInt f;
3659: for (f = 0; f < numFields; ++f)
3660: if (b->field == fields[f]) break;
3661: if (f == numFields) continue;
3662: fieldNew = f;
3663: }
3664: PetscCall(DSBoundaryDuplicate_Internal(b, &bNew));
3665: bNew->field = fieldNew < 0 ? b->field : fieldNew;
3666: *lastnext = bNew;
3667: lastnext = &bNew->next;
3668: }
3669: PetscFunctionReturn(PETSC_SUCCESS);
3670: }
3672: /*@
3673: PetscDSDestroyBoundary - Remove all `DMBoundary` objects from the `PetscDS`
3675: Not Collective
3677: Input Parameter:
3678: . ds - The `PetscDS` object
3680: Level: intermediate
3682: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`
3683: @*/
3684: PetscErrorCode PetscDSDestroyBoundary(PetscDS ds)
3685: {
3686: DSBoundary next = ds->boundary;
3688: PetscFunctionBegin;
3689: while (next) {
3690: DSBoundary b = next;
3692: next = b->next;
3693: PetscCall(PetscWeakFormDestroy(&b->wf));
3694: PetscCall(PetscFree(b->name));
3695: PetscCall(PetscFree(b->lname));
3696: PetscCall(PetscFree(b->values));
3697: PetscCall(PetscFree(b->comps));
3698: PetscCall(PetscFree(b));
3699: }
3700: PetscFunctionReturn(PETSC_SUCCESS);
3701: }
3703: /*@
3704: PetscDSSelectDiscretizations - Copy discretizations to the new `PetscDS` with different field layout
3706: Not Collective
3708: Input Parameters:
3709: + prob - The `PetscDS` object
3710: . numFields - Number of new fields
3711: . fields - Old field number for each new field
3712: . minDegree - Minimum degree for a discretization, or `PETSC_DETERMINE` for no limit
3713: - maxDegree - Maximum degree for a discretization, or `PETSC_DETERMINE` for no limit
3715: Output Parameter:
3716: . newprob - The `PetscDS` copy
3718: Level: intermediate
3720: .seealso: `PetscDS`, `PetscDSSelectEquations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3721: @*/
3722: PetscErrorCode PetscDSSelectDiscretizations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscInt minDegree, PetscInt maxDegree, PetscDS newprob)
3723: {
3724: PetscInt Nf, Nfn, fn;
3726: PetscFunctionBegin;
3728: if (fields) PetscAssertPointer(fields, 3);
3730: PetscCall(PetscDSGetNumFields(prob, &Nf));
3731: PetscCall(PetscDSGetNumFields(newprob, &Nfn));
3732: numFields = numFields < 0 ? Nf : numFields;
3733: for (fn = 0; fn < numFields; ++fn) {
3734: const PetscInt f = fields ? fields[fn] : fn;
3735: PetscObject disc;
3736: PetscClassId id;
3738: if (f >= Nf) continue;
3739: PetscCall(PetscDSGetDiscretization(prob, f, &disc));
3740: PetscCallContinue(PetscObjectGetClassId(disc, &id));
3741: if (id == PETSCFE_CLASSID) {
3742: PetscFE fe;
3744: PetscCall(PetscFELimitDegree((PetscFE)disc, minDegree, maxDegree, &fe));
3745: PetscCall(PetscDSSetDiscretization(newprob, fn, (PetscObject)fe));
3746: PetscCall(PetscFEDestroy(&fe));
3747: } else {
3748: PetscCall(PetscDSSetDiscretization(newprob, fn, disc));
3749: }
3750: }
3751: PetscFunctionReturn(PETSC_SUCCESS);
3752: }
3754: /*@
3755: PetscDSSelectEquations - Copy pointwise function pointers to the new `PetscDS` with different field layout
3757: Not Collective
3759: Input Parameters:
3760: + prob - The `PetscDS` object
3761: . numFields - Number of new fields
3762: - fields - Old field number for each new field
3764: Output Parameter:
3765: . newprob - The `PetscDS` copy
3767: Level: intermediate
3769: .seealso: `PetscDS`, `PetscDSSelectDiscretizations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3770: @*/
3771: PetscErrorCode PetscDSSelectEquations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
3772: {
3773: PetscInt Nf, Nfn, fn, gn;
3775: PetscFunctionBegin;
3777: if (fields) PetscAssertPointer(fields, 3);
3779: PetscCall(PetscDSGetNumFields(prob, &Nf));
3780: PetscCall(PetscDSGetNumFields(newprob, &Nfn));
3781: PetscCheck(numFields <= Nfn, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_SIZ, "Number of fields %" PetscInt_FMT " to transfer must not be greater than the total number of fields %" PetscInt_FMT, numFields, Nfn);
3782: for (fn = 0; fn < numFields; ++fn) {
3783: const PetscInt f = fields ? fields[fn] : fn;
3784: PetscPointFn *obj;
3785: PetscPointFn *f0, *f1;
3786: PetscBdPointFn *f0Bd, *f1Bd;
3787: PetscRiemannFn *r;
3789: if (f >= Nf) continue;
3790: PetscCall(PetscDSGetObjective(prob, f, &obj));
3791: PetscCall(PetscDSGetResidual(prob, f, &f0, &f1));
3792: PetscCall(PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd));
3793: PetscCall(PetscDSGetRiemannSolver(prob, f, &r));
3794: PetscCall(PetscDSSetObjective(newprob, fn, obj));
3795: PetscCall(PetscDSSetResidual(newprob, fn, f0, f1));
3796: PetscCall(PetscDSSetBdResidual(newprob, fn, f0Bd, f1Bd));
3797: PetscCall(PetscDSSetRiemannSolver(newprob, fn, r));
3798: for (gn = 0; gn < numFields; ++gn) {
3799: const PetscInt g = fields ? fields[gn] : gn;
3800: PetscPointJacFn *g0, *g1, *g2, *g3;
3801: PetscPointJacFn *g0p, *g1p, *g2p, *g3p;
3802: PetscBdPointJacFn *g0Bd, *g1Bd, *g2Bd, *g3Bd;
3804: if (g >= Nf) continue;
3805: PetscCall(PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3));
3806: PetscCall(PetscDSGetJacobianPreconditioner(prob, f, g, &g0p, &g1p, &g2p, &g3p));
3807: PetscCall(PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd));
3808: PetscCall(PetscDSSetJacobian(newprob, fn, gn, g0, g1, g2, g3));
3809: PetscCall(PetscDSSetJacobianPreconditioner(newprob, fn, gn, g0p, g1p, g2p, g3p));
3810: PetscCall(PetscDSSetBdJacobian(newprob, fn, gn, g0Bd, g1Bd, g2Bd, g3Bd));
3811: }
3812: }
3813: PetscFunctionReturn(PETSC_SUCCESS);
3814: }
3816: /*@
3817: PetscDSCopyEquations - Copy all pointwise function pointers to another `PetscDS`
3819: Not Collective
3821: Input Parameter:
3822: . prob - The `PetscDS` object
3824: Output Parameter:
3825: . newprob - The `PetscDS` copy
3827: Level: intermediate
3829: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3830: @*/
3831: PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob)
3832: {
3833: PetscWeakForm wf, newwf;
3834: PetscInt Nf, Ng;
3836: PetscFunctionBegin;
3839: PetscCall(PetscDSGetNumFields(prob, &Nf));
3840: PetscCall(PetscDSGetNumFields(newprob, &Ng));
3841: PetscCheck(Nf == Ng, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %" PetscInt_FMT " != %" PetscInt_FMT, Nf, Ng);
3842: PetscCall(PetscDSGetWeakForm(prob, &wf));
3843: PetscCall(PetscDSGetWeakForm(newprob, &newwf));
3844: PetscCall(PetscWeakFormCopy(wf, newwf));
3845: PetscFunctionReturn(PETSC_SUCCESS);
3846: }
3848: /*@
3849: PetscDSCopyConstants - Copy all constants set with `PetscDSSetConstants()` to another `PetscDS`
3851: Not Collective
3853: Input Parameter:
3854: . prob - The `PetscDS` object
3856: Output Parameter:
3857: . newprob - The `PetscDS` copy
3859: Level: intermediate
3861: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3862: @*/
3863: PetscErrorCode PetscDSCopyConstants(PetscDS prob, PetscDS newprob)
3864: {
3865: PetscInt Nc;
3866: const PetscScalar *constants;
3868: PetscFunctionBegin;
3871: PetscCall(PetscDSGetConstants(prob, &Nc, &constants));
3872: PetscCall(PetscDSSetConstants(newprob, Nc, (PetscScalar *)constants));
3873: PetscFunctionReturn(PETSC_SUCCESS);
3874: }
3876: /*@
3877: PetscDSCopyExactSolutions - Copy all exact solutions set with `PetscDSSetExactSolution()` and `PetscDSSetExactSolutionTimeDerivative()` to another `PetscDS`
3879: Not Collective
3881: Input Parameter:
3882: . ds - The `PetscDS` object
3884: Output Parameter:
3885: . newds - The `PetscDS` copy
3887: Level: intermediate
3889: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSCopyBounds()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3890: @*/
3891: PetscErrorCode PetscDSCopyExactSolutions(PetscDS ds, PetscDS newds)
3892: {
3893: PetscSimplePointFn *sol;
3894: void *ctx;
3895: PetscInt Nf, f;
3897: PetscFunctionBegin;
3900: PetscCall(PetscDSGetNumFields(ds, &Nf));
3901: for (f = 0; f < Nf; ++f) {
3902: PetscCall(PetscDSGetExactSolution(ds, f, &sol, &ctx));
3903: PetscCall(PetscDSSetExactSolution(newds, f, sol, ctx));
3904: PetscCall(PetscDSGetExactSolutionTimeDerivative(ds, f, &sol, &ctx));
3905: PetscCall(PetscDSSetExactSolutionTimeDerivative(newds, f, sol, ctx));
3906: }
3907: PetscFunctionReturn(PETSC_SUCCESS);
3908: }
3910: /*@
3911: PetscDSCopyBounds - Copy lower and upper solution bounds set with `PetscDSSetLowerBound()` and `PetscDSSetLowerBound()` to another `PetscDS`
3913: Not Collective
3915: Input Parameter:
3916: . ds - The `PetscDS` object
3918: Output Parameter:
3919: . newds - The `PetscDS` copy
3921: Level: intermediate
3923: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSCopyExactSolutions()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3924: @*/
3925: PetscErrorCode PetscDSCopyBounds(PetscDS ds, PetscDS newds)
3926: {
3927: PetscSimplePointFn *bound;
3928: void *ctx;
3929: PetscInt Nf, f;
3931: PetscFunctionBegin;
3934: PetscCall(PetscDSGetNumFields(ds, &Nf));
3935: for (f = 0; f < Nf; ++f) {
3936: PetscCall(PetscDSGetLowerBound(ds, f, &bound, &ctx));
3937: PetscCall(PetscDSSetLowerBound(newds, f, bound, ctx));
3938: PetscCall(PetscDSGetUpperBound(ds, f, &bound, &ctx));
3939: PetscCall(PetscDSSetUpperBound(newds, f, bound, ctx));
3940: }
3941: PetscFunctionReturn(PETSC_SUCCESS);
3942: }
3944: /*@
3945: PetscDSCopy - Copy the contents of a `PetscDS` into another `PetscDS` on a new `DM`.
3947: Collective
3949: Input Parameters:
3950: + ds - the source `PetscDS`
3951: . minDegree - the minimum polynomial degree to consider when selecting discretizations, or `PETSC_DETERMINE`
3952: . maxDegree - the maximum polynomial degree to consider when selecting discretizations, or `PETSC_DETERMINE`
3953: - dmNew - the target `DM` used to resolve boundary condition labels for the copied boundaries
3955: Output Parameter:
3956: . dsNew - the destination `PetscDS`
3958: Level: developer
3960: Note:
3961: This copies constants, exact solutions, bounds, discretizations, equations, field contexts,
3962: cohesive flags, jet degrees, and boundary conditions.
3964: .seealso: `PetscDS`, `PetscDSCopyEquations()`, `PetscDSCopyConstants()`, `PetscDSCopyExactSolutions()`, `PetscDSCopyBounds()`, `PetscDSCopyBoundary()`
3965: @*/
3966: PetscErrorCode PetscDSCopy(PetscDS ds, PetscInt minDegree, PetscInt maxDegree, DM dmNew, PetscDS dsNew)
3967: {
3968: DSBoundary b;
3969: PetscInt cdim, Nf, f, d;
3970: PetscBool isCohesive;
3971: void *ctx;
3973: PetscFunctionBegin;
3974: PetscCall(PetscDSCopyConstants(ds, dsNew));
3975: PetscCall(PetscDSCopyExactSolutions(ds, dsNew));
3976: PetscCall(PetscDSCopyBounds(ds, dsNew));
3977: PetscCall(PetscDSSelectDiscretizations(ds, PETSC_DETERMINE, NULL, minDegree, maxDegree, dsNew));
3978: PetscCall(PetscDSCopyEquations(ds, dsNew));
3979: PetscCall(PetscDSGetNumFields(ds, &Nf));
3980: for (f = 0; f < Nf; ++f) {
3981: PetscCall(PetscDSGetContext(ds, f, &ctx));
3982: PetscCall(PetscDSSetContext(dsNew, f, ctx));
3983: PetscCall(PetscDSGetCohesive(ds, f, &isCohesive));
3984: PetscCall(PetscDSSetCohesive(dsNew, f, isCohesive));
3985: PetscCall(PetscDSGetJetDegree(ds, f, &d));
3986: PetscCall(PetscDSSetJetDegree(dsNew, f, d));
3987: }
3988: if (Nf) {
3989: PetscCall(PetscDSGetCoordinateDimension(ds, &cdim));
3990: PetscCall(PetscDSSetCoordinateDimension(dsNew, cdim));
3991: }
3992: PetscCall(PetscDSCopyBoundary(ds, PETSC_DETERMINE, NULL, dsNew));
3993: for (b = dsNew->boundary; b; b = b->next) {
3994: PetscCall(DMGetLabel(dmNew, b->lname, &b->label));
3995: /* Do not check if label exists here, since p4est calls this for the reference tree which does not have the labels */
3996: //PetscCheck(b->label,PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Label %s missing in new DM", name);
3997: }
3998: PetscFunctionReturn(PETSC_SUCCESS);
3999: }
4001: /*@
4002: PetscDSGetHeightSubspace - Get the `PetscDS` for the trace subspace at a given height in the mesh.
4004: Not Collective
4006: Input Parameters:
4007: + prob - the `PetscDS`
4008: - height - the height (0 for the ambient cell, 1 for faces, etc.)
4010: Output Parameter:
4011: . subprob - the `PetscDS` for the trace subspace; `prob` itself is returned when `height` is 0
4013: Level: developer
4015: Note:
4016: Only `PetscFE` discretizations are currently supported.
4018: .seealso: `PetscDS`, `PetscFE`, `PetscFEGetHeightSubspace()`, `PetscDSGetSpatialDimension()`
4019: @*/
4020: PetscErrorCode PetscDSGetHeightSubspace(PetscDS prob, PetscInt height, PetscDS *subprob)
4021: {
4022: PetscInt dim, Nf, f;
4024: PetscFunctionBegin;
4026: PetscAssertPointer(subprob, 3);
4027: if (height == 0) {
4028: *subprob = prob;
4029: PetscFunctionReturn(PETSC_SUCCESS);
4030: }
4031: PetscCall(PetscDSGetNumFields(prob, &Nf));
4032: PetscCall(PetscDSGetSpatialDimension(prob, &dim));
4033: PetscCheck(height <= dim, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_OUTOFRANGE, "DS can only handle height in [0, %" PetscInt_FMT "], not %" PetscInt_FMT, dim, height);
4034: if (!prob->subprobs) PetscCall(PetscCalloc1(dim, &prob->subprobs));
4035: if (!prob->subprobs[height - 1]) {
4036: PetscInt cdim;
4038: PetscCall(PetscDSCreate(PetscObjectComm((PetscObject)prob), &prob->subprobs[height - 1]));
4039: PetscCall(PetscDSGetCoordinateDimension(prob, &cdim));
4040: PetscCall(PetscDSSetCoordinateDimension(prob->subprobs[height - 1], cdim));
4041: for (f = 0; f < Nf; ++f) {
4042: PetscFE subfe;
4043: PetscObject obj;
4044: PetscClassId id;
4046: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
4047: PetscCall(PetscObjectGetClassId(obj, &id));
4048: PetscCheck(id == PETSCFE_CLASSID, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unsupported discretization type for field %" PetscInt_FMT, f);
4049: PetscCall(PetscFEGetHeightSubspace((PetscFE)obj, height, &subfe));
4050: PetscCall(PetscDSSetDiscretization(prob->subprobs[height - 1], f, (PetscObject)subfe));
4051: }
4052: }
4053: *subprob = prob->subprobs[height - 1];
4054: PetscFunctionReturn(PETSC_SUCCESS);
4055: }
4057: /*@
4058: PetscDSPermuteQuadPoint - Permute a quadrature point index according to a cell orientation.
4060: Not Collective
4062: Input Parameters:
4063: + ds - the `PetscDS`
4064: . ornt - the cell orientation, in `[-Na, Na)` where `Na` is half the number of arrangements for the cell type
4065: . field - the field number whose quadrature is used
4066: - q - the input quadrature point index in `[0, Nq)`
4068: Output Parameter:
4069: . qperm - the permuted quadrature point index
4071: Level: developer
4073: .seealso: `PetscDS`, `PetscQuadrature`, `PetscQuadratureComputePermutations()`, `DMPolytopeTypeGetNumArrangements()`
4074: @*/
4075: PetscErrorCode PetscDSPermuteQuadPoint(PetscDS ds, PetscInt ornt, PetscInt field, PetscInt q, PetscInt *qperm)
4076: {
4077: IS permIS;
4078: PetscQuadrature quad;
4079: DMPolytopeType ct;
4080: const PetscInt *perm;
4081: PetscInt Na, Nq;
4083: PetscFunctionBeginHot;
4084: PetscCall(PetscFEGetQuadrature((PetscFE)ds->disc[field], &quad));
4085: PetscCall(PetscQuadratureGetData(quad, NULL, NULL, &Nq, NULL, NULL));
4086: PetscCall(PetscQuadratureGetCellType(quad, &ct));
4087: PetscCheck(q >= 0 && q < Nq, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Quadrature point %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", q, Nq);
4088: Na = DMPolytopeTypeGetNumArrangements(ct) / 2;
4089: 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);
4090: if (!ds->quadPerm[(PetscInt)ct]) PetscCall(PetscQuadratureComputePermutations(quad, NULL, &ds->quadPerm[(PetscInt)ct]));
4091: permIS = ds->quadPerm[(PetscInt)ct][ornt + Na];
4092: PetscCall(ISGetIndices(permIS, &perm));
4093: *qperm = perm[q];
4094: PetscCall(ISRestoreIndices(permIS, &perm));
4095: PetscFunctionReturn(PETSC_SUCCESS);
4096: }
4098: PetscErrorCode PetscDSGetDiscType_Internal(PetscDS ds, PetscInt f, PetscDiscType *disctype)
4099: {
4100: PetscObject obj;
4101: PetscClassId id;
4102: PetscInt Nf;
4104: PetscFunctionBegin;
4106: PetscAssertPointer(disctype, 3);
4107: *disctype = PETSC_DISC_NONE;
4108: PetscCall(PetscDSGetNumFields(ds, &Nf));
4109: PetscCheck(f < Nf, PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_SIZ, "Field %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, Nf);
4110: PetscCall(PetscDSGetDiscretization(ds, f, &obj));
4111: if (obj) {
4112: PetscCall(PetscObjectGetClassId(obj, &id));
4113: if (id == PETSCFE_CLASSID) *disctype = PETSC_DISC_FE;
4114: else *disctype = PETSC_DISC_FV;
4115: }
4116: PetscFunctionReturn(PETSC_SUCCESS);
4117: }
4119: static PetscErrorCode PetscDSDestroy_Basic(PetscDS ds)
4120: {
4121: PetscFunctionBegin;
4122: PetscCall(PetscFree(ds->data));
4123: PetscFunctionReturn(PETSC_SUCCESS);
4124: }
4126: static PetscErrorCode PetscDSInitialize_Basic(PetscDS ds)
4127: {
4128: PetscFunctionBegin;
4129: ds->ops->setfromoptions = NULL;
4130: ds->ops->setup = NULL;
4131: ds->ops->view = NULL;
4132: ds->ops->destroy = PetscDSDestroy_Basic;
4133: PetscFunctionReturn(PETSC_SUCCESS);
4134: }
4136: /*MC
4137: PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions
4139: Level: intermediate
4141: .seealso: `PetscDSType`, `PetscDSCreate()`, `PetscDSSetType()`
4142: M*/
4144: PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS ds)
4145: {
4146: PetscDS_Basic *b;
4148: PetscFunctionBegin;
4150: PetscCall(PetscNew(&b));
4151: ds->data = b;
4153: PetscCall(PetscDSInitialize_Basic(ds));
4154: PetscFunctionReturn(PETSC_SUCCESS);
4155: }