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: /*@C
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()`, `PetscDSRegisterDestroy()`
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;
162: PetscInt c, i;
164: if (b->field != f) continue;
165: PetscCall(PetscViewerASCIIPushTab(viewer));
166: PetscCall(PetscViewerASCIIPrintf(viewer, "Boundary %s (%s) %s\n", b->name, b->lname, DMBoundaryConditionTypes[b->type]));
167: if (!b->Nc) {
168: PetscCall(PetscViewerASCIIPrintf(viewer, " all components\n"));
169: } else {
170: PetscCall(PetscViewerASCIIPrintf(viewer, " components: "));
171: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
172: for (c = 0; c < b->Nc; ++c) {
173: if (c > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", "));
174: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->comps[c]));
175: }
176: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
177: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
178: }
179: PetscCall(PetscViewerASCIIPrintf(viewer, " values: "));
180: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_FALSE));
181: for (i = 0; i < b->Nv; ++i) {
182: if (i > 0) PetscCall(PetscViewerASCIIPrintf(viewer, ", "));
183: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT, b->values[i]));
184: }
185: PetscCall(PetscViewerASCIIPrintf(viewer, "\n"));
186: PetscCall(PetscViewerASCIIUseTabs(viewer, PETSC_TRUE));
187: #if defined(__clang__)
188: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat-pedantic")
189: #elif defined(__GNUC__) || defined(__GNUG__)
190: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_BEGIN("-Wformat")
191: #endif
192: if (b->func) {
193: PetscCall(PetscDLAddr(b->func, &name));
194: if (name) PetscCall(PetscViewerASCIIPrintf(viewer, " func: %s\n", name));
195: else PetscCall(PetscViewerASCIIPrintf(viewer, " func: %p\n", b->func));
196: PetscCall(PetscFree(name));
197: }
198: if (b->func_t) {
199: PetscCall(PetscDLAddr(b->func_t, &name));
200: if (name) PetscCall(PetscViewerASCIIPrintf(viewer, " func_t: %s\n", name));
201: else PetscCall(PetscViewerASCIIPrintf(viewer, " func_t: %p\n", b->func_t));
202: PetscCall(PetscFree(name));
203: }
204: PETSC_PRAGMA_DIAGNOSTIC_IGNORED_END()
205: PetscCall(PetscWeakFormView(b->wf, viewer));
206: PetscCall(PetscViewerASCIIPopTab(viewer));
207: }
208: }
209: PetscCall(PetscDSGetConstants(ds, &numConstants, &constants));
210: if (numConstants) {
211: PetscCall(PetscViewerASCIIPrintf(viewer, "%" PetscInt_FMT " constants\n", numConstants));
212: PetscCall(PetscViewerASCIIPushTab(viewer));
213: for (f = 0; f < numConstants; ++f) PetscCall(PetscViewerASCIIPrintf(viewer, "%g\n", (double)PetscRealPart(constants[f])));
214: PetscCall(PetscViewerASCIIPopTab(viewer));
215: }
216: PetscCall(PetscWeakFormView(ds->wf, viewer));
217: PetscCall(PetscViewerASCIIPopTab(viewer));
218: PetscFunctionReturn(PETSC_SUCCESS);
219: }
221: /*@
222: PetscDSViewFromOptions - View a `PetscDS` based on values in the options database
224: Collective
226: Input Parameters:
227: + A - the `PetscDS` object
228: . obj - Optional object that provides the options prefix used in the search of the options database
229: - name - command line option
231: Level: intermediate
233: .seealso: `PetscDSType`, `PetscDS`, `PetscDSView()`, `PetscObjectViewFromOptions()`, `PetscDSCreate()`
234: @*/
235: PetscErrorCode PetscDSViewFromOptions(PetscDS A, PetscObject obj, const char name[])
236: {
237: PetscFunctionBegin;
239: PetscCall(PetscObjectViewFromOptions((PetscObject)A, obj, name));
240: PetscFunctionReturn(PETSC_SUCCESS);
241: }
243: /*@
244: PetscDSView - Views a `PetscDS`
246: Collective
248: Input Parameters:
249: + prob - the `PetscDS` object to view
250: - v - the viewer
252: Level: developer
254: .seealso: `PetscDSType`, `PetscDS`, `PetscViewer`, `PetscDSDestroy()`, `PetscDSViewFromOptions()`
255: @*/
256: PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v)
257: {
258: PetscBool isascii;
260: PetscFunctionBegin;
262: if (!v) PetscCall(PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject)prob), &v));
264: PetscCall(PetscObjectTypeCompare((PetscObject)v, PETSCVIEWERASCII, &isascii));
265: if (isascii) PetscCall(PetscDSView_Ascii(prob, v));
266: PetscTryTypeMethod(prob, view, v);
267: PetscFunctionReturn(PETSC_SUCCESS);
268: }
270: /*@
271: PetscDSSetFromOptions - sets parameters in a `PetscDS` from the options database
273: Collective
275: Input Parameter:
276: . prob - the `PetscDS` object to set options for
278: Options Database Keys:
279: + -petscds_type <type> - Set the `PetscDS` type
280: . -petscds_view <view opt> - View the `PetscDS`
281: . -petscds_jac_pre - Turn formation of a separate Jacobian preconditioner on or off
282: . -bc_<name> <ids> - Specify a list of label ids for a boundary condition
283: - -bc_<name>_comp <comps> - Specify a list of field components to constrain for a boundary condition
285: Level: intermediate
287: .seealso: `PetscDS`, `PetscDSView()`
288: @*/
289: PetscErrorCode PetscDSSetFromOptions(PetscDS prob)
290: {
291: DSBoundary b;
292: const char *defaultType;
293: char name[256];
294: PetscBool flg;
296: PetscFunctionBegin;
298: if (!((PetscObject)prob)->type_name) {
299: defaultType = PETSCDSBASIC;
300: } else {
301: defaultType = ((PetscObject)prob)->type_name;
302: }
303: PetscCall(PetscDSRegisterAll());
305: PetscObjectOptionsBegin((PetscObject)prob);
306: for (b = prob->boundary; b; b = b->next) {
307: char optname[1024];
308: PetscInt ids[1024], len = 1024;
309: PetscBool flg;
311: PetscCall(PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name));
312: PetscCall(PetscMemzero(ids, sizeof(ids)));
313: PetscCall(PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg));
314: if (flg) {
315: b->Nv = len;
316: PetscCall(PetscFree(b->values));
317: PetscCall(PetscMalloc1(len, &b->values));
318: PetscCall(PetscArraycpy(b->values, ids, len));
319: PetscCall(PetscWeakFormRewriteKeys(b->wf, b->label, len, b->values));
320: }
321: len = 1024;
322: PetscCall(PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name));
323: PetscCall(PetscMemzero(ids, sizeof(ids)));
324: PetscCall(PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg));
325: if (flg) {
326: b->Nc = len;
327: PetscCall(PetscFree(b->comps));
328: PetscCall(PetscMalloc1(len, &b->comps));
329: PetscCall(PetscArraycpy(b->comps, ids, len));
330: }
331: }
332: PetscCall(PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg));
333: if (flg) {
334: PetscCall(PetscDSSetType(prob, name));
335: } else if (!((PetscObject)prob)->type_name) {
336: PetscCall(PetscDSSetType(prob, defaultType));
337: }
338: PetscCall(PetscOptionsBool("-petscds_jac_pre", "Discrete System", "PetscDSUseJacobianPreconditioner", prob->useJacPre, &prob->useJacPre, &flg));
339: PetscCall(PetscOptionsBool("-petscds_force_quad", "Discrete System", "PetscDSSetForceQuad", prob->forceQuad, &prob->forceQuad, &flg));
340: PetscCall(PetscOptionsInt("-petscds_print_integrate", "Discrete System", "", prob->printIntegrate, &prob->printIntegrate, NULL));
341: PetscTryTypeMethod(prob, setfromoptions);
342: /* process any options handlers added with PetscObjectAddOptionsHandler() */
343: PetscCall(PetscObjectProcessOptionsHandlers((PetscObject)prob, PetscOptionsObject));
344: PetscOptionsEnd();
345: if (prob->Nf) PetscCall(PetscDSViewFromOptions(prob, NULL, "-petscds_view"));
346: PetscFunctionReturn(PETSC_SUCCESS);
347: }
349: /*@
350: PetscDSSetUp - Construct data structures for the `PetscDS`
352: Collective
354: Input Parameter:
355: . prob - the `PetscDS` object to setup
357: Level: developer
359: .seealso: `PetscDS`, `PetscDSView()`, `PetscDSDestroy()`
360: @*/
361: PetscErrorCode PetscDSSetUp(PetscDS prob)
362: {
363: const PetscInt Nf = prob->Nf;
364: PetscBool hasH = PETSC_FALSE;
365: PetscInt maxOrder[4] = {-2, -2, -2, -2};
366: PetscInt dim, dimEmbed, NbMax = 0, NcMax = 0, NqMax = 0, NsMax = 1, f;
368: PetscFunctionBegin;
370: if (prob->setup) PetscFunctionReturn(PETSC_SUCCESS);
371: /* Calculate sizes */
372: PetscCall(PetscDSGetSpatialDimension(prob, &dim));
373: PetscCall(PetscDSGetCoordinateDimension(prob, &dimEmbed));
374: prob->totDim = prob->totComp = 0;
375: PetscCall(PetscMalloc2(Nf, &prob->Nc, Nf, &prob->Nb));
376: PetscCall(PetscCalloc2(Nf + 1, &prob->off, Nf + 1, &prob->offDer));
377: 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]));
378: PetscCall(PetscMalloc2(Nf, &prob->T, Nf, &prob->Tf));
379: if (prob->forceQuad) {
380: // Note: This assumes we have one kind of cell at each dimension.
381: // We can fix this by having quadrature hold the celltype
382: PetscQuadrature maxQuad[4] = {NULL, NULL, NULL, NULL};
384: for (f = 0; f < Nf; ++f) {
385: PetscObject obj;
386: PetscClassId id;
387: PetscQuadrature q = NULL, fq = NULL;
388: PetscInt dim = -1, order = -1, forder = -1;
390: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
391: if (!obj) continue;
392: PetscCall(PetscObjectGetClassId(obj, &id));
393: if (id == PETSCFE_CLASSID) {
394: PetscFE fe = (PetscFE)obj;
396: PetscCall(PetscFEGetQuadrature(fe, &q));
397: PetscCall(PetscFEGetFaceQuadrature(fe, &fq));
398: } else if (id == PETSCFV_CLASSID) {
399: PetscFV fv = (PetscFV)obj;
401: PetscCall(PetscFVGetQuadrature(fv, &q));
402: }
403: if (q) {
404: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
405: PetscCall(PetscQuadratureGetOrder(q, &order));
406: if (order > maxOrder[dim]) {
407: maxOrder[dim] = order;
408: maxQuad[dim] = q;
409: }
410: }
411: if (fq) {
412: PetscCall(PetscQuadratureGetData(fq, &dim, NULL, NULL, NULL, NULL));
413: PetscCall(PetscQuadratureGetOrder(fq, &forder));
414: if (forder > maxOrder[dim]) {
415: maxOrder[dim] = forder;
416: maxQuad[dim] = fq;
417: }
418: }
419: }
420: for (f = 0; f < Nf; ++f) {
421: PetscObject obj;
422: PetscClassId id;
423: PetscQuadrature q;
424: PetscInt dim;
426: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
427: if (!obj) continue;
428: PetscCall(PetscObjectGetClassId(obj, &id));
429: if (id == PETSCFE_CLASSID) {
430: PetscFE fe = (PetscFE)obj;
432: PetscCall(PetscFEGetQuadrature(fe, &q));
433: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
434: PetscCall(PetscFESetQuadrature(fe, maxQuad[dim]));
435: PetscCall(PetscFESetFaceQuadrature(fe, dim ? maxQuad[dim - 1] : NULL));
436: } else if (id == PETSCFV_CLASSID) {
437: PetscFV fv = (PetscFV)obj;
439: PetscCall(PetscFVGetQuadrature(fv, &q));
440: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
441: PetscCall(PetscFVSetQuadrature(fv, maxQuad[dim]));
442: }
443: }
444: }
445: for (f = 0; f < Nf; ++f) {
446: PetscObject obj;
447: PetscClassId id;
448: PetscQuadrature q = NULL;
449: PetscInt Nq = 0, Nb, Nc;
451: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
452: if (prob->jetDegree[f] > 1) hasH = PETSC_TRUE;
453: if (!obj) {
454: /* Empty mesh */
455: Nb = Nc = 0;
456: prob->T[f] = prob->Tf[f] = NULL;
457: } else {
458: PetscCall(PetscObjectGetClassId(obj, &id));
459: if (id == PETSCFE_CLASSID) {
460: PetscFE fe = (PetscFE)obj;
462: PetscCall(PetscFEGetQuadrature(fe, &q));
463: {
464: PetscQuadrature fq;
465: PetscInt dim, order;
467: PetscCall(PetscQuadratureGetData(q, &dim, NULL, NULL, NULL, NULL));
468: PetscCall(PetscQuadratureGetOrder(q, &order));
469: if (maxOrder[dim] < 0) maxOrder[dim] = order;
470: 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]);
471: PetscCall(PetscFEGetFaceQuadrature(fe, &fq));
472: if (fq) {
473: PetscCall(PetscQuadratureGetData(fq, &dim, NULL, NULL, NULL, NULL));
474: PetscCall(PetscQuadratureGetOrder(fq, &order));
475: if (maxOrder[dim] < 0) maxOrder[dim] = order;
476: 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]);
477: }
478: }
479: PetscCall(PetscFEGetDimension(fe, &Nb));
480: PetscCall(PetscFEGetNumComponents(fe, &Nc));
481: PetscCall(PetscFEGetCellTabulation(fe, prob->jetDegree[f], &prob->T[f]));
482: PetscCall(PetscFEGetFaceTabulation(fe, prob->jetDegree[f], &prob->Tf[f]));
483: } else if (id == PETSCFV_CLASSID) {
484: PetscFV fv = (PetscFV)obj;
486: PetscCall(PetscFVGetQuadrature(fv, &q));
487: PetscCall(PetscFVGetNumComponents(fv, &Nc));
488: Nb = Nc;
489: PetscCall(PetscFVGetCellTabulation(fv, &prob->T[f]));
490: /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */
491: } else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %" PetscInt_FMT, f);
492: }
493: prob->Nc[f] = Nc;
494: prob->Nb[f] = Nb;
495: prob->off[f + 1] = Nc + prob->off[f];
496: prob->offDer[f + 1] = Nc * dim + prob->offDer[f];
497: prob->offCohesive[0][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[0][f];
498: prob->offDerCohesive[0][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[0][f];
499: prob->offCohesive[1][f] = (prob->cohesive[f] ? 0 : Nc) + prob->offCohesive[0][f];
500: prob->offDerCohesive[1][f] = (prob->cohesive[f] ? 0 : Nc) * dimEmbed + prob->offDerCohesive[0][f];
501: prob->offCohesive[2][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) + prob->offCohesive[2][f];
502: prob->offDerCohesive[2][f + 1] = (prob->cohesive[f] ? Nc : Nc * 2) * dimEmbed + prob->offDerCohesive[2][f];
503: if (q) PetscCall(PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL));
504: NqMax = PetscMax(NqMax, Nq);
505: NbMax = PetscMax(NbMax, Nb);
506: NcMax = PetscMax(NcMax, Nc);
507: prob->totDim += Nb;
508: prob->totComp += Nc;
509: /* There are two faces for all fields on a cohesive cell, except for cohesive fields */
510: if (prob->isCohesive && !prob->cohesive[f]) prob->totDim += Nb;
511: }
512: prob->offCohesive[1][Nf] = prob->offCohesive[0][Nf];
513: prob->offDerCohesive[1][Nf] = prob->offDerCohesive[0][Nf];
514: /* Allocate works space */
515: NsMax = 2; /* A non-cohesive discretizations can be used on a cohesive cell, so we need this extra workspace for all DS */
516: 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));
517: PetscCall(PetscMalloc5(dimEmbed, &prob->x, NbMax * NcMax, &prob->basisReal, NbMax * NcMax * dimEmbed, &prob->basisDerReal, NbMax * NcMax, &prob->testReal, NbMax * NcMax * dimEmbed, &prob->testDerReal));
518: 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,
519: &prob->g2, NsMax * NsMax * NqMax * NcMax * NcMax * dimEmbed * dimEmbed, &prob->g3));
520: PetscTryTypeMethod(prob, setup);
521: prob->setup = PETSC_TRUE;
522: PetscFunctionReturn(PETSC_SUCCESS);
523: }
525: static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob)
526: {
527: PetscFunctionBegin;
528: PetscCall(PetscFree2(prob->Nc, prob->Nb));
529: PetscCall(PetscFree2(prob->off, prob->offDer));
530: PetscCall(PetscFree6(prob->offCohesive[0], prob->offCohesive[1], prob->offCohesive[2], prob->offDerCohesive[0], prob->offDerCohesive[1], prob->offDerCohesive[2]));
531: PetscCall(PetscFree2(prob->T, prob->Tf));
532: PetscCall(PetscFree3(prob->u, prob->u_t, prob->u_x));
533: PetscCall(PetscFree5(prob->x, prob->basisReal, prob->basisDerReal, prob->testReal, prob->testDerReal));
534: PetscCall(PetscFree6(prob->f0, prob->f1, prob->g0, prob->g1, prob->g2, prob->g3));
535: PetscFunctionReturn(PETSC_SUCCESS);
536: }
538: static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew)
539: {
540: PetscObject *tmpd;
541: PetscBool *tmpi;
542: PetscInt *tmpk;
543: PetscBool *tmpc;
544: PetscPointFn **tmpup;
545: PetscSimplePointFn **tmpexactSol, **tmpexactSol_t, **tmplowerBound, **tmpupperBound;
546: void **tmpexactCtx, **tmpexactCtx_t, **tmplowerCtx, **tmpupperCtx;
547: void **tmpctx;
548: PetscInt Nf = prob->Nf, f;
550: PetscFunctionBegin;
551: if (Nf >= NfNew) PetscFunctionReturn(PETSC_SUCCESS);
552: prob->setup = PETSC_FALSE;
553: PetscCall(PetscDSDestroyStructs_Static(prob));
554: PetscCall(PetscMalloc4(NfNew, &tmpd, NfNew, &tmpi, NfNew, &tmpc, NfNew, &tmpk));
555: for (f = 0; f < Nf; ++f) {
556: tmpd[f] = prob->disc[f];
557: tmpi[f] = prob->implicit[f];
558: tmpc[f] = prob->cohesive[f];
559: tmpk[f] = prob->jetDegree[f];
560: }
561: for (f = Nf; f < NfNew; ++f) {
562: tmpd[f] = NULL;
563: tmpi[f] = PETSC_TRUE, tmpc[f] = PETSC_FALSE;
564: tmpk[f] = 1;
565: }
566: PetscCall(PetscFree4(prob->disc, prob->implicit, prob->cohesive, prob->jetDegree));
567: PetscCall(PetscWeakFormSetNumFields(prob->wf, NfNew));
568: prob->Nf = NfNew;
569: prob->disc = tmpd;
570: prob->implicit = tmpi;
571: prob->cohesive = tmpc;
572: prob->jetDegree = tmpk;
573: PetscCall(PetscCalloc2(NfNew, &tmpup, NfNew, &tmpctx));
574: for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f];
575: for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f];
576: for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL;
577: for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL;
578: PetscCall(PetscFree2(prob->update, prob->ctx));
579: prob->update = tmpup;
580: prob->ctx = tmpctx;
581: PetscCall(PetscCalloc4(NfNew, &tmpexactSol, NfNew, &tmpexactCtx, NfNew, &tmpexactSol_t, NfNew, &tmpexactCtx_t));
582: PetscCall(PetscCalloc4(NfNew, &tmplowerBound, NfNew, &tmplowerCtx, NfNew, &tmpupperBound, NfNew, &tmpupperCtx));
583: for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f];
584: for (f = 0; f < Nf; ++f) tmpexactCtx[f] = prob->exactCtx[f];
585: for (f = 0; f < Nf; ++f) tmpexactSol_t[f] = prob->exactSol_t[f];
586: for (f = 0; f < Nf; ++f) tmpexactCtx_t[f] = prob->exactCtx_t[f];
587: for (f = 0; f < Nf; ++f) tmplowerBound[f] = prob->lowerBound[f];
588: for (f = 0; f < Nf; ++f) tmplowerCtx[f] = prob->lowerCtx[f];
589: for (f = 0; f < Nf; ++f) tmpupperBound[f] = prob->upperBound[f];
590: for (f = 0; f < Nf; ++f) tmpupperCtx[f] = prob->upperCtx[f];
591: for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL;
592: for (f = Nf; f < NfNew; ++f) tmpexactCtx[f] = NULL;
593: for (f = Nf; f < NfNew; ++f) tmpexactSol_t[f] = NULL;
594: for (f = Nf; f < NfNew; ++f) tmpexactCtx_t[f] = NULL;
595: for (f = Nf; f < NfNew; ++f) tmplowerBound[f] = NULL;
596: for (f = Nf; f < NfNew; ++f) tmplowerCtx[f] = NULL;
597: for (f = Nf; f < NfNew; ++f) tmpupperBound[f] = NULL;
598: for (f = Nf; f < NfNew; ++f) tmpupperCtx[f] = NULL;
599: PetscCall(PetscFree4(prob->exactSol, prob->exactCtx, prob->exactSol_t, prob->exactCtx_t));
600: PetscCall(PetscFree4(prob->lowerBound, prob->lowerCtx, prob->upperBound, prob->upperCtx));
601: prob->exactSol = tmpexactSol;
602: prob->exactCtx = tmpexactCtx;
603: prob->exactSol_t = tmpexactSol_t;
604: prob->exactCtx_t = tmpexactCtx_t;
605: prob->lowerBound = tmplowerBound;
606: prob->lowerCtx = tmplowerCtx;
607: prob->upperBound = tmpupperBound;
608: prob->upperCtx = tmpupperCtx;
609: PetscFunctionReturn(PETSC_SUCCESS);
610: }
612: /*@
613: PetscDSDestroy - Destroys a `PetscDS` object
615: Collective
617: Input Parameter:
618: . ds - the `PetscDS` object to destroy
620: Level: developer
622: .seealso: `PetscDSView()`
623: @*/
624: PetscErrorCode PetscDSDestroy(PetscDS *ds)
625: {
626: PetscInt f;
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, d;
640: PetscCall(PetscDSGetSpatialDimension(*ds, &dim));
641: for (d = 0; d < dim; ++d) PetscCall(PetscDSDestroy(&(*ds)->subprobs[d]));
642: }
643: PetscCall(PetscFree((*ds)->subprobs));
644: PetscCall(PetscDSDestroyStructs_Static(*ds));
645: for (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: PetscInt f;
898: PetscFunctionBegin;
900: PetscAssertPointer(numCohesive, 2);
901: *numCohesive = 0;
902: for (f = 0; f < ds->Nf; ++f) *numCohesive += ds->cohesive[f] ? 1 : 0;
903: PetscFunctionReturn(PETSC_SUCCESS);
904: }
906: /*@
907: PetscDSGetCohesive - Returns the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell
909: Not Collective
911: Input Parameters:
912: + ds - The `PetscDS` object
913: - f - The field index
915: Output Parameter:
916: . isCohesive - The flag
918: Level: developer
920: .seealso: `PetscDS`, `PetscDSSetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
921: @*/
922: PetscErrorCode PetscDSGetCohesive(PetscDS ds, PetscInt f, PetscBool *isCohesive)
923: {
924: PetscFunctionBegin;
926: PetscAssertPointer(isCohesive, 3);
927: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
928: *isCohesive = ds->cohesive[f];
929: PetscFunctionReturn(PETSC_SUCCESS);
930: }
932: /*@
933: PetscDSSetCohesive - Set the flag indicating that a field is cohesive, meaning it is defined on the interior of a cohesive cell
935: Not Collective
937: Input Parameters:
938: + ds - The `PetscDS` object
939: . f - The field index
940: - isCohesive - The flag for a cohesive field
942: Level: developer
944: .seealso: `PetscDS`, `PetscDSGetCohesive()`, `PetscDSIsCohesive()`, `PetscDSCreate()`
945: @*/
946: PetscErrorCode PetscDSSetCohesive(PetscDS ds, PetscInt f, PetscBool isCohesive)
947: {
948: PetscInt i;
950: PetscFunctionBegin;
952: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
953: ds->cohesive[f] = isCohesive;
954: ds->isCohesive = PETSC_FALSE;
955: for (i = 0; i < ds->Nf; ++i) ds->isCohesive = ds->isCohesive || ds->cohesive[f] ? PETSC_TRUE : PETSC_FALSE;
956: PetscFunctionReturn(PETSC_SUCCESS);
957: }
959: /*@
960: PetscDSGetTotalDimension - Returns the total size of the approximation space for this system
962: Not Collective
964: Input Parameter:
965: . prob - The `PetscDS` object
967: Output Parameter:
968: . dim - The total problem dimension
970: Level: beginner
972: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
973: @*/
974: PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim)
975: {
976: PetscFunctionBegin;
978: PetscCall(PetscDSSetUp(prob));
979: PetscAssertPointer(dim, 2);
980: *dim = prob->totDim;
981: PetscFunctionReturn(PETSC_SUCCESS);
982: }
984: /*@
985: PetscDSGetTotalComponents - Returns the total number of components in this system
987: Not Collective
989: Input Parameter:
990: . prob - The `PetscDS` object
992: Output Parameter:
993: . Nc - The total number of components
995: Level: beginner
997: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
998: @*/
999: PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc)
1000: {
1001: PetscFunctionBegin;
1003: PetscCall(PetscDSSetUp(prob));
1004: PetscAssertPointer(Nc, 2);
1005: *Nc = prob->totComp;
1006: PetscFunctionReturn(PETSC_SUCCESS);
1007: }
1009: /*@
1010: PetscDSGetDiscretization - Returns the discretization object for the given field
1012: Not Collective
1014: Input Parameters:
1015: + prob - The `PetscDS` object
1016: - f - The field number
1018: Output Parameter:
1019: . disc - The discretization object, this can be a `PetscFE` or a `PetscFV`
1021: Level: beginner
1023: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1024: @*/
1025: PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc)
1026: {
1027: PetscFunctionBeginHot;
1029: PetscAssertPointer(disc, 3);
1030: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
1031: *disc = prob->disc[f];
1032: PetscFunctionReturn(PETSC_SUCCESS);
1033: }
1035: /*@
1036: PetscDSSetDiscretization - Sets the discretization object for the given field
1038: Not Collective
1040: Input Parameters:
1041: + prob - The `PetscDS` object
1042: . f - The field number
1043: - disc - The discretization object, this can be a `PetscFE` or a `PetscFV`
1045: Level: beginner
1047: .seealso: `PetscDS`, `PetscFE`, `PetscFV`, `PetscDSGetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1048: @*/
1049: PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc)
1050: {
1051: PetscFunctionBegin;
1053: if (disc) PetscAssertPointer(disc, 3);
1054: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1055: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
1056: PetscCall(PetscObjectDereference(prob->disc[f]));
1057: prob->disc[f] = disc;
1058: PetscCall(PetscObjectReference(disc));
1059: if (disc) {
1060: PetscClassId id;
1062: PetscCall(PetscObjectGetClassId(disc, &id));
1063: if (id == PETSCFE_CLASSID) {
1064: PetscCall(PetscDSSetImplicit(prob, f, PETSC_TRUE));
1065: } else if (id == PETSCFV_CLASSID) {
1066: PetscCall(PetscDSSetImplicit(prob, f, PETSC_FALSE));
1067: }
1068: PetscCall(PetscDSSetJetDegree(prob, f, 1));
1069: }
1070: PetscFunctionReturn(PETSC_SUCCESS);
1071: }
1073: /*@
1074: PetscDSGetWeakForm - Returns the weak form object from within the `PetscDS`
1076: Not Collective
1078: Input Parameter:
1079: . ds - The `PetscDS` object
1081: Output Parameter:
1082: . wf - The weak form object
1084: Level: beginner
1086: .seealso: `PetscWeakForm`, `PetscDSSetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1087: @*/
1088: PetscErrorCode PetscDSGetWeakForm(PetscDS ds, PetscWeakForm *wf)
1089: {
1090: PetscFunctionBegin;
1092: PetscAssertPointer(wf, 2);
1093: *wf = ds->wf;
1094: PetscFunctionReturn(PETSC_SUCCESS);
1095: }
1097: /*@
1098: PetscDSSetWeakForm - Sets the weak form object to be used by the `PetscDS`
1100: Not Collective
1102: Input Parameters:
1103: + ds - The `PetscDS` object
1104: - wf - The weak form object
1106: Level: beginner
1108: .seealso: `PetscWeakForm`, `PetscDSGetWeakForm()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1109: @*/
1110: PetscErrorCode PetscDSSetWeakForm(PetscDS ds, PetscWeakForm wf)
1111: {
1112: PetscFunctionBegin;
1115: PetscCall(PetscObjectDereference((PetscObject)ds->wf));
1116: ds->wf = wf;
1117: PetscCall(PetscObjectReference((PetscObject)wf));
1118: PetscCall(PetscWeakFormSetNumFields(wf, ds->Nf));
1119: PetscFunctionReturn(PETSC_SUCCESS);
1120: }
1122: /*@
1123: PetscDSAddDiscretization - Adds a discretization object
1125: Not Collective
1127: Input Parameters:
1128: + prob - The `PetscDS` object
1129: - disc - The discretization object, this can be a `PetscFE` or `PetscFV`
1131: Level: beginner
1133: .seealso: `PetscWeakForm`, `PetscFE`, `PetscFV`, `PetscDSGetDiscretization()`, `PetscDSSetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1134: @*/
1135: PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc)
1136: {
1137: PetscFunctionBegin;
1138: PetscCall(PetscDSSetDiscretization(prob, prob->Nf, disc));
1139: PetscFunctionReturn(PETSC_SUCCESS);
1140: }
1142: /*@
1143: PetscDSGetQuadrature - Returns the quadrature, which must agree for all fields in the `PetscDS`
1145: Not Collective
1147: Input Parameter:
1148: . prob - The `PetscDS` object
1150: Output Parameter:
1151: . q - The quadrature object
1153: Level: intermediate
1155: .seealso: `PetscDS`, `PetscQuadrature`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1156: @*/
1157: PetscErrorCode PetscDSGetQuadrature(PetscDS prob, PetscQuadrature *q)
1158: {
1159: PetscObject obj;
1160: PetscClassId id;
1162: PetscFunctionBegin;
1163: *q = NULL;
1164: if (!prob->Nf) PetscFunctionReturn(PETSC_SUCCESS);
1165: PetscCall(PetscDSGetDiscretization(prob, 0, &obj));
1166: PetscCall(PetscObjectGetClassId(obj, &id));
1167: if (id == PETSCFE_CLASSID) PetscCall(PetscFEGetQuadrature((PetscFE)obj, q));
1168: else if (id == PETSCFV_CLASSID) PetscCall(PetscFVGetQuadrature((PetscFV)obj, q));
1169: else SETERRQ(PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
1170: PetscFunctionReturn(PETSC_SUCCESS);
1171: }
1173: /*@
1174: PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for `TSARKIMEX`
1176: Not Collective
1178: Input Parameters:
1179: + prob - The `PetscDS` object
1180: - f - The field number
1182: Output Parameter:
1183: . implicit - The flag indicating what kind of solve to use for this field
1185: Level: developer
1187: .seealso: `TSARKIMEX`, `PetscDS`, `PetscDSSetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1188: @*/
1189: PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit)
1190: {
1191: PetscFunctionBegin;
1193: PetscAssertPointer(implicit, 3);
1194: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
1195: *implicit = prob->implicit[f];
1196: PetscFunctionReturn(PETSC_SUCCESS);
1197: }
1199: /*@
1200: PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for `TSARKIMEX`
1202: Not Collective
1204: Input Parameters:
1205: + prob - The `PetscDS` object
1206: . f - The field number
1207: - implicit - The flag indicating what kind of solve to use for this field
1209: Level: developer
1211: .seealso: `TSARKIMEX`, `PetscDSGetImplicit()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1212: @*/
1213: PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit)
1214: {
1215: PetscFunctionBegin;
1217: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
1218: prob->implicit[f] = implicit;
1219: PetscFunctionReturn(PETSC_SUCCESS);
1220: }
1222: /*@
1223: PetscDSGetJetDegree - Returns the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate.
1225: Not Collective
1227: Input Parameters:
1228: + ds - The `PetscDS` object
1229: - f - The field number
1231: Output Parameter:
1232: . k - The highest derivative we need to tabulate
1234: Level: developer
1236: .seealso: `PetscDS`, `PetscDSSetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1237: @*/
1238: PetscErrorCode PetscDSGetJetDegree(PetscDS ds, PetscInt f, PetscInt *k)
1239: {
1240: PetscFunctionBegin;
1242: PetscAssertPointer(k, 3);
1243: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1244: *k = ds->jetDegree[f];
1245: PetscFunctionReturn(PETSC_SUCCESS);
1246: }
1248: /*@
1249: PetscDSSetJetDegree - Set the highest derivative for this field equation, or the k-jet that the discretization needs to tabulate.
1251: Not Collective
1253: Input Parameters:
1254: + ds - The `PetscDS` object
1255: . f - The field number
1256: - k - The highest derivative we need to tabulate
1258: Level: developer
1260: .seealso: `PetscDS`, `PetscDSGetJetDegree()`, `PetscDSSetDiscretization()`, `PetscDSAddDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
1261: @*/
1262: PetscErrorCode PetscDSSetJetDegree(PetscDS ds, PetscInt f, PetscInt k)
1263: {
1264: PetscFunctionBegin;
1266: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1267: ds->jetDegree[f] = k;
1268: PetscFunctionReturn(PETSC_SUCCESS);
1269: }
1271: /*@C
1272: PetscDSGetObjective - Get the pointwise objective function for a given test field that was provided with `PetscDSSetObjective()`
1274: Not Collective
1276: Input Parameters:
1277: + ds - The `PetscDS`
1278: - f - The test field number
1280: Output Parameter:
1281: . obj - integrand for the test function term, see `PetscPointFn`
1283: Level: intermediate
1285: Note:
1286: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi\,\mathrm{obj}(u, u_t, \nabla u, x, t)$
1288: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSSetObjective()`, `PetscDSGetResidual()`
1289: @*/
1290: PetscErrorCode PetscDSGetObjective(PetscDS ds, PetscInt f, PetscPointFn **obj)
1291: {
1292: PetscPointFn **tmp;
1293: PetscInt n;
1295: PetscFunctionBegin;
1297: PetscAssertPointer(obj, 3);
1298: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1299: PetscCall(PetscWeakFormGetObjective(ds->wf, NULL, 0, f, 0, &n, &tmp));
1300: *obj = tmp ? tmp[0] : NULL;
1301: PetscFunctionReturn(PETSC_SUCCESS);
1302: }
1304: /*@C
1305: PetscDSSetObjective - Set the pointwise objective function for a given test field
1307: Not Collective
1309: Input Parameters:
1310: + ds - The `PetscDS`
1311: . f - The test field number
1312: - obj - integrand for the test function term, see `PetscPointFn`
1314: Level: intermediate
1316: Note:
1317: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi\,\mathrm{obj}(u, u_t, \nabla u, x, t)$
1319: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSGetObjective()`, `PetscDSSetResidual()`
1320: @*/
1321: PetscErrorCode PetscDSSetObjective(PetscDS ds, PetscInt f, PetscPointFn *obj)
1322: {
1323: PetscFunctionBegin;
1326: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1327: PetscCall(PetscWeakFormSetIndexObjective(ds->wf, NULL, 0, f, 0, 0, obj));
1328: PetscFunctionReturn(PETSC_SUCCESS);
1329: }
1331: /*@C
1332: PetscDSGetResidual - Get the pointwise residual function for a given test field
1334: Not Collective
1336: Input Parameters:
1337: + ds - The `PetscDS`
1338: - f - The test field number
1340: Output Parameters:
1341: + f0 - integrand for the test function term, see `PetscPointFn`
1342: - f1 - integrand for the test function gradient term, see `PetscPointFn`
1344: Level: intermediate
1346: Note:
1347: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t)$
1349: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSSetResidual()`
1350: @*/
1351: PetscErrorCode PetscDSGetResidual(PetscDS ds, PetscInt f, PetscPointFn **f0, PetscPointFn **f1)
1352: {
1353: PetscPointFn **tmp0, **tmp1;
1354: PetscInt n0, n1;
1356: PetscFunctionBegin;
1358: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1359: PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1));
1360: *f0 = tmp0 ? tmp0[0] : NULL;
1361: *f1 = tmp1 ? tmp1[0] : NULL;
1362: PetscFunctionReturn(PETSC_SUCCESS);
1363: }
1365: /*@C
1366: PetscDSSetResidual - Set the pointwise residual function for a given test field
1368: Not Collective
1370: Input Parameters:
1371: + ds - The `PetscDS`
1372: . f - The test field number
1373: . f0 - integrand for the test function term, see `PetscPointFn`
1374: - f1 - integrand for the test function gradient term, see `PetscPointFn`
1376: Level: intermediate
1378: Note:
1379: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t)$
1381: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSGetResidual()`
1382: @*/
1383: PetscErrorCode PetscDSSetResidual(PetscDS ds, PetscInt f, PetscPointFn *f0, PetscPointFn *f1)
1384: {
1385: PetscFunctionBegin;
1389: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1390: PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1));
1391: PetscFunctionReturn(PETSC_SUCCESS);
1392: }
1394: /*@C
1395: PetscDSGetRHSResidual - Get the pointwise RHS residual function for explicit timestepping for a given test field
1397: Not Collective
1399: Input Parameters:
1400: + ds - The `PetscDS`
1401: - f - The test field number
1403: Output Parameters:
1404: + f0 - integrand for the test function term, see `PetscPointFn`
1405: - f1 - integrand for the test function gradient term, see `PetscPointFn`
1407: Level: intermediate
1409: Note:
1410: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t)$
1412: .seealso: `PetscPointFn`, `PetscDS`, `PetscDSSetRHSResidual()`
1413: @*/
1414: PetscErrorCode PetscDSGetRHSResidual(PetscDS ds, PetscInt f, PetscPointFn **f0, PetscPointFn **f1)
1415: {
1416: PetscPointFn **tmp0, **tmp1;
1417: PetscInt n0, n1;
1419: PetscFunctionBegin;
1421: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1422: PetscCall(PetscWeakFormGetResidual(ds->wf, NULL, 0, f, 100, &n0, &tmp0, &n1, &tmp1));
1423: *f0 = tmp0 ? tmp0[0] : NULL;
1424: *f1 = tmp1 ? tmp1[0] : NULL;
1425: PetscFunctionReturn(PETSC_SUCCESS);
1426: }
1428: /*@C
1429: PetscDSSetRHSResidual - Set the pointwise residual function for explicit timestepping for a given test field
1431: Not Collective
1433: Input Parameters:
1434: + ds - The `PetscDS`
1435: . f - The test field number
1436: . f0 - integrand for the test function term, see `PetscPointFn`
1437: - f1 - integrand for the test function gradient term, see `PetscPointFn`
1439: Level: intermediate
1441: Note:
1442: We are using a first order FEM model for the weak form\: $ \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t)$
1444: .seealso: `PetscDS`, `PetscDSGetResidual()`
1445: @*/
1446: PetscErrorCode PetscDSSetRHSResidual(PetscDS ds, PetscInt f, PetscPointFn *f0, PetscPointFn *f1)
1447: {
1448: PetscFunctionBegin;
1452: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1453: PetscCall(PetscWeakFormSetIndexResidual(ds->wf, NULL, 0, f, 100, 0, f0, 0, f1));
1454: PetscFunctionReturn(PETSC_SUCCESS);
1455: }
1457: /*@
1458: PetscDSHasJacobian - Checks that the Jacobian functions have been set
1460: Not Collective
1462: Input Parameter:
1463: . ds - The `PetscDS`
1465: Output Parameter:
1466: . hasJac - flag that indicates the pointwise function for the Jacobian has been set
1468: Level: intermediate
1470: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1471: @*/
1472: PetscErrorCode PetscDSHasJacobian(PetscDS ds, PetscBool *hasJac)
1473: {
1474: PetscFunctionBegin;
1476: PetscCall(PetscWeakFormHasJacobian(ds->wf, hasJac));
1477: PetscFunctionReturn(PETSC_SUCCESS);
1478: }
1480: /*@C
1481: PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field
1483: Not Collective
1485: Input Parameters:
1486: + ds - The `PetscDS`
1487: . f - The test field number
1488: - g - The field number
1490: Output Parameters:
1491: + g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1492: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1493: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1494: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1496: Level: intermediate
1498: Note:
1499: We are using a first order FEM model for the weak form\:
1501: $$
1502: \int_\Omega \phi\, g_0(u, u_t, \nabla u, x, t) \psi + \phi\, {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi
1503: + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1504: $$
1506: .seealso: `PetscDS`, `PetscDSSetJacobian()`, `PetscPointJacFn`
1507: @*/
1508: PetscErrorCode PetscDSGetJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn **g0, PetscPointJacFn **g1, PetscPointJacFn **g2, PetscPointJacFn **g3)
1509: {
1510: PetscPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
1511: PetscInt n0, n1, n2, n3;
1513: PetscFunctionBegin;
1515: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1516: PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf);
1517: PetscCall(PetscWeakFormGetJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1518: *g0 = tmp0 ? tmp0[0] : NULL;
1519: *g1 = tmp1 ? tmp1[0] : NULL;
1520: *g2 = tmp2 ? tmp2[0] : NULL;
1521: *g3 = tmp3 ? tmp3[0] : NULL;
1522: PetscFunctionReturn(PETSC_SUCCESS);
1523: }
1525: /*@C
1526: PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields
1528: Not Collective
1530: Input Parameters:
1531: + ds - The `PetscDS`
1532: . f - The test field number
1533: . g - The field number
1534: . g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1535: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1536: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1537: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1539: Level: intermediate
1541: Note:
1542: We are using a first order FEM model for the weak form\:
1544: $$
1545: \int_\Omega \phi\, g_0(u, u_t, \nabla u, x, t) \psi + \phi\, {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi
1546: + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1547: $$
1549: .seealso: `PetscDS`, `PetscDSGetJacobian()`, `PetscPointJacFn`
1550: @*/
1551: PetscErrorCode PetscDSSetJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn *g0, PetscPointJacFn *g1, PetscPointJacFn *g2, PetscPointJacFn *g3)
1552: {
1553: PetscFunctionBegin;
1559: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1560: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1561: PetscCall(PetscWeakFormSetIndexJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1562: PetscFunctionReturn(PETSC_SUCCESS);
1563: }
1565: /*@
1566: PetscDSUseJacobianPreconditioner - Set whether to construct a Jacobian preconditioner
1568: Not Collective
1570: Input Parameters:
1571: + prob - The `PetscDS`
1572: - useJacPre - flag that enables construction of a Jacobian preconditioner
1574: Level: intermediate
1576: Developer Note:
1577: Should be called `PetscDSSetUseJacobianPreconditioner()`
1579: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1580: @*/
1581: PetscErrorCode PetscDSUseJacobianPreconditioner(PetscDS prob, PetscBool useJacPre)
1582: {
1583: PetscFunctionBegin;
1585: prob->useJacPre = useJacPre;
1586: PetscFunctionReturn(PETSC_SUCCESS);
1587: }
1589: /*@
1590: PetscDSHasJacobianPreconditioner - Checks if a Jacobian matrix for constructing a preconditioner has been set
1592: Not Collective
1594: Input Parameter:
1595: . ds - The `PetscDS`
1597: Output Parameter:
1598: . hasJacPre - the flag
1600: Level: intermediate
1602: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`
1603: @*/
1604: PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS ds, PetscBool *hasJacPre)
1605: {
1606: PetscFunctionBegin;
1608: *hasJacPre = PETSC_FALSE;
1609: if (!ds->useJacPre) PetscFunctionReturn(PETSC_SUCCESS);
1610: PetscCall(PetscWeakFormHasJacobianPreconditioner(ds->wf, hasJacPre));
1611: PetscFunctionReturn(PETSC_SUCCESS);
1612: }
1614: /*@C
1615: PetscDSGetJacobianPreconditioner - Get the pointwise Jacobian function for given test and basis field that constructs the matrix used
1616: to compute the preconditioner. If this is missing, the system matrix is used to build the preconditioner.
1618: Not Collective
1620: Input Parameters:
1621: + ds - The `PetscDS`
1622: . f - The test field number
1623: - g - The field number
1625: Output Parameters:
1626: + g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1627: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1628: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1629: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1631: Level: intermediate
1633: Note:
1634: We are using a first order FEM model for the weak form\:
1636: $$
1637: \int_\Omega \phi\, g_0(u, u_t, \nabla u, x, t) \psi + \phi\, {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi
1638: + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1639: $$
1641: Developer Note:
1642: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
1644: .seealso: `PetscDS`, `PetscDSSetJacobianPreconditioner()`, `PetscDSGetJacobian()`, `PetscPointJacFn`
1645: @*/
1646: PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn **g0, PetscPointJacFn **g1, PetscPointJacFn **g2, PetscPointJacFn **g3)
1647: {
1648: PetscPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
1649: PetscInt n0, n1, n2, n3;
1651: PetscFunctionBegin;
1653: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1654: PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf);
1655: PetscCall(PetscWeakFormGetJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1656: *g0 = tmp0 ? tmp0[0] : NULL;
1657: *g1 = tmp1 ? tmp1[0] : NULL;
1658: *g2 = tmp2 ? tmp2[0] : NULL;
1659: *g3 = tmp3 ? tmp3[0] : NULL;
1660: PetscFunctionReturn(PETSC_SUCCESS);
1661: }
1663: /*@C
1664: PetscDSSetJacobianPreconditioner - Set the pointwise Jacobian function for given test and basis fields that constructs the matrix used
1665: to compute the preconditioner. If this is missing, the system matrix is used to build the preconditioner.
1667: Not Collective
1669: Input Parameters:
1670: + ds - The `PetscDS`
1671: . f - The test field number
1672: . g - The field number
1673: . g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1674: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1675: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1676: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1678: Level: intermediate
1680: Note:
1681: We are using a first order FEM model for the weak form\:
1683: $$
1684: \int_\Omega \phi\, g_0(u, u_t, \nabla u, x, t) \psi + \phi\, {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi
1685: + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1686: $$
1688: Developer Note:
1689: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
1691: .seealso: `PetscDS`, `PetscDSGetJacobianPreconditioner()`, `PetscDSSetJacobian()`, `PetscPointJacFn`
1692: @*/
1693: PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn *g0, PetscPointJacFn *g1, PetscPointJacFn *g2, PetscPointJacFn *g3)
1694: {
1695: PetscFunctionBegin;
1701: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1702: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1703: PetscCall(PetscWeakFormSetIndexJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1704: PetscFunctionReturn(PETSC_SUCCESS);
1705: }
1707: /*@
1708: PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, $dF/du_t$, has been set
1710: Not Collective
1712: Input Parameter:
1713: . ds - The `PetscDS`
1715: Output Parameter:
1716: . hasDynJac - flag that pointwise function for dynamic Jacobian has been set
1718: Level: intermediate
1720: .seealso: `PetscDS`, `PetscDSGetDynamicJacobian()`, `PetscDSSetDynamicJacobian()`, `PetscDSGetJacobian()`
1721: @*/
1722: PetscErrorCode PetscDSHasDynamicJacobian(PetscDS ds, PetscBool *hasDynJac)
1723: {
1724: PetscFunctionBegin;
1726: PetscCall(PetscWeakFormHasDynamicJacobian(ds->wf, hasDynJac));
1727: PetscFunctionReturn(PETSC_SUCCESS);
1728: }
1730: /*@C
1731: PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, $dF/du_t$, function for given test and basis field
1733: Not Collective
1735: Input Parameters:
1736: + ds - The `PetscDS`
1737: . f - The test field number
1738: - g - The field number
1740: Output Parameters:
1741: + g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1742: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1743: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1744: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1746: Level: intermediate
1748: Note:
1749: We are using a first order FEM model for the weak form\:
1751: $$
1752: \int_\Omega \phi\, g_0(u, u_t, \nabla u, x, t) \psi + \phi\, {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi
1753: + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1754: $$
1756: .seealso: `PetscDS`, `PetscDSSetJacobian()`, `PetscDSSetDynamicJacobian()`, `PetscPointJacFn`
1757: @*/
1758: PetscErrorCode PetscDSGetDynamicJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn **g0, PetscPointJacFn **g1, PetscPointJacFn **g2, PetscPointJacFn **g3)
1759: {
1760: PetscPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
1761: PetscInt n0, n1, n2, n3;
1763: PetscFunctionBegin;
1765: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1766: PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf);
1767: PetscCall(PetscWeakFormGetDynamicJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
1768: *g0 = tmp0 ? tmp0[0] : NULL;
1769: *g1 = tmp1 ? tmp1[0] : NULL;
1770: *g2 = tmp2 ? tmp2[0] : NULL;
1771: *g3 = tmp3 ? tmp3[0] : NULL;
1772: PetscFunctionReturn(PETSC_SUCCESS);
1773: }
1775: /*@C
1776: PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, $dF/du_t$, function for given test and basis fields
1778: Not Collective
1780: Input Parameters:
1781: + ds - The `PetscDS`
1782: . f - The test field number
1783: . g - The field number
1784: . g0 - integrand for the test and basis function term, see `PetscPointJacFn`
1785: . g1 - integrand for the test function and basis function gradient term, see `PetscPointJacFn`
1786: . g2 - integrand for the test function gradient and basis function term, see `PetscPointJacFn`
1787: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscPointJacFn`
1789: Level: intermediate
1791: Note:
1792: We are using a first order FEM model for the weak form\:
1794: $$
1795: \int_\Omega \phi\, g_0(u, u_t, \nabla u, x, t) \psi + \phi\, {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi
1796: + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1797: $$
1799: .seealso: `PetscDS`, `PetscDSGetDynamicJacobian()`, `PetscDSGetJacobian()`, `PetscPointJacFn`
1800: @*/
1801: PetscErrorCode PetscDSSetDynamicJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscPointJacFn *g0, PetscPointJacFn *g1, PetscPointJacFn *g2, PetscPointJacFn *g3)
1802: {
1803: PetscFunctionBegin;
1809: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1810: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
1811: PetscCall(PetscWeakFormSetIndexDynamicJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
1812: PetscFunctionReturn(PETSC_SUCCESS);
1813: }
1815: /*@C
1816: PetscDSGetRiemannSolver - Returns the Riemann solver for the given field
1818: Not Collective
1820: Input Parameters:
1821: + ds - The `PetscDS` object
1822: - f - The field number
1824: Output Parameter:
1825: . r - Riemann solver, see `PetscRiemannFn`
1827: Level: intermediate
1829: .seealso: `PetscDS`, `PetscRiemannFn`, `PetscDSSetRiemannSolver()`
1830: @*/
1831: PetscErrorCode PetscDSGetRiemannSolver(PetscDS ds, PetscInt f, PetscRiemannFn **r)
1832: {
1833: PetscRiemannFn **tmp;
1834: PetscInt n;
1836: PetscFunctionBegin;
1838: PetscAssertPointer(r, 3);
1839: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1840: PetscCall(PetscWeakFormGetRiemannSolver(ds->wf, NULL, 0, f, 0, &n, &tmp));
1841: *r = tmp ? tmp[0] : NULL;
1842: PetscFunctionReturn(PETSC_SUCCESS);
1843: }
1845: /*@C
1846: PetscDSSetRiemannSolver - Sets the Riemann solver for the given field
1848: Not Collective
1850: Input Parameters:
1851: + ds - The `PetscDS` object
1852: . f - The field number
1853: - r - Riemann solver, see `PetscRiemannFn`
1855: Level: intermediate
1857: .seealso: `PetscDS`, `PetscRiemannFn`, `PetscDSGetRiemannSolver()`
1858: @*/
1859: PetscErrorCode PetscDSSetRiemannSolver(PetscDS ds, PetscInt f, PetscRiemannFn *r)
1860: {
1861: PetscFunctionBegin;
1864: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1865: PetscCall(PetscWeakFormSetIndexRiemannSolver(ds->wf, NULL, 0, f, 0, 0, r));
1866: PetscFunctionReturn(PETSC_SUCCESS);
1867: }
1869: /*@C
1870: PetscDSGetUpdate - Get the pointwise update function for a given field
1872: Not Collective
1874: Input Parameters:
1875: + ds - The `PetscDS`
1876: - f - The field number
1878: Output Parameter:
1879: . update - update function, see `PetscPointFn`
1881: Level: intermediate
1883: .seealso: `PetscDS`, `PetscPointFn`, `PetscDSSetUpdate()`, `PetscDSSetResidual()`
1884: @*/
1885: PetscErrorCode PetscDSGetUpdate(PetscDS ds, PetscInt f, PetscPointFn **update)
1886: {
1887: PetscFunctionBegin;
1889: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1890: if (update) {
1891: PetscAssertPointer(update, 3);
1892: *update = ds->update[f];
1893: }
1894: PetscFunctionReturn(PETSC_SUCCESS);
1895: }
1897: /*@C
1898: PetscDSSetUpdate - Set the pointwise update function for a given field
1900: Not Collective
1902: Input Parameters:
1903: + ds - The `PetscDS`
1904: . f - The field number
1905: - update - update function, see `PetscPointFn`
1907: Level: intermediate
1909: .seealso: `PetscDS`, `PetscPointFn`, `PetscDSGetResidual()`
1910: @*/
1911: PetscErrorCode PetscDSSetUpdate(PetscDS ds, PetscInt f, PetscPointFn *update)
1912: {
1913: PetscFunctionBegin;
1916: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1917: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
1918: ds->update[f] = update;
1919: PetscFunctionReturn(PETSC_SUCCESS);
1920: }
1922: /*@C
1923: PetscDSGetContext - Returns the context that was passed by `PetscDSSetContext()`
1925: Not Collective
1927: Input Parameters:
1928: + ds - The `PetscDS`
1929: . f - The field number
1930: - ctx - the context
1932: Level: intermediate
1934: Fortran Notes:
1935: This only works when the context is a Fortran derived type or a `PetscObject`. Define `ctx` with
1936: .vb
1937: type(tUsertype), pointer :: ctx
1938: .ve
1940: .seealso: `PetscDS`, `PetscPointFn`, `PetscDSSetContext()`
1941: @*/
1942: PetscErrorCode PetscDSGetContext(PetscDS ds, PetscInt f, PetscCtxRt ctx)
1943: {
1944: PetscFunctionBegin;
1946: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
1947: PetscAssertPointer(ctx, 3);
1948: *(void **)ctx = ds->ctx[f];
1949: PetscFunctionReturn(PETSC_SUCCESS);
1950: }
1952: /*@C
1953: PetscDSSetContext - Sets the context that is passed back to some of the pointwise function callbacks used by this `PetscDS`
1955: Not Collective
1957: Input Parameters:
1958: + ds - The `PetscDS`
1959: . f - The field number
1960: - ctx - the context
1962: Level: intermediate
1964: .seealso: `PetscDS`, `PetscPointFn`, `PetscDSGetContext()`
1965: @*/
1966: PetscErrorCode PetscDSSetContext(PetscDS ds, PetscInt f, PetscCtx ctx)
1967: {
1968: PetscFunctionBegin;
1970: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
1971: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
1972: ds->ctx[f] = ctx;
1973: PetscFunctionReturn(PETSC_SUCCESS);
1974: }
1976: /*@C
1977: PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field
1979: Not Collective
1981: Input Parameters:
1982: + ds - The PetscDS
1983: - f - The test field number
1985: Output Parameters:
1986: + f0 - boundary integrand for the test function term, see `PetscBdPointFn`
1987: - f1 - boundary integrand for the test function gradient term, see `PetscBdPointFn`
1989: Level: intermediate
1991: Note:
1992: We are using a first order FEM model for the weak form\:
1994: $$
1995: \int_\Gamma \phi {\vec f}_0(u, u_t, \nabla u, x, t) \cdot \hat n + \nabla\phi \cdot {\overleftrightarrow f}_1(u, u_t, \nabla u, x, t) \cdot \hat n
1996: $$
1998: .seealso: `PetscDS`, `PetscBdPointFn`, `PetscDSSetBdResidual()`
1999: @*/
2000: PetscErrorCode PetscDSGetBdResidual(PetscDS ds, PetscInt f, PetscBdPointFn **f0, PetscBdPointFn **f1)
2001: {
2002: PetscBdPointFn **tmp0, **tmp1;
2003: PetscInt n0, n1;
2005: PetscFunctionBegin;
2007: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
2008: PetscCall(PetscWeakFormGetBdResidual(ds->wf, NULL, 0, f, 0, &n0, &tmp0, &n1, &tmp1));
2009: *f0 = tmp0 ? tmp0[0] : NULL;
2010: *f1 = tmp1 ? tmp1[0] : NULL;
2011: PetscFunctionReturn(PETSC_SUCCESS);
2012: }
2014: /*@C
2015: PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field
2017: Not Collective
2019: Input Parameters:
2020: + ds - The `PetscDS`
2021: . f - The test field number
2022: . f0 - boundary integrand for the test function term, see `PetscBdPointFn`
2023: - f1 - boundary integrand for the test function gradient term, see `PetscBdPointFn`
2025: Level: intermediate
2027: Note:
2028: We are using a first order FEM model for the weak form\:
2030: $$
2031: \int_\Gamma \phi {\vec f}_0(u, u_t, \nabla u, x, t) \cdot \hat n + \nabla\phi \cdot {\overleftrightarrow f}_1(u, u_t, \nabla u, x, t) \cdot \hat n
2032: $$
2034: .seealso: `PetscDS`, `PetscBdPointFn`, `PetscDSGetBdResidual()`
2035: @*/
2036: PetscErrorCode PetscDSSetBdResidual(PetscDS ds, PetscInt f, PetscBdPointFn *f0, PetscBdPointFn *f1)
2037: {
2038: PetscFunctionBegin;
2040: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2041: PetscCall(PetscWeakFormSetIndexBdResidual(ds->wf, NULL, 0, f, 0, 0, f0, 0, f1));
2042: PetscFunctionReturn(PETSC_SUCCESS);
2043: }
2045: /*@
2046: PetscDSHasBdJacobian - Indicates that boundary Jacobian functions have been set
2048: Not Collective
2050: Input Parameter:
2051: . ds - The `PetscDS`
2053: Output Parameter:
2054: . hasBdJac - flag that pointwise function for the boundary Jacobian has been set
2056: Level: intermediate
2058: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2059: @*/
2060: PetscErrorCode PetscDSHasBdJacobian(PetscDS ds, PetscBool *hasBdJac)
2061: {
2062: PetscFunctionBegin;
2064: PetscAssertPointer(hasBdJac, 2);
2065: PetscCall(PetscWeakFormHasBdJacobian(ds->wf, hasBdJac));
2066: PetscFunctionReturn(PETSC_SUCCESS);
2067: }
2069: /*@C
2070: PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field
2072: Not Collective
2074: Input Parameters:
2075: + ds - The `PetscDS`
2076: . f - The test field number
2077: - g - The field number
2079: Output Parameters:
2080: + g0 - integrand for the test and basis function term, see `PetscBdPointJacFn`
2081: . g1 - integrand for the test function and basis function gradient term, see `PetscBdPointJacFn`
2082: . g2 - integrand for the test function gradient and basis function term, see `PetscBdPointJacFn`
2083: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscBdPointJacFn`
2085: Level: intermediate
2087: Note:
2088: We are using a first order FEM model for the weak form\:
2090: $$
2091: \int_\Gamma \phi\, {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi\, {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi
2092: + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi
2093: $$
2095: .seealso: `PetscDS`, `PetscBdPointJacFn`, `PetscDSSetBdJacobian()`
2096: @*/
2097: PetscErrorCode PetscDSGetBdJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscBdPointJacFn **g0, PetscBdPointJacFn **g1, PetscBdPointJacFn **g2, PetscBdPointJacFn **g3)
2098: {
2099: PetscBdPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
2100: PetscInt n0, n1, n2, n3;
2102: PetscFunctionBegin;
2104: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
2105: PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf);
2106: PetscCall(PetscWeakFormGetBdJacobian(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2107: *g0 = tmp0 ? tmp0[0] : NULL;
2108: *g1 = tmp1 ? tmp1[0] : NULL;
2109: *g2 = tmp2 ? tmp2[0] : NULL;
2110: *g3 = tmp3 ? tmp3[0] : NULL;
2111: PetscFunctionReturn(PETSC_SUCCESS);
2112: }
2114: /*@C
2115: PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field
2117: Not Collective
2119: Input Parameters:
2120: + ds - The PetscDS
2121: . f - The test field number
2122: . g - The field number
2123: . g0 - integrand for the test and basis function term, see `PetscBdPointJacFn`
2124: . g1 - integrand for the test function and basis function gradient term, see `PetscBdPointJacFn`
2125: . g2 - integrand for the test function gradient and basis function term, see `PetscBdPointJacFn`
2126: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscBdPointJacFn`
2128: Level: intermediate
2130: Note:
2131: We are using a first order FEM model for the weak form\:
2133: $$
2134: \int_\Gamma \phi\, {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi\, {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi
2135: + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi
2136: $$
2138: .seealso: `PetscDS`, `PetscBdPointJacFn`, `PetscDSGetBdJacobian()`
2139: @*/
2140: PetscErrorCode PetscDSSetBdJacobian(PetscDS ds, PetscInt f, PetscInt g, PetscBdPointJacFn *g0, PetscBdPointJacFn *g1, PetscBdPointJacFn *g2, PetscBdPointJacFn *g3)
2141: {
2142: PetscFunctionBegin;
2148: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2149: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2150: PetscCall(PetscWeakFormSetIndexBdJacobian(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2151: PetscFunctionReturn(PETSC_SUCCESS);
2152: }
2154: /*@
2155: PetscDSHasBdJacobianPreconditioner - Signals that boundary Jacobian preconditioner functions have been set with `PetscDSSetBdJacobianPreconditioner()`
2157: Not Collective
2159: Input Parameter:
2160: . ds - The `PetscDS`
2162: Output Parameter:
2163: . hasBdJacPre - flag that pointwise function for the boundary Jacobian matrix to construct the preconditioner has been set
2165: Level: intermediate
2167: Developer Note:
2168: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
2170: .seealso: `PetscDS`, `PetscDSHasJacobian()`, `PetscDSSetBdJacobian()`, `PetscDSGetBdJacobian()`
2171: @*/
2172: PetscErrorCode PetscDSHasBdJacobianPreconditioner(PetscDS ds, PetscBool *hasBdJacPre)
2173: {
2174: PetscFunctionBegin;
2176: PetscAssertPointer(hasBdJacPre, 2);
2177: PetscCall(PetscWeakFormHasBdJacobianPreconditioner(ds->wf, hasBdJacPre));
2178: PetscFunctionReturn(PETSC_SUCCESS);
2179: }
2181: /*@C
2182: PetscDSGetBdJacobianPreconditioner - Get the pointwise boundary Jacobian function for given test and basis field that constructs the
2183: matrix used to construct the preconditioner
2185: Not Collective; No Fortran Support
2187: Input Parameters:
2188: + ds - The `PetscDS`
2189: . f - The test field number
2190: - g - The field number
2192: Output Parameters:
2193: + g0 - integrand for the test and basis function term, see `PetscBdPointJacFn`
2194: . g1 - integrand for the test function and basis function gradient term, see `PetscBdPointJacFn`
2195: . g2 - integrand for the test function gradient and basis function term, see `PetscBdPointJacFn`
2196: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscBdPointJacFn`
2198: Level: intermediate
2200: Note:
2201: We are using a first order FEM model for the weak form\:
2203: $$
2204: \int_\Gamma \phi\, {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi\, {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi
2205: + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi
2206: $$
2208: Developer Note:
2209: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
2211: .seealso: `PetscDS`, `PetscBdPointJacFn`, `PetscDSSetBdJacobianPreconditioner()`
2212: @*/
2213: PetscErrorCode PetscDSGetBdJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, PetscBdPointJacFn **g0, PetscBdPointJacFn **g1, PetscBdPointJacFn **g2, PetscBdPointJacFn **g3)
2214: {
2215: PetscBdPointJacFn **tmp0, **tmp1, **tmp2, **tmp3;
2216: PetscInt n0, n1, n2, n3;
2218: PetscFunctionBegin;
2220: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
2221: PetscCheck(!(g < 0) && !(g >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", g, ds->Nf);
2222: PetscCall(PetscWeakFormGetBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, &n0, &tmp0, &n1, &tmp1, &n2, &tmp2, &n3, &tmp3));
2223: *g0 = tmp0 ? tmp0[0] : NULL;
2224: *g1 = tmp1 ? tmp1[0] : NULL;
2225: *g2 = tmp2 ? tmp2[0] : NULL;
2226: *g3 = tmp3 ? tmp3[0] : NULL;
2227: PetscFunctionReturn(PETSC_SUCCESS);
2228: }
2230: /*@C
2231: PetscDSSetBdJacobianPreconditioner - Set the pointwise boundary Jacobian preconditioner function for given test and basis field that constructs the
2232: matrix used to construct the preconditioner
2234: Not Collective; No Fortran Support
2236: Input Parameters:
2237: + ds - The `PetscDS`
2238: . f - The test field number
2239: . g - The field number
2240: . g0 - integrand for the test and basis function term, see `PetscBdPointJacFn`
2241: . g1 - integrand for the test function and basis function gradient term, see `PetscBdPointJacFn`
2242: . g2 - integrand for the test function gradient and basis function term, see `PetscBdPointJacFn`
2243: - g3 - integrand for the test function gradient and basis function gradient term, see `PetscBdPointJacFn`
2245: Level: intermediate
2247: Note:
2248: We are using a first order FEM model for the weak form\:
2250: $$
2251: \int_\Gamma \phi\, {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi\, {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi
2252: + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi
2253: $$
2255: Developer Note:
2256: The name is confusing since the function computes a matrix used to construct the preconditioner, not a preconditioner.
2258: .seealso: `PetscDS`, `PetscBdPointJacFn`, `PetscDSGetBdJacobianPreconditioner()`
2259: @*/
2260: PetscErrorCode PetscDSSetBdJacobianPreconditioner(PetscDS ds, PetscInt f, PetscInt g, PetscBdPointJacFn *g0, PetscBdPointJacFn *g1, PetscBdPointJacFn *g2, PetscBdPointJacFn *g3)
2261: {
2262: PetscFunctionBegin;
2268: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2269: PetscCheck(g >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", g);
2270: PetscCall(PetscWeakFormSetIndexBdJacobianPreconditioner(ds->wf, NULL, 0, f, g, 0, 0, g0, 0, g1, 0, g2, 0, g3));
2271: PetscFunctionReturn(PETSC_SUCCESS);
2272: }
2274: /*@C
2275: PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field
2277: Not Collective
2279: Input Parameters:
2280: + prob - The `PetscDS`
2281: - f - The test field number
2283: Output Parameters:
2284: + sol - exact solution function for the test field, see `PetscPointExactSolutionFn`
2285: - ctx - exact solution context
2287: Level: intermediate
2289: .seealso: `PetscDS`, `PetscPointExactSolutionFn`, `PetscDSSetExactSolution()`, `PetscDSGetExactSolutionTimeDerivative()`
2290: @*/
2291: PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscPointExactSolutionFn **sol, void **ctx)
2292: {
2293: PetscFunctionBegin;
2295: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
2296: if (sol) {
2297: PetscAssertPointer(sol, 3);
2298: *sol = prob->exactSol[f];
2299: }
2300: if (ctx) {
2301: PetscAssertPointer(ctx, 4);
2302: *ctx = prob->exactCtx[f];
2303: }
2304: PetscFunctionReturn(PETSC_SUCCESS);
2305: }
2307: /*@C
2308: PetscDSSetExactSolution - Set the pointwise exact solution function for a given test field
2310: Not Collective
2312: Input Parameters:
2313: + prob - The `PetscDS`
2314: . f - The test field number
2315: . sol - solution function for the test fields, see `PetscPointExactSolutionFn`
2316: - ctx - solution context or `NULL`
2318: Level: intermediate
2320: .seealso: `PetscDS`, `PetscPointExactSolutionFn`, `PetscDSGetExactSolution()`
2321: @*/
2322: PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscPointExactSolutionFn *sol, PetscCtx ctx)
2323: {
2324: PetscFunctionBegin;
2326: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2327: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
2328: if (sol) {
2330: prob->exactSol[f] = sol;
2331: }
2332: if (ctx) {
2334: prob->exactCtx[f] = ctx;
2335: }
2336: PetscFunctionReturn(PETSC_SUCCESS);
2337: }
2339: /*@C
2340: PetscDSGetExactSolutionTimeDerivative - Get the pointwise time derivative of the exact solution function for a given test field
2342: Not Collective
2344: Input Parameters:
2345: + prob - The `PetscDS`
2346: - f - The test field number
2348: Output Parameters:
2349: + sol - time derivative of the exact solution for the test field, see `PetscPointExactSolutionFn`
2350: - ctx - the exact solution context
2352: Level: intermediate
2354: .seealso: `PetscDS`, `PetscPointExactSolutionFn`, `PetscDSSetExactSolutionTimeDerivative()`, `PetscDSGetExactSolution()`
2355: @*/
2356: PetscErrorCode PetscDSGetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscPointExactSolutionFn **sol, void **ctx)
2357: {
2358: PetscFunctionBegin;
2360: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
2361: if (sol) {
2362: PetscAssertPointer(sol, 3);
2363: *sol = prob->exactSol_t[f];
2364: }
2365: if (ctx) {
2366: PetscAssertPointer(ctx, 4);
2367: *ctx = prob->exactCtx_t[f];
2368: }
2369: PetscFunctionReturn(PETSC_SUCCESS);
2370: }
2372: /*@C
2373: PetscDSSetExactSolutionTimeDerivative - Set the pointwise time derivative of the exact solution function for a given test field
2375: Not Collective
2377: Input Parameters:
2378: + prob - The `PetscDS`
2379: . f - The test field number
2380: . sol - time derivative of the solution function for the test fields, see `PetscPointExactSolutionFn`
2381: - ctx - the solution context or `NULL`
2383: Level: intermediate
2385: .seealso: `PetscDS`, `PetscPointExactSolutionFn`, `PetscDSGetExactSolutionTimeDerivative()`, `PetscDSSetExactSolution()`
2386: @*/
2387: PetscErrorCode PetscDSSetExactSolutionTimeDerivative(PetscDS prob, PetscInt f, PetscPointExactSolutionFn *sol, PetscCtx ctx)
2388: {
2389: PetscFunctionBegin;
2391: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2392: PetscCall(PetscDSEnlarge_Static(prob, f + 1));
2393: if (sol) {
2395: prob->exactSol_t[f] = sol;
2396: }
2397: if (ctx) {
2399: prob->exactCtx_t[f] = ctx;
2400: }
2401: PetscFunctionReturn(PETSC_SUCCESS);
2402: }
2404: /*@C
2405: PetscDSGetLowerBound - Get the pointwise lower bound function for a given field
2407: Not Collective
2409: Input Parameters:
2410: + ds - The PetscDS
2411: - f - The field number
2413: Output Parameters:
2414: + lb - lower bound function for the field, see `PetscPointBoundFn`
2415: - ctx - lower bound context that was set with `PetscDSSetLowerBound()`
2417: Level: intermediate
2419: .seealso: `PetscDS`, `PetscPointBoundFn`, `PetscDSSetLowerBound()`, `PetscDSGetUpperBound()`, `PetscDSGetExactSolution()`
2420: @*/
2421: PetscErrorCode PetscDSGetLowerBound(PetscDS ds, PetscInt f, PetscPointBoundFn **lb, void **ctx)
2422: {
2423: PetscFunctionBegin;
2425: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
2426: if (lb) {
2427: PetscAssertPointer(lb, 3);
2428: *lb = ds->lowerBound[f];
2429: }
2430: if (ctx) {
2431: PetscAssertPointer(ctx, 4);
2432: *ctx = ds->lowerCtx[f];
2433: }
2434: PetscFunctionReturn(PETSC_SUCCESS);
2435: }
2437: /*@C
2438: PetscDSSetLowerBound - Set the pointwise lower bound function for a given field
2440: Not Collective
2442: Input Parameters:
2443: + ds - The `PetscDS`
2444: . f - The field number
2445: . lb - lower bound function for the test fields, see `PetscPointBoundFn`
2446: - ctx - lower bound context or `NULL` which will be passed to `lb`
2448: Level: intermediate
2450: .seealso: `PetscDS`, `PetscPointBoundFn`, `PetscDSGetLowerBound()`, `PetscDSGetUpperBound()`, `PetscDSGetExactSolution()`
2451: @*/
2452: PetscErrorCode PetscDSSetLowerBound(PetscDS ds, PetscInt f, PetscPointBoundFn *lb, PetscCtx ctx)
2453: {
2454: PetscFunctionBegin;
2456: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2457: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
2458: if (lb) {
2460: ds->lowerBound[f] = lb;
2461: }
2462: if (ctx) {
2464: ds->lowerCtx[f] = ctx;
2465: }
2466: PetscFunctionReturn(PETSC_SUCCESS);
2467: }
2469: /*@C
2470: PetscDSGetUpperBound - Get the pointwise upper bound function for a given field
2472: Not Collective
2474: Input Parameters:
2475: + ds - The `PetscDS`
2476: - f - The field number
2478: Output Parameters:
2479: + ub - upper bound function for the field, see `PetscPointBoundFn`
2480: - ctx - upper bound context that was set with `PetscDSSetUpperBound()`
2482: Level: intermediate
2484: .seealso: `PetscDS`, `PetscPointBoundFn`, `PetscDSSetUpperBound()`, `PetscDSGetLowerBound()`, `PetscDSGetExactSolution()`
2485: @*/
2486: PetscErrorCode PetscDSGetUpperBound(PetscDS ds, PetscInt f, PetscPointBoundFn **ub, void **ctx)
2487: {
2488: PetscFunctionBegin;
2490: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
2491: if (ub) {
2492: PetscAssertPointer(ub, 3);
2493: *ub = ds->upperBound[f];
2494: }
2495: if (ctx) {
2496: PetscAssertPointer(ctx, 4);
2497: *ctx = ds->upperCtx[f];
2498: }
2499: PetscFunctionReturn(PETSC_SUCCESS);
2500: }
2502: /*@C
2503: PetscDSSetUpperBound - Set the pointwise upper bound function for a given field
2505: Not Collective
2507: Input Parameters:
2508: + ds - The `PetscDS`
2509: . f - The field number
2510: . ub - upper bound function for the test fields, see `PetscPointBoundFn`
2511: - ctx - context or `NULL` that will be passed to `ub`
2513: Level: intermediate
2515: .seealso: `PetscDS`, `PetscPointBoundFn`, `PetscDSGetUpperBound()`, `PetscDSGetLowerBound()`, `PetscDSGetExactSolution()`
2516: @*/
2517: PetscErrorCode PetscDSSetUpperBound(PetscDS ds, PetscInt f, PetscPointBoundFn *ub, PetscCtx ctx)
2518: {
2519: PetscFunctionBegin;
2521: PetscCheck(f >= 0, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be non-negative", f);
2522: PetscCall(PetscDSEnlarge_Static(ds, f + 1));
2523: if (ub) {
2525: ds->upperBound[f] = ub;
2526: }
2527: if (ctx) {
2529: ds->upperCtx[f] = ctx;
2530: }
2531: PetscFunctionReturn(PETSC_SUCCESS);
2532: }
2534: /*@C
2535: PetscDSGetConstants - Returns the array of constants passed to point functions from a `PetscDS` object
2537: Not Collective
2539: Input Parameter:
2540: . ds - The `PetscDS` object
2542: Output Parameters:
2543: + numConstants - The number of constants, or pass in `NULL` if not required
2544: - constants - The array of constants, `NULL` if there are none
2546: Level: intermediate
2548: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSCreate()`
2549: @*/
2550: PetscErrorCode PetscDSGetConstants(PetscDS ds, PeOp PetscInt *numConstants, PeOp const PetscScalar *constants[])
2551: {
2552: PetscFunctionBegin;
2554: if (numConstants) {
2555: PetscAssertPointer(numConstants, 2);
2556: *numConstants = ds->numConstants;
2557: }
2558: if (constants) {
2559: PetscAssertPointer(constants, 3);
2560: *constants = ds->constants;
2561: }
2562: PetscFunctionReturn(PETSC_SUCCESS);
2563: }
2565: /*@C
2566: PetscDSSetConstants - Set the array of constants passed to point functions from a `PetscDS`
2568: Not Collective
2570: Input Parameters:
2571: + ds - The `PetscDS` object
2572: . numConstants - The number of constants
2573: - constants - The array of constants, `NULL` if there are none
2575: Level: intermediate
2577: .seealso: `PetscDS`, `PetscDSGetConstants()`, `PetscDSCreate()`
2578: @*/
2579: PetscErrorCode PetscDSSetConstants(PetscDS ds, PetscInt numConstants, PetscScalar constants[])
2580: {
2581: PetscFunctionBegin;
2583: if (numConstants != ds->numConstants) {
2584: PetscCall(PetscFree(ds->constants));
2585: ds->numConstants = numConstants;
2586: PetscCall(PetscMalloc1(ds->numConstants + ds->numFuncConstants, &ds->constants));
2587: }
2588: if (ds->numConstants) {
2589: PetscAssertPointer(constants, 3);
2590: PetscCall(PetscArraycpy(ds->constants, constants, ds->numConstants));
2591: }
2592: PetscFunctionReturn(PETSC_SUCCESS);
2593: }
2595: /*@C
2596: PetscDSSetIntegrationParameters - Set the parameters for a particular integration
2598: Not Collective
2600: Input Parameters:
2601: + ds - The `PetscDS` object
2602: . fieldI - The test field for a given point function, or `PETSC_DETERMINE`
2603: - fieldJ - The basis field for a given point function, or `PETSC_DETERMINE`
2605: Level: intermediate
2607: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSGetConstants()`, `PetscDSCreate()`
2608: @*/
2609: PetscErrorCode PetscDSSetIntegrationParameters(PetscDS ds, PetscInt fieldI, PetscInt fieldJ)
2610: {
2611: PetscFunctionBegin;
2613: ds->constants[ds->numConstants] = fieldI;
2614: ds->constants[ds->numConstants + 1] = fieldJ;
2615: PetscFunctionReturn(PETSC_SUCCESS);
2616: }
2618: /*@C
2619: PetscDSSetCellParameters - Set the parameters for a particular cell
2621: Not Collective
2623: Input Parameters:
2624: + ds - The `PetscDS` object
2625: - volume - The cell volume
2627: Level: intermediate
2629: .seealso: `PetscDS`, `PetscDSSetConstants()`, `PetscDSGetConstants()`, `PetscDSCreate()`
2630: @*/
2631: PetscErrorCode PetscDSSetCellParameters(PetscDS ds, PetscReal volume)
2632: {
2633: PetscFunctionBegin;
2635: ds->constants[ds->numConstants + 2] = volume;
2636: PetscFunctionReturn(PETSC_SUCCESS);
2637: }
2639: /*@
2640: PetscDSGetFieldIndex - Returns the index of the given field
2642: Not Collective
2644: Input Parameters:
2645: + prob - The `PetscDS` object
2646: - disc - The discretization object
2648: Output Parameter:
2649: . f - The field number
2651: Level: beginner
2653: .seealso: `PetscDS`, `PetscGetDiscretization()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2654: @*/
2655: PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f)
2656: {
2657: PetscInt g;
2659: PetscFunctionBegin;
2661: PetscAssertPointer(f, 3);
2662: *f = -1;
2663: for (g = 0; g < prob->Nf; ++g) {
2664: if (disc == prob->disc[g]) break;
2665: }
2666: PetscCheck(g != prob->Nf, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS.");
2667: *f = g;
2668: PetscFunctionReturn(PETSC_SUCCESS);
2669: }
2671: /*@
2672: PetscDSGetFieldSize - Returns the size of the given field in the full space basis
2674: Not Collective
2676: Input Parameters:
2677: + prob - The `PetscDS` object
2678: - f - The field number
2680: Output Parameter:
2681: . size - The size
2683: Level: beginner
2685: .seealso: `PetscDS`, `PetscDSGetFieldOffset()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2686: @*/
2687: PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size)
2688: {
2689: PetscFunctionBegin;
2691: PetscAssertPointer(size, 3);
2692: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
2693: PetscCall(PetscDSSetUp(prob));
2694: *size = prob->Nb[f];
2695: PetscFunctionReturn(PETSC_SUCCESS);
2696: }
2698: /*@
2699: PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis
2701: Not Collective
2703: Input Parameters:
2704: + prob - The `PetscDS` object
2705: - f - The field number
2707: Output Parameter:
2708: . off - The offset
2710: Level: beginner
2712: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2713: @*/
2714: PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off)
2715: {
2716: PetscInt size, g;
2718: PetscFunctionBegin;
2720: PetscAssertPointer(off, 3);
2721: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
2722: *off = 0;
2723: for (g = 0; g < f; ++g) {
2724: PetscCall(PetscDSGetFieldSize(prob, g, &size));
2725: *off += size;
2726: }
2727: PetscFunctionReturn(PETSC_SUCCESS);
2728: }
2730: /*@
2731: PetscDSGetFieldOffsetCohesive - Returns the offset of the given field in the full space basis on a cohesive cell
2733: Not Collective
2735: Input Parameters:
2736: + ds - The `PetscDS` object
2737: - f - The field number
2739: Output Parameter:
2740: . off - The offset
2742: Level: beginner
2744: .seealso: `PetscDS`, `PetscDSGetFieldSize()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2745: @*/
2746: PetscErrorCode PetscDSGetFieldOffsetCohesive(PetscDS ds, PetscInt f, PetscInt *off)
2747: {
2748: PetscInt size, g;
2750: PetscFunctionBegin;
2752: PetscAssertPointer(off, 3);
2753: PetscCheck(!(f < 0) && !(f >= ds->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, ds->Nf);
2754: *off = 0;
2755: for (g = 0; g < f; ++g) {
2756: PetscBool cohesive;
2758: PetscCall(PetscDSGetCohesive(ds, g, &cohesive));
2759: PetscCall(PetscDSGetFieldSize(ds, g, &size));
2760: *off += cohesive ? size : size * 2;
2761: }
2762: PetscFunctionReturn(PETSC_SUCCESS);
2763: }
2765: /*@
2766: PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point
2768: Not Collective
2770: Input Parameter:
2771: . prob - The `PetscDS` object
2773: Output Parameter:
2774: . dimensions - The number of dimensions
2776: Level: beginner
2778: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2779: @*/
2780: PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[])
2781: {
2782: PetscFunctionBegin;
2784: PetscCall(PetscDSSetUp(prob));
2785: PetscAssertPointer(dimensions, 2);
2786: *dimensions = prob->Nb;
2787: PetscFunctionReturn(PETSC_SUCCESS);
2788: }
2790: /*@
2791: PetscDSGetComponents - Returns the number of components for each field on an evaluation point
2793: Not Collective
2795: Input Parameter:
2796: . prob - The `PetscDS` object
2798: Output Parameter:
2799: . components - The number of components
2801: Level: beginner
2803: .seealso: `PetscDS`, `PetscDSGetComponentOffsets()`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2804: @*/
2805: PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[])
2806: {
2807: PetscFunctionBegin;
2809: PetscCall(PetscDSSetUp(prob));
2810: PetscAssertPointer(components, 2);
2811: *components = prob->Nc;
2812: PetscFunctionReturn(PETSC_SUCCESS);
2813: }
2815: /*@
2816: PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point
2818: Not Collective
2820: Input Parameters:
2821: + prob - The `PetscDS` object
2822: - f - The field number
2824: Output Parameter:
2825: . off - The offset
2827: Level: beginner
2829: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2830: @*/
2831: PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off)
2832: {
2833: PetscFunctionBegin;
2835: PetscAssertPointer(off, 3);
2836: PetscCheck(!(f < 0) && !(f >= prob->Nf), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, prob->Nf);
2837: PetscCall(PetscDSSetUp(prob));
2838: *off = prob->off[f];
2839: PetscFunctionReturn(PETSC_SUCCESS);
2840: }
2842: /*@
2843: PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point
2845: Not Collective
2847: Input Parameter:
2848: . prob - The `PetscDS` object
2850: Output Parameter:
2851: . offsets - The offsets
2853: Level: beginner
2855: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2856: @*/
2857: PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[])
2858: {
2859: PetscFunctionBegin;
2861: PetscAssertPointer(offsets, 2);
2862: PetscCall(PetscDSSetUp(prob));
2863: *offsets = prob->off;
2864: PetscFunctionReturn(PETSC_SUCCESS);
2865: }
2867: /*@
2868: PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point
2870: Not Collective
2872: Input Parameter:
2873: . prob - The `PetscDS` object
2875: Output Parameter:
2876: . offsets - The offsets
2878: Level: beginner
2880: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2881: @*/
2882: PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[])
2883: {
2884: PetscFunctionBegin;
2886: PetscAssertPointer(offsets, 2);
2887: PetscCall(PetscDSSetUp(prob));
2888: *offsets = prob->offDer;
2889: PetscFunctionReturn(PETSC_SUCCESS);
2890: }
2892: /*@
2893: PetscDSGetComponentOffsetsCohesive - Returns the offset of each field on an evaluation point
2895: Not Collective
2897: Input Parameters:
2898: + ds - The `PetscDS` object
2899: - s - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive
2901: Output Parameter:
2902: . offsets - The offsets
2904: Level: beginner
2906: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2907: @*/
2908: PetscErrorCode PetscDSGetComponentOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
2909: {
2910: PetscFunctionBegin;
2912: PetscAssertPointer(offsets, 3);
2913: PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS");
2914: PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s);
2915: PetscCall(PetscDSSetUp(ds));
2916: *offsets = ds->offCohesive[s];
2917: PetscFunctionReturn(PETSC_SUCCESS);
2918: }
2920: /*@
2921: PetscDSGetComponentDerivativeOffsetsCohesive - Returns the offset of each field derivative on an evaluation point
2923: Not Collective
2925: Input Parameters:
2926: + ds - The `PetscDS` object
2927: - s - The cohesive side, 0 for negative, 1 for positive, 2 for cohesive
2929: Output Parameter:
2930: . offsets - The offsets
2932: Level: beginner
2934: .seealso: `PetscDS`, `PetscDSGetNumFields()`, `PetscDSCreate()`
2935: @*/
2936: PetscErrorCode PetscDSGetComponentDerivativeOffsetsCohesive(PetscDS ds, PetscInt s, PetscInt *offsets[])
2937: {
2938: PetscFunctionBegin;
2940: PetscAssertPointer(offsets, 3);
2941: PetscCheck(ds->isCohesive, PETSC_COMM_SELF, PETSC_ERR_ARG_WRONGSTATE, "Cohesive offsets are only valid for a cohesive DS");
2942: PetscCheck(!(s < 0) && !(s > 2), PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Cohesive side %" PetscInt_FMT " is not in [0, 2]", s);
2943: PetscCall(PetscDSSetUp(ds));
2944: *offsets = ds->offDerCohesive[s];
2945: PetscFunctionReturn(PETSC_SUCCESS);
2946: }
2948: /*@C
2949: PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization
2951: Not Collective
2953: Input Parameter:
2954: . prob - The `PetscDS` object
2956: Output Parameter:
2957: . T - The basis function and derivatives tabulation at quadrature points for each field, see `PetscTabulation` for its details
2959: Level: intermediate
2961: Note:
2962: The tabulation is only valid so long as the `PetscDS` has not be destroyed. There is no `PetscDSRestoreTabulation()` in C.
2964: Fortran Note:
2965: Use the declaration
2966: .vb
2967: PetscTabulation, pointer :: tab(:)
2968: .ve
2969: and access the values using, for example,
2970: .vb
2971: tab(i)%ptr%K
2972: tab(i)%ptr%T(j)%ptr
2973: .ve
2974: where $ i = 1, 2, ..., Nf $ and $ j = 1, 2, ..., tab(i)%ptr%K+1 $.
2976: Use `PetscDSRestoreTabulation()` to restore the array
2978: Developer Note:
2979: The Fortran language syntax does not directly support arrays of pointers, the '%ptr' notation allows mimicking their use in Fortran.
2981: .seealso: `PetscDS`, `PetscTabulation`, `PetscDSCreate()`
2982: @*/
2983: PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscTabulation *T[]) PeNS
2984: {
2985: PetscFunctionBegin;
2987: PetscAssertPointer(T, 2);
2988: PetscCall(PetscDSSetUp(prob));
2989: *T = prob->T;
2990: PetscFunctionReturn(PETSC_SUCCESS);
2991: }
2993: /*@C
2994: PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces
2996: Not Collective
2998: Input Parameter:
2999: . prob - The `PetscDS` object
3001: Output Parameter:
3002: . Tf - The basis function and derivative tabulation on each local face at quadrature points for each field
3004: Level: intermediate
3006: Note:
3007: The tabulation is only valid so long as the `PetscDS` has not be destroyed. There is no `PetscDSRestoreFaceTabulation()` in C.
3009: .seealso: `PetscTabulation`, `PetscDS`, `PetscDSGetTabulation()`, `PetscDSCreate()`
3010: @*/
3011: PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscTabulation *Tf[])
3012: {
3013: PetscFunctionBegin;
3015: PetscAssertPointer(Tf, 2);
3016: PetscCall(PetscDSSetUp(prob));
3017: *Tf = prob->Tf;
3018: PetscFunctionReturn(PETSC_SUCCESS);
3019: }
3021: PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar *u[], PetscScalar *u_t[], PetscScalar *u_x[])
3022: {
3023: PetscFunctionBegin;
3025: PetscCall(PetscDSSetUp(prob));
3026: if (u) {
3027: PetscAssertPointer(u, 2);
3028: *u = prob->u;
3029: }
3030: if (u_t) {
3031: PetscAssertPointer(u_t, 3);
3032: *u_t = prob->u_t;
3033: }
3034: if (u_x) {
3035: PetscAssertPointer(u_x, 4);
3036: *u_x = prob->u_x;
3037: }
3038: PetscFunctionReturn(PETSC_SUCCESS);
3039: }
3041: PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar *f0[], PetscScalar *f1[], PetscScalar *g0[], PetscScalar *g1[], PetscScalar *g2[], PetscScalar *g3[])
3042: {
3043: PetscFunctionBegin;
3045: PetscCall(PetscDSSetUp(prob));
3046: if (f0) {
3047: PetscAssertPointer(f0, 2);
3048: *f0 = prob->f0;
3049: }
3050: if (f1) {
3051: PetscAssertPointer(f1, 3);
3052: *f1 = prob->f1;
3053: }
3054: if (g0) {
3055: PetscAssertPointer(g0, 4);
3056: *g0 = prob->g0;
3057: }
3058: if (g1) {
3059: PetscAssertPointer(g1, 5);
3060: *g1 = prob->g1;
3061: }
3062: if (g2) {
3063: PetscAssertPointer(g2, 6);
3064: *g2 = prob->g2;
3065: }
3066: if (g3) {
3067: PetscAssertPointer(g3, 7);
3068: *g3 = prob->g3;
3069: }
3070: PetscFunctionReturn(PETSC_SUCCESS);
3071: }
3073: PetscErrorCode PetscDSGetWorkspace(PetscDS prob, PetscReal **x, PetscScalar **basisReal, PetscScalar **basisDerReal, PetscScalar **testReal, PetscScalar **testDerReal)
3074: {
3075: PetscFunctionBegin;
3077: PetscCall(PetscDSSetUp(prob));
3078: if (x) {
3079: PetscAssertPointer(x, 2);
3080: *x = prob->x;
3081: }
3082: if (basisReal) {
3083: PetscAssertPointer(basisReal, 3);
3084: *basisReal = prob->basisReal;
3085: }
3086: if (basisDerReal) {
3087: PetscAssertPointer(basisDerReal, 4);
3088: *basisDerReal = prob->basisDerReal;
3089: }
3090: if (testReal) {
3091: PetscAssertPointer(testReal, 5);
3092: *testReal = prob->testReal;
3093: }
3094: if (testDerReal) {
3095: PetscAssertPointer(testDerReal, 6);
3096: *testDerReal = prob->testDerReal;
3097: }
3098: PetscFunctionReturn(PETSC_SUCCESS);
3099: }
3101: /*@C
3102: PetscDSAddBoundary - Add a boundary condition to the model.
3104: Collective
3106: Input Parameters:
3107: + ds - The `PetscDS` object
3108: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3109: . name - The name for the boundary condition
3110: . label - The label defining constrained points
3111: . Nv - The number of `DMLabel` values for constrained points
3112: . values - An array of label values for constrained points
3113: . field - The field to constrain
3114: . Nc - The number of constrained field components (0 will constrain all fields)
3115: . comps - An array of constrained component numbers
3116: . bcFunc - A pointwise function giving boundary values
3117: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or `NULL`
3118: - ctx - An optional user context for `bcFunc`
3120: Output Parameter:
3121: . bd - The boundary number
3123: Options Database Keys:
3124: + -bc_<boundary name> <num> - Overrides the boundary ids
3125: - -bc_<boundary name>_comp <num> - Overrides the boundary components
3127: Level: developer
3129: Note:
3130: Both `bcFunc` and `bcFunc_t` will depend on the boundary condition type. If the type if `DM_BC_ESSENTIAL`, then the calling sequence is\:
3131: .vb
3132: void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])
3133: .ve
3135: If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value, then the calling sequence is\:
3136: .vb
3137: void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3138: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3139: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3140: PetscReal time, const PetscReal x[], PetscScalar bcval[])
3141: .ve
3142: + dim - the coordinate dimension
3143: . Nf - the number of fields
3144: . uOff - the offset into u[] and u_t[] for each field
3145: . uOff_x - the offset into u_x[] for each field
3146: . u - each field evaluated at the current point
3147: . u_t - the time derivative of each field evaluated at the current point
3148: . u_x - the gradient of each field evaluated at the current point
3149: . aOff - the offset into a[] and a_t[] for each auxiliary field
3150: . aOff_x - the offset into a_x[] for each auxiliary field
3151: . a - each auxiliary field evaluated at the current point
3152: . a_t - the time derivative of each auxiliary field evaluated at the current point
3153: . a_x - the gradient of auxiliary each field evaluated at the current point
3154: . t - current time
3155: . x - coordinates of the current point
3156: . numConstants - number of constant parameters
3157: . constants - constant parameters
3158: - bcval - output values at the current point
3160: Notes:
3161: The pointwise functions are used to provide boundary values for essential boundary
3162: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3163: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3164: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3166: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundaryByName()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3167: @*/
3168: 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)
3169: {
3170: DSBoundary head = ds->boundary, b;
3171: PetscInt n = 0;
3172: const char *lname;
3174: PetscFunctionBegin;
3177: PetscAssertPointer(name, 3);
3182: 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);
3183: if (Nc > 0) {
3184: PetscInt *fcomps;
3185: PetscInt c;
3187: PetscCall(PetscDSGetComponents(ds, &fcomps));
3188: 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);
3189: for (c = 0; c < Nc; ++c) {
3190: 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);
3191: }
3192: }
3193: PetscCall(PetscNew(&b));
3194: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3195: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf));
3196: PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf));
3197: PetscCall(PetscMalloc1(Nv, &b->values));
3198: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3199: PetscCall(PetscMalloc1(Nc, &b->comps));
3200: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3201: PetscCall(PetscObjectGetName((PetscObject)label, &lname));
3202: PetscCall(PetscStrallocpy(lname, (char **)&b->lname));
3203: b->type = type;
3204: b->label = label;
3205: b->Nv = Nv;
3206: b->field = field;
3207: b->Nc = Nc;
3208: b->func = bcFunc;
3209: b->func_t = bcFunc_t;
3210: b->ctx = ctx;
3211: b->next = NULL;
3212: /* Append to linked list so that we can preserve the order */
3213: if (!head) ds->boundary = b;
3214: while (head) {
3215: if (!head->next) {
3216: head->next = b;
3217: head = b;
3218: }
3219: head = head->next;
3220: ++n;
3221: }
3222: if (bd) {
3223: PetscAssertPointer(bd, 13);
3224: *bd = n;
3225: }
3226: PetscFunctionReturn(PETSC_SUCCESS);
3227: }
3229: // PetscClangLinter pragma ignore: -fdoc-section-header-unknown
3230: /*@C
3231: PetscDSAddBoundaryByName - Add a boundary condition to the model.
3233: Collective
3235: Input Parameters:
3236: + ds - The `PetscDS` object
3237: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3238: . name - The boundary condition name
3239: . lname - The name of the label defining constrained points
3240: . Nv - The number of `DMLabel` values for constrained points
3241: . values - An array of label values for constrained points
3242: . field - The field to constrain
3243: . Nc - The number of constrained field components (0 will constrain all fields)
3244: . comps - An array of constrained component numbers
3245: . bcFunc - A pointwise function giving boundary values
3246: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or `NULL`
3247: - ctx - An optional user context for `bcFunc`
3249: Output Parameter:
3250: . bd - The boundary number
3252: Options Database Keys:
3253: + -bc_<boundary name> <num> - Overrides the boundary ids
3254: - -bc_<boundary name>_comp <num> - Overrides the boundary components
3256: Calling Sequence of `bcFunc` and `bcFunc_t`:
3257: If the type is `DM_BC_ESSENTIAL`
3258: .vb
3259: void bcFunc(PetscInt dim, PetscReal time, const PetscReal x[], PetscInt Nc, PetscScalar bcval[])
3260: .ve
3261: If the type is `DM_BC_ESSENTIAL_FIELD` or other _FIELD value,
3262: .vb
3263: void bcFunc(PetscInt dim, PetscInt Nf, PetscInt NfAux,
3264: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
3265: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
3266: PetscReal time, const PetscReal x[], PetscScalar bcval[])
3267: .ve
3268: + dim - the coordinate dimension
3269: . Nf - the number of fields
3270: . uOff - the offset into `u`[] and `u_t`[] for each field
3271: . uOff_x - the offset into `u_x`[] for each field
3272: . u - each field evaluated at the current point
3273: . u_t - the time derivative of each field evaluated at the current point
3274: . u_x - the gradient of each field evaluated at the current point
3275: . aOff - the offset into `a`[] and `a_t`[] for each auxiliary field
3276: . aOff_x - the offset into `a_x`[] for each auxiliary field
3277: . a - each auxiliary field evaluated at the current point
3278: . a_t - the time derivative of each auxiliary field evaluated at the current point
3279: . a_x - the gradient of auxiliary each field evaluated at the current point
3280: . t - current time
3281: . x - coordinates of the current point
3282: . numConstants - number of constant parameters
3283: . constants - constant parameters
3284: - bcval - output values at the current point
3286: Level: developer
3288: Notes:
3289: The pointwise functions are used to provide boundary values for essential boundary
3290: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3291: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3292: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3294: This function should only be used with `DMFOREST` currently, since labels cannot be defined before the underlying `DMPLEX` is built.
3296: .seealso: `PetscDS`, `PetscWeakForm`, `DMLabel`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSSetResidual()`, `PetscDSSetBdResidual()`
3297: @*/
3298: 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)
3299: {
3300: DSBoundary head = ds->boundary, b;
3301: PetscInt n = 0;
3303: PetscFunctionBegin;
3306: PetscAssertPointer(name, 3);
3307: PetscAssertPointer(lname, 4);
3311: PetscCall(PetscNew(&b));
3312: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3313: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &b->wf));
3314: PetscCall(PetscWeakFormSetNumFields(b->wf, ds->Nf));
3315: PetscCall(PetscMalloc1(Nv, &b->values));
3316: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3317: PetscCall(PetscMalloc1(Nc, &b->comps));
3318: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3319: PetscCall(PetscStrallocpy(lname, (char **)&b->lname));
3320: b->type = type;
3321: b->label = NULL;
3322: b->Nv = Nv;
3323: b->field = field;
3324: b->Nc = Nc;
3325: b->func = bcFunc;
3326: b->func_t = bcFunc_t;
3327: b->ctx = ctx;
3328: b->next = NULL;
3329: /* Append to linked list so that we can preserve the order */
3330: if (!head) ds->boundary = b;
3331: while (head) {
3332: if (!head->next) {
3333: head->next = b;
3334: head = b;
3335: }
3336: head = head->next;
3337: ++n;
3338: }
3339: if (bd) {
3340: PetscAssertPointer(bd, 13);
3341: *bd = n;
3342: }
3343: PetscFunctionReturn(PETSC_SUCCESS);
3344: }
3346: /*@C
3347: PetscDSUpdateBoundary - Change a boundary condition for the model.
3349: Input Parameters:
3350: + ds - The `PetscDS` object
3351: . bd - The boundary condition number
3352: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3353: . name - The boundary condition name
3354: . label - The label defining constrained points
3355: . Nv - The number of `DMLabel` ids for constrained points
3356: . values - An array of ids for constrained points
3357: . field - The field to constrain
3358: . Nc - The number of constrained field components
3359: . comps - An array of constrained component numbers
3360: . bcFunc - A pointwise function giving boundary values
3361: . bcFunc_t - A pointwise function giving the time derivative of the boundary values, or `NULL`
3362: - ctx - An optional user context for `bcFunc`
3364: Level: developer
3366: Notes:
3367: The pointwise functions are used to provide boundary values for essential boundary
3368: conditions. In FEM, they are acting upon by dual basis functionals to generate FEM
3369: coefficients which are fixed. Natural boundary conditions signal to PETSc that boundary
3370: integrals should be performed, using the kernels from `PetscDSSetBdResidual()`.
3372: The boundary condition number is the order in which it was registered. The user can get the number of boundary conditions from `PetscDSGetNumBoundary()`.
3373: See `PetscDSAddBoundary()` for a description of the calling sequences for the callbacks.
3375: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`, `PetscDSGetNumBoundary()`, `DMLabel`
3376: @*/
3377: 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)
3378: {
3379: DSBoundary b = ds->boundary;
3380: PetscInt n = 0;
3382: PetscFunctionBegin;
3384: while (b) {
3385: if (n == bd) break;
3386: b = b->next;
3387: ++n;
3388: }
3389: PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n);
3390: if (name) {
3391: PetscCall(PetscFree(b->name));
3392: PetscCall(PetscStrallocpy(name, (char **)&b->name));
3393: }
3394: b->type = type;
3395: if (label) {
3396: const char *name;
3398: b->label = label;
3399: PetscCall(PetscFree(b->lname));
3400: PetscCall(PetscObjectGetName((PetscObject)label, &name));
3401: PetscCall(PetscStrallocpy(name, (char **)&b->lname));
3402: }
3403: if (Nv >= 0) {
3404: b->Nv = Nv;
3405: PetscCall(PetscFree(b->values));
3406: PetscCall(PetscMalloc1(Nv, &b->values));
3407: if (Nv) PetscCall(PetscArraycpy(b->values, values, Nv));
3408: }
3409: if (field >= 0) b->field = field;
3410: if (Nc >= 0) {
3411: b->Nc = Nc;
3412: PetscCall(PetscFree(b->comps));
3413: PetscCall(PetscMalloc1(Nc, &b->comps));
3414: if (Nc) PetscCall(PetscArraycpy(b->comps, comps, Nc));
3415: }
3416: if (bcFunc) b->func = bcFunc;
3417: if (bcFunc_t) b->func_t = bcFunc_t;
3418: if (ctx) b->ctx = ctx;
3419: PetscFunctionReturn(PETSC_SUCCESS);
3420: }
3422: /*@
3423: PetscDSGetNumBoundary - Get the number of registered boundary conditions
3425: Input Parameter:
3426: . ds - The `PetscDS` object
3428: Output Parameter:
3429: . numBd - The number of boundary conditions
3431: Level: intermediate
3433: .seealso: `PetscDS`, `PetscDSAddBoundary()`, `PetscDSGetBoundary()`
3434: @*/
3435: PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd)
3436: {
3437: DSBoundary b = ds->boundary;
3439: PetscFunctionBegin;
3441: PetscAssertPointer(numBd, 2);
3442: *numBd = 0;
3443: while (b) {
3444: ++(*numBd);
3445: b = b->next;
3446: }
3447: PetscFunctionReturn(PETSC_SUCCESS);
3448: }
3450: /*@C
3451: PetscDSGetBoundary - Gets a boundary condition from the model
3453: Input Parameters:
3454: + ds - The `PetscDS` object
3455: - bd - The boundary condition number
3457: Output Parameters:
3458: + wf - The `PetscWeakForm` holding the pointwise functions
3459: . type - The type of condition, e.g. `DM_BC_ESSENTIAL`/`DM_BC_ESSENTIAL_FIELD` (Dirichlet), or `DM_BC_NATURAL` (Neumann)
3460: . name - The boundary condition name
3461: . label - The label defining constrained points
3462: . Nv - The number of `DMLabel` ids for constrained points
3463: . values - An array of ids for constrained points
3464: . field - The field to constrain
3465: . Nc - The number of constrained field components
3466: . comps - An array of constrained component numbers
3467: . func - A pointwise function giving boundary values
3468: . func_t - A pointwise function giving the time derivative of the boundary values
3469: - ctx - An optional user context for `bcFunc`
3471: Options Database Keys:
3472: + -bc_<boundary name> <num> - Overrides the boundary ids
3473: - -bc_<boundary name>_comp <num> - Overrides the boundary components
3475: Level: developer
3477: .seealso: `PetscDS`, `PetscWeakForm`, `DMBoundaryConditionType`, `PetscDSAddBoundary()`, `DMLabel`
3478: @*/
3479: 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, void **ctx)
3480: {
3481: DSBoundary b = ds->boundary;
3482: PetscInt n = 0;
3484: PetscFunctionBegin;
3486: while (b) {
3487: if (n == bd) break;
3488: b = b->next;
3489: ++n;
3490: }
3491: PetscCheck(b, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", bd, n);
3492: if (wf) {
3493: PetscAssertPointer(wf, 3);
3494: *wf = b->wf;
3495: }
3496: if (type) {
3497: PetscAssertPointer(type, 4);
3498: *type = b->type;
3499: }
3500: if (name) {
3501: PetscAssertPointer(name, 5);
3502: *name = b->name;
3503: }
3504: if (label) {
3505: PetscAssertPointer(label, 6);
3506: *label = b->label;
3507: }
3508: if (Nv) {
3509: PetscAssertPointer(Nv, 7);
3510: *Nv = b->Nv;
3511: }
3512: if (values) {
3513: PetscAssertPointer(values, 8);
3514: *values = b->values;
3515: }
3516: if (field) {
3517: PetscAssertPointer(field, 9);
3518: *field = b->field;
3519: }
3520: if (Nc) {
3521: PetscAssertPointer(Nc, 10);
3522: *Nc = b->Nc;
3523: }
3524: if (comps) {
3525: PetscAssertPointer(comps, 11);
3526: *comps = b->comps;
3527: }
3528: if (func) {
3529: PetscAssertPointer(func, 12);
3530: *func = b->func;
3531: }
3532: if (func_t) {
3533: PetscAssertPointer(func_t, 13);
3534: *func_t = b->func_t;
3535: }
3536: if (ctx) {
3537: PetscAssertPointer(ctx, 14);
3538: *ctx = b->ctx;
3539: }
3540: PetscFunctionReturn(PETSC_SUCCESS);
3541: }
3543: /*@
3544: PetscDSUpdateBoundaryLabels - Update `DMLabel` in each boundary condition using the label name and the input `DM`
3546: Not Collective
3548: Input Parameters:
3549: + ds - The source `PetscDS` object
3550: - dm - The `DM` holding labels
3552: Level: intermediate
3554: .seealso: `PetscDS`, `DMBoundary`, `DM`, `PetscDSCopyBoundary()`, `PetscDSCreate()`, `DMGetLabel()`
3555: @*/
3556: PetscErrorCode PetscDSUpdateBoundaryLabels(PetscDS ds, DM dm)
3557: {
3558: DSBoundary b;
3560: PetscFunctionBegin;
3563: for (b = ds->boundary; b; b = b->next) {
3564: if (b->lname) PetscCall(DMGetLabel(dm, b->lname, &b->label));
3565: }
3566: PetscFunctionReturn(PETSC_SUCCESS);
3567: }
3569: static PetscErrorCode DSBoundaryDuplicate_Internal(DSBoundary b, DSBoundary *bNew)
3570: {
3571: PetscFunctionBegin;
3572: PetscCall(PetscNew(bNew));
3573: PetscCall(PetscWeakFormCreate(PETSC_COMM_SELF, &(*bNew)->wf));
3574: PetscCall(PetscWeakFormCopy(b->wf, (*bNew)->wf));
3575: PetscCall(PetscStrallocpy(b->name, (char **)&((*bNew)->name)));
3576: PetscCall(PetscStrallocpy(b->lname, (char **)&((*bNew)->lname)));
3577: (*bNew)->type = b->type;
3578: (*bNew)->label = b->label;
3579: (*bNew)->Nv = b->Nv;
3580: PetscCall(PetscMalloc1(b->Nv, &(*bNew)->values));
3581: PetscCall(PetscArraycpy((*bNew)->values, b->values, b->Nv));
3582: (*bNew)->field = b->field;
3583: (*bNew)->Nc = b->Nc;
3584: PetscCall(PetscMalloc1(b->Nc, &(*bNew)->comps));
3585: PetscCall(PetscArraycpy((*bNew)->comps, b->comps, b->Nc));
3586: (*bNew)->func = b->func;
3587: (*bNew)->func_t = b->func_t;
3588: (*bNew)->ctx = b->ctx;
3589: PetscFunctionReturn(PETSC_SUCCESS);
3590: }
3592: /*@
3593: PetscDSCopyBoundary - Copy all boundary condition objects to the new `PetscDS`
3595: Not Collective
3597: Input Parameters:
3598: + ds - The source `PetscDS` object
3599: . numFields - The number of selected fields, or `PETSC_DEFAULT` for all fields
3600: - fields - The selected fields, or `NULL` for all fields
3602: Output Parameter:
3603: . newds - The target `PetscDS`, now with a copy of the boundary conditions
3605: Level: intermediate
3607: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3608: @*/
3609: PetscErrorCode PetscDSCopyBoundary(PetscDS ds, PetscInt numFields, const PetscInt fields[], PetscDS newds)
3610: {
3611: DSBoundary b, *lastnext;
3613: PetscFunctionBegin;
3616: if (ds == newds) PetscFunctionReturn(PETSC_SUCCESS);
3617: PetscCall(PetscDSDestroyBoundary(newds));
3618: lastnext = &newds->boundary;
3619: for (b = ds->boundary; b; b = b->next) {
3620: DSBoundary bNew;
3621: PetscInt fieldNew = -1;
3623: if (numFields > 0 && fields) {
3624: PetscInt f;
3626: for (f = 0; f < numFields; ++f)
3627: if (b->field == fields[f]) break;
3628: if (f == numFields) continue;
3629: fieldNew = f;
3630: }
3631: PetscCall(DSBoundaryDuplicate_Internal(b, &bNew));
3632: bNew->field = fieldNew < 0 ? b->field : fieldNew;
3633: *lastnext = bNew;
3634: lastnext = &bNew->next;
3635: }
3636: PetscFunctionReturn(PETSC_SUCCESS);
3637: }
3639: /*@
3640: PetscDSDestroyBoundary - Remove all `DMBoundary` objects from the `PetscDS`
3642: Not Collective
3644: Input Parameter:
3645: . ds - The `PetscDS` object
3647: Level: intermediate
3649: .seealso: `PetscDS`, `DMBoundary`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`
3650: @*/
3651: PetscErrorCode PetscDSDestroyBoundary(PetscDS ds)
3652: {
3653: DSBoundary next = ds->boundary;
3655: PetscFunctionBegin;
3656: while (next) {
3657: DSBoundary b = next;
3659: next = b->next;
3660: PetscCall(PetscWeakFormDestroy(&b->wf));
3661: PetscCall(PetscFree(b->name));
3662: PetscCall(PetscFree(b->lname));
3663: PetscCall(PetscFree(b->values));
3664: PetscCall(PetscFree(b->comps));
3665: PetscCall(PetscFree(b));
3666: }
3667: PetscFunctionReturn(PETSC_SUCCESS);
3668: }
3670: /*@
3671: PetscDSSelectDiscretizations - Copy discretizations to the new `PetscDS` with different field layout
3673: Not Collective
3675: Input Parameters:
3676: + prob - The `PetscDS` object
3677: . numFields - Number of new fields
3678: . fields - Old field number for each new field
3679: . minDegree - Minimum degree for a discretization, or `PETSC_DETERMINE` for no limit
3680: - maxDegree - Maximum degree for a discretization, or `PETSC_DETERMINE` for no limit
3682: Output Parameter:
3683: . newprob - The `PetscDS` copy
3685: Level: intermediate
3687: .seealso: `PetscDS`, `PetscDSSelectEquations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3688: @*/
3689: PetscErrorCode PetscDSSelectDiscretizations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscInt minDegree, PetscInt maxDegree, PetscDS newprob)
3690: {
3691: PetscInt Nf, Nfn, fn;
3693: PetscFunctionBegin;
3695: if (fields) PetscAssertPointer(fields, 3);
3697: PetscCall(PetscDSGetNumFields(prob, &Nf));
3698: PetscCall(PetscDSGetNumFields(newprob, &Nfn));
3699: numFields = numFields < 0 ? Nf : numFields;
3700: for (fn = 0; fn < numFields; ++fn) {
3701: const PetscInt f = fields ? fields[fn] : fn;
3702: PetscObject disc;
3703: PetscClassId id;
3705: if (f >= Nf) continue;
3706: PetscCall(PetscDSGetDiscretization(prob, f, &disc));
3707: PetscCallContinue(PetscObjectGetClassId(disc, &id));
3708: if (id == PETSCFE_CLASSID) {
3709: PetscFE fe;
3711: PetscCall(PetscFELimitDegree((PetscFE)disc, minDegree, maxDegree, &fe));
3712: PetscCall(PetscDSSetDiscretization(newprob, fn, (PetscObject)fe));
3713: PetscCall(PetscFEDestroy(&fe));
3714: } else {
3715: PetscCall(PetscDSSetDiscretization(newprob, fn, disc));
3716: }
3717: }
3718: PetscFunctionReturn(PETSC_SUCCESS);
3719: }
3721: /*@
3722: PetscDSSelectEquations - Copy pointwise function pointers to the new `PetscDS` with different field layout
3724: Not Collective
3726: Input Parameters:
3727: + prob - The `PetscDS` object
3728: . numFields - Number of new fields
3729: - fields - Old field number for each new field
3731: Output Parameter:
3732: . newprob - The `PetscDS` copy
3734: Level: intermediate
3736: .seealso: `PetscDS`, `PetscDSSelectDiscretizations()`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3737: @*/
3738: PetscErrorCode PetscDSSelectEquations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
3739: {
3740: PetscInt Nf, Nfn, fn, gn;
3742: PetscFunctionBegin;
3744: if (fields) PetscAssertPointer(fields, 3);
3746: PetscCall(PetscDSGetNumFields(prob, &Nf));
3747: PetscCall(PetscDSGetNumFields(newprob, &Nfn));
3748: 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);
3749: for (fn = 0; fn < numFields; ++fn) {
3750: const PetscInt f = fields ? fields[fn] : fn;
3751: PetscPointFn *obj;
3752: PetscPointFn *f0, *f1;
3753: PetscBdPointFn *f0Bd, *f1Bd;
3754: PetscRiemannFn *r;
3756: if (f >= Nf) continue;
3757: PetscCall(PetscDSGetObjective(prob, f, &obj));
3758: PetscCall(PetscDSGetResidual(prob, f, &f0, &f1));
3759: PetscCall(PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd));
3760: PetscCall(PetscDSGetRiemannSolver(prob, f, &r));
3761: PetscCall(PetscDSSetObjective(newprob, fn, obj));
3762: PetscCall(PetscDSSetResidual(newprob, fn, f0, f1));
3763: PetscCall(PetscDSSetBdResidual(newprob, fn, f0Bd, f1Bd));
3764: PetscCall(PetscDSSetRiemannSolver(newprob, fn, r));
3765: for (gn = 0; gn < numFields; ++gn) {
3766: const PetscInt g = fields ? fields[gn] : gn;
3767: PetscPointJacFn *g0, *g1, *g2, *g3;
3768: PetscPointJacFn *g0p, *g1p, *g2p, *g3p;
3769: PetscBdPointJacFn *g0Bd, *g1Bd, *g2Bd, *g3Bd;
3771: if (g >= Nf) continue;
3772: PetscCall(PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3));
3773: PetscCall(PetscDSGetJacobianPreconditioner(prob, f, g, &g0p, &g1p, &g2p, &g3p));
3774: PetscCall(PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd));
3775: PetscCall(PetscDSSetJacobian(newprob, fn, gn, g0, g1, g2, g3));
3776: PetscCall(PetscDSSetJacobianPreconditioner(newprob, fn, gn, g0p, g1p, g2p, g3p));
3777: PetscCall(PetscDSSetBdJacobian(newprob, fn, gn, g0Bd, g1Bd, g2Bd, g3Bd));
3778: }
3779: }
3780: PetscFunctionReturn(PETSC_SUCCESS);
3781: }
3783: /*@
3784: PetscDSCopyEquations - Copy all pointwise function pointers to another `PetscDS`
3786: Not Collective
3788: Input Parameter:
3789: . prob - The `PetscDS` object
3791: Output Parameter:
3792: . newprob - The `PetscDS` copy
3794: Level: intermediate
3796: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3797: @*/
3798: PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob)
3799: {
3800: PetscWeakForm wf, newwf;
3801: PetscInt Nf, Ng;
3803: PetscFunctionBegin;
3806: PetscCall(PetscDSGetNumFields(prob, &Nf));
3807: PetscCall(PetscDSGetNumFields(newprob, &Ng));
3808: PetscCheck(Nf == Ng, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %" PetscInt_FMT " != %" PetscInt_FMT, Nf, Ng);
3809: PetscCall(PetscDSGetWeakForm(prob, &wf));
3810: PetscCall(PetscDSGetWeakForm(newprob, &newwf));
3811: PetscCall(PetscWeakFormCopy(wf, newwf));
3812: PetscFunctionReturn(PETSC_SUCCESS);
3813: }
3815: /*@
3816: PetscDSCopyConstants - Copy all constants set with `PetscDSSetConstants()` to another `PetscDS`
3818: Not Collective
3820: Input Parameter:
3821: . prob - The `PetscDS` object
3823: Output Parameter:
3824: . newprob - The `PetscDS` copy
3826: Level: intermediate
3828: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3829: @*/
3830: PetscErrorCode PetscDSCopyConstants(PetscDS prob, PetscDS newprob)
3831: {
3832: PetscInt Nc;
3833: const PetscScalar *constants;
3835: PetscFunctionBegin;
3838: PetscCall(PetscDSGetConstants(prob, &Nc, &constants));
3839: PetscCall(PetscDSSetConstants(newprob, Nc, (PetscScalar *)constants));
3840: PetscFunctionReturn(PETSC_SUCCESS);
3841: }
3843: /*@
3844: PetscDSCopyExactSolutions - Copy all exact solutions set with `PetscDSSetExactSolution()` and `PetscDSSetExactSolutionTimeDerivative()` to another `PetscDS`
3846: Not Collective
3848: Input Parameter:
3849: . ds - The `PetscDS` object
3851: Output Parameter:
3852: . newds - The `PetscDS` copy
3854: Level: intermediate
3856: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSCopyBounds()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3857: @*/
3858: PetscErrorCode PetscDSCopyExactSolutions(PetscDS ds, PetscDS newds)
3859: {
3860: PetscSimplePointFn *sol;
3861: void *ctx;
3862: PetscInt Nf, f;
3864: PetscFunctionBegin;
3867: PetscCall(PetscDSGetNumFields(ds, &Nf));
3868: for (f = 0; f < Nf; ++f) {
3869: PetscCall(PetscDSGetExactSolution(ds, f, &sol, &ctx));
3870: PetscCall(PetscDSSetExactSolution(newds, f, sol, ctx));
3871: PetscCall(PetscDSGetExactSolutionTimeDerivative(ds, f, &sol, &ctx));
3872: PetscCall(PetscDSSetExactSolutionTimeDerivative(newds, f, sol, ctx));
3873: }
3874: PetscFunctionReturn(PETSC_SUCCESS);
3875: }
3877: /*@
3878: PetscDSCopyBounds - Copy lower and upper solution bounds set with `PetscDSSetLowerBound()` and `PetscDSSetLowerBound()` to another `PetscDS`
3880: Not Collective
3882: Input Parameter:
3883: . ds - The `PetscDS` object
3885: Output Parameter:
3886: . newds - The `PetscDS` copy
3888: Level: intermediate
3890: .seealso: `PetscDS`, `PetscDSCopyBoundary()`, `PetscDSCopyEquations()`, `PetscDSCopyExactSolutions()`, `PetscDSSetResidual()`, `PetscDSSetJacobian()`, `PetscDSSetRiemannSolver()`, `PetscDSSetBdResidual()`, `PetscDSSetBdJacobian()`, `PetscDSCreate()`
3891: @*/
3892: PetscErrorCode PetscDSCopyBounds(PetscDS ds, PetscDS newds)
3893: {
3894: PetscSimplePointFn *bound;
3895: void *ctx;
3896: PetscInt Nf, f;
3898: PetscFunctionBegin;
3901: PetscCall(PetscDSGetNumFields(ds, &Nf));
3902: for (f = 0; f < Nf; ++f) {
3903: PetscCall(PetscDSGetLowerBound(ds, f, &bound, &ctx));
3904: PetscCall(PetscDSSetLowerBound(newds, f, bound, ctx));
3905: PetscCall(PetscDSGetUpperBound(ds, f, &bound, &ctx));
3906: PetscCall(PetscDSSetUpperBound(newds, f, bound, ctx));
3907: }
3908: PetscFunctionReturn(PETSC_SUCCESS);
3909: }
3911: PetscErrorCode PetscDSCopy(PetscDS ds, PetscInt minDegree, PetscInt maxDegree, DM dmNew, PetscDS dsNew)
3912: {
3913: DSBoundary b;
3914: PetscInt cdim, Nf, f, d;
3915: PetscBool isCohesive;
3916: void *ctx;
3918: PetscFunctionBegin;
3919: PetscCall(PetscDSCopyConstants(ds, dsNew));
3920: PetscCall(PetscDSCopyExactSolutions(ds, dsNew));
3921: PetscCall(PetscDSCopyBounds(ds, dsNew));
3922: PetscCall(PetscDSSelectDiscretizations(ds, PETSC_DETERMINE, NULL, minDegree, maxDegree, dsNew));
3923: PetscCall(PetscDSCopyEquations(ds, dsNew));
3924: PetscCall(PetscDSGetNumFields(ds, &Nf));
3925: for (f = 0; f < Nf; ++f) {
3926: PetscCall(PetscDSGetContext(ds, f, &ctx));
3927: PetscCall(PetscDSSetContext(dsNew, f, ctx));
3928: PetscCall(PetscDSGetCohesive(ds, f, &isCohesive));
3929: PetscCall(PetscDSSetCohesive(dsNew, f, isCohesive));
3930: PetscCall(PetscDSGetJetDegree(ds, f, &d));
3931: PetscCall(PetscDSSetJetDegree(dsNew, f, d));
3932: }
3933: if (Nf) {
3934: PetscCall(PetscDSGetCoordinateDimension(ds, &cdim));
3935: PetscCall(PetscDSSetCoordinateDimension(dsNew, cdim));
3936: }
3937: PetscCall(PetscDSCopyBoundary(ds, PETSC_DETERMINE, NULL, dsNew));
3938: for (b = dsNew->boundary; b; b = b->next) {
3939: PetscCall(DMGetLabel(dmNew, b->lname, &b->label));
3940: /* Do not check if label exists here, since p4est calls this for the reference tree which does not have the labels */
3941: //PetscCheck(b->label,PETSC_COMM_SELF, PETSC_ERR_ARG_WRONG, "Label %s missing in new DM", name);
3942: }
3943: PetscFunctionReturn(PETSC_SUCCESS);
3944: }
3946: PetscErrorCode PetscDSGetHeightSubspace(PetscDS prob, PetscInt height, PetscDS *subprob)
3947: {
3948: PetscInt dim, Nf, f;
3950: PetscFunctionBegin;
3952: PetscAssertPointer(subprob, 3);
3953: if (height == 0) {
3954: *subprob = prob;
3955: PetscFunctionReturn(PETSC_SUCCESS);
3956: }
3957: PetscCall(PetscDSGetNumFields(prob, &Nf));
3958: PetscCall(PetscDSGetSpatialDimension(prob, &dim));
3959: PetscCheck(height <= dim, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_OUTOFRANGE, "DS can only handle height in [0, %" PetscInt_FMT "], not %" PetscInt_FMT, dim, height);
3960: if (!prob->subprobs) PetscCall(PetscCalloc1(dim, &prob->subprobs));
3961: if (!prob->subprobs[height - 1]) {
3962: PetscInt cdim;
3964: PetscCall(PetscDSCreate(PetscObjectComm((PetscObject)prob), &prob->subprobs[height - 1]));
3965: PetscCall(PetscDSGetCoordinateDimension(prob, &cdim));
3966: PetscCall(PetscDSSetCoordinateDimension(prob->subprobs[height - 1], cdim));
3967: for (f = 0; f < Nf; ++f) {
3968: PetscFE subfe;
3969: PetscObject obj;
3970: PetscClassId id;
3972: PetscCall(PetscDSGetDiscretization(prob, f, &obj));
3973: PetscCall(PetscObjectGetClassId(obj, &id));
3974: PetscCheck(id == PETSCFE_CLASSID, PetscObjectComm((PetscObject)prob), PETSC_ERR_ARG_WRONG, "Unsupported discretization type for field %" PetscInt_FMT, f);
3975: PetscCall(PetscFEGetHeightSubspace((PetscFE)obj, height, &subfe));
3976: PetscCall(PetscDSSetDiscretization(prob->subprobs[height - 1], f, (PetscObject)subfe));
3977: }
3978: }
3979: *subprob = prob->subprobs[height - 1];
3980: PetscFunctionReturn(PETSC_SUCCESS);
3981: }
3983: PetscErrorCode PetscDSPermuteQuadPoint(PetscDS ds, PetscInt ornt, PetscInt field, PetscInt q, PetscInt *qperm)
3984: {
3985: IS permIS;
3986: PetscQuadrature quad;
3987: DMPolytopeType ct;
3988: const PetscInt *perm;
3989: PetscInt Na, Nq;
3991: PetscFunctionBeginHot;
3992: PetscCall(PetscFEGetQuadrature((PetscFE)ds->disc[field], &quad));
3993: PetscCall(PetscQuadratureGetData(quad, NULL, NULL, &Nq, NULL, NULL));
3994: PetscCall(PetscQuadratureGetCellType(quad, &ct));
3995: PetscCheck(q >= 0 && q < Nq, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Quadrature point %" PetscInt_FMT " is not in [0, %" PetscInt_FMT ")", q, Nq);
3996: Na = DMPolytopeTypeGetNumArrangements(ct) / 2;
3997: 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);
3998: if (!ds->quadPerm[(PetscInt)ct]) PetscCall(PetscQuadratureComputePermutations(quad, NULL, &ds->quadPerm[(PetscInt)ct]));
3999: permIS = ds->quadPerm[(PetscInt)ct][ornt + Na];
4000: PetscCall(ISGetIndices(permIS, &perm));
4001: *qperm = perm[q];
4002: PetscCall(ISRestoreIndices(permIS, &perm));
4003: PetscFunctionReturn(PETSC_SUCCESS);
4004: }
4006: PetscErrorCode PetscDSGetDiscType_Internal(PetscDS ds, PetscInt f, PetscDiscType *disctype)
4007: {
4008: PetscObject obj;
4009: PetscClassId id;
4010: PetscInt Nf;
4012: PetscFunctionBegin;
4014: PetscAssertPointer(disctype, 3);
4015: *disctype = PETSC_DISC_NONE;
4016: PetscCall(PetscDSGetNumFields(ds, &Nf));
4017: PetscCheck(f < Nf, PetscObjectComm((PetscObject)ds), PETSC_ERR_ARG_SIZ, "Field %" PetscInt_FMT " must be in [0, %" PetscInt_FMT ")", f, Nf);
4018: PetscCall(PetscDSGetDiscretization(ds, f, &obj));
4019: if (obj) {
4020: PetscCall(PetscObjectGetClassId(obj, &id));
4021: if (id == PETSCFE_CLASSID) *disctype = PETSC_DISC_FE;
4022: else *disctype = PETSC_DISC_FV;
4023: }
4024: PetscFunctionReturn(PETSC_SUCCESS);
4025: }
4027: static PetscErrorCode PetscDSDestroy_Basic(PetscDS ds)
4028: {
4029: PetscFunctionBegin;
4030: PetscCall(PetscFree(ds->data));
4031: PetscFunctionReturn(PETSC_SUCCESS);
4032: }
4034: static PetscErrorCode PetscDSInitialize_Basic(PetscDS ds)
4035: {
4036: PetscFunctionBegin;
4037: ds->ops->setfromoptions = NULL;
4038: ds->ops->setup = NULL;
4039: ds->ops->view = NULL;
4040: ds->ops->destroy = PetscDSDestroy_Basic;
4041: PetscFunctionReturn(PETSC_SUCCESS);
4042: }
4044: /*MC
4045: PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions
4047: Level: intermediate
4049: .seealso: `PetscDSType`, `PetscDSCreate()`, `PetscDSSetType()`
4050: M*/
4052: PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS ds)
4053: {
4054: PetscDS_Basic *b;
4056: PetscFunctionBegin;
4058: PetscCall(PetscNew(&b));
4059: ds->data = b;
4061: PetscCall(PetscDSInitialize_Basic(ds));
4062: PetscFunctionReturn(PETSC_SUCCESS);
4063: }