Actual source code: dtds.c
petsc-3.12.5 2020-03-29
1: #include <petsc/private/petscdsimpl.h>
3: PetscClassId PETSCDS_CLASSID = 0;
5: PetscFunctionList PetscDSList = NULL;
6: PetscBool PetscDSRegisterAllCalled = PETSC_FALSE;
8: /* A PetscDS (Discrete System) encodes a set of equations posed in a discrete space, which represents a set of
9: nonlinear continuum equations. The equations can have multiple fields, each field having a different
10: discretization. In addition, different pieces of the domain can have different field combinations and equations.
12: The DS provides the user a description of the approximation space on any given cell. It also gives pointwise
13: functions representing the equations.
15: Each field is associated with a label, marking the cells on which it is supported. Note that a field can be
16: supported on the closure of a cell not in the label due to overlap of the boundary of neighboring cells. The DM
17: then creates a DS for each set of cells with identical approximation spaces. When assembling, the user asks for
18: the space associated with a given cell. DMPlex uses the labels associated with each DS in the default integration loop.
19: */
21: /*@C
22: PetscDSRegister - Adds a new PetscDS implementation
24: Not Collective
26: Input Parameters:
27: + name - The name of a new user-defined creation routine
28: - create_func - The creation routine itself
30: Notes:
31: PetscDSRegister() may be called multiple times to add several user-defined PetscDSs
33: Sample usage:
34: .vb
35: PetscDSRegister("my_ds", MyPetscDSCreate);
36: .ve
38: Then, your PetscDS type can be chosen with the procedural interface via
39: .vb
40: PetscDSCreate(MPI_Comm, PetscDS *);
41: PetscDSSetType(PetscDS, "my_ds");
42: .ve
43: or at runtime via the option
44: .vb
45: -petscds_type my_ds
46: .ve
48: Level: advanced
50: Not available from Fortran
52: .seealso: PetscDSRegisterAll(), PetscDSRegisterDestroy()
54: @*/
55: PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS))
56: {
60: PetscFunctionListAdd(&PetscDSList, sname, function);
61: return(0);
62: }
64: /*@C
65: PetscDSSetType - Builds a particular PetscDS
67: Collective on prob
69: Input Parameters:
70: + prob - The PetscDS object
71: - name - The kind of system
73: Options Database Key:
74: . -petscds_type <type> - Sets the PetscDS type; use -help for a list of available types
76: Level: intermediate
78: Not available from Fortran
80: .seealso: PetscDSGetType(), PetscDSCreate()
81: @*/
82: PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name)
83: {
84: PetscErrorCode (*r)(PetscDS);
85: PetscBool match;
90: PetscObjectTypeCompare((PetscObject) prob, name, &match);
91: if (match) return(0);
93: PetscDSRegisterAll();
94: PetscFunctionListFind(PetscDSList, name, &r);
95: if (!r) SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscDS type: %s", name);
97: if (prob->ops->destroy) {
98: (*prob->ops->destroy)(prob);
99: prob->ops->destroy = NULL;
100: }
101: (*r)(prob);
102: PetscObjectChangeTypeName((PetscObject) prob, name);
103: return(0);
104: }
106: /*@C
107: PetscDSGetType - Gets the PetscDS type name (as a string) from the object.
109: Not Collective
111: Input Parameter:
112: . prob - The PetscDS
114: Output Parameter:
115: . name - The PetscDS type name
117: Level: intermediate
119: Not available from Fortran
121: .seealso: PetscDSSetType(), PetscDSCreate()
122: @*/
123: PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name)
124: {
130: PetscDSRegisterAll();
131: *name = ((PetscObject) prob)->type_name;
132: return(0);
133: }
135: static PetscErrorCode PetscDSView_Ascii(PetscDS prob, PetscViewer viewer)
136: {
137: PetscViewerFormat format;
138: const PetscScalar *constants;
139: PetscInt numConstants, f;
140: PetscErrorCode ierr;
143: PetscViewerGetFormat(viewer, &format);
144: PetscViewerASCIIPrintf(viewer, "Discrete System with %d fields\n", prob->Nf);
145: PetscViewerASCIIPushTab(viewer);
146: PetscViewerASCIIPrintf(viewer, " cell total dim %D total comp %D\n", prob->totDim, prob->totComp);
147: if (prob->isHybrid) {PetscViewerASCIIPrintf(viewer, " hybrid cell\n");}
148: for (f = 0; f < prob->Nf; ++f) {
149: DSBoundary b;
150: PetscObject obj;
151: PetscClassId id;
152: PetscQuadrature q;
153: const char *name;
154: PetscInt Nc, Nq, Nqc;
156: PetscDSGetDiscretization(prob, f, &obj);
157: PetscObjectGetClassId(obj, &id);
158: PetscObjectGetName(obj, &name);
159: PetscViewerASCIIPrintf(viewer, "Field %s", name ? name : "<unknown>");
160: PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
161: if (id == PETSCFE_CLASSID) {
162: PetscFEGetNumComponents((PetscFE) obj, &Nc);
163: PetscFEGetQuadrature((PetscFE) obj, &q);
164: PetscViewerASCIIPrintf(viewer, " FEM");
165: } else if (id == PETSCFV_CLASSID) {
166: PetscFVGetNumComponents((PetscFV) obj, &Nc);
167: PetscFVGetQuadrature((PetscFV) obj, &q);
168: PetscViewerASCIIPrintf(viewer, " FVM");
169: }
170: else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %D", f);
171: if (Nc > 1) {PetscViewerASCIIPrintf(viewer, "%D components", Nc);}
172: else {PetscViewerASCIIPrintf(viewer, "%D component ", Nc);}
173: if (prob->implicit[f]) {PetscViewerASCIIPrintf(viewer, " (implicit)");}
174: else {PetscViewerASCIIPrintf(viewer, " (explicit)");}
175: if (q) {
176: PetscQuadratureGetData(q, NULL, &Nqc, &Nq, NULL, NULL);
177: PetscViewerASCIIPrintf(viewer, " (Nq %D Nqc %D)", Nq, Nqc);
178: }
179: PetscViewerASCIIPrintf(viewer, "\n");
180: PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
181: PetscViewerASCIIPushTab(viewer);
182: if (id == PETSCFE_CLASSID) {PetscFEView((PetscFE) obj, viewer);}
183: else if (id == PETSCFV_CLASSID) {PetscFVView((PetscFV) obj, viewer);}
184: PetscViewerASCIIPopTab(viewer);
186: for (b = prob->boundary; b; b = b->next) {
187: PetscInt c, i;
189: if (b->field != f) continue;
190: PetscViewerASCIIPushTab(viewer);
191: PetscViewerASCIIPrintf(viewer, "Boundary %s (%s) %s\n", b->name, b->labelname, DMBoundaryConditionTypes[b->type]);
192: if (!b->numcomps) {
193: PetscViewerASCIIPrintf(viewer, " all components\n");
194: } else {
195: PetscViewerASCIIPrintf(viewer, " components: ");
196: PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
197: for (c = 0; c < b->numcomps; ++c) {
198: if (c > 0) {PetscViewerASCIIPrintf(viewer, ", ");}
199: PetscViewerASCIIPrintf(viewer, "%D", b->comps[c]);
200: }
201: PetscViewerASCIIPrintf(viewer, "\n");
202: PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
203: }
204: PetscViewerASCIIPrintf(viewer, " ids: ");
205: PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
206: for (i = 0; i < b->numids; ++i) {
207: if (i > 0) {PetscViewerASCIIPrintf(viewer, ", ");}
208: PetscViewerASCIIPrintf(viewer, "%D", b->ids[i]);
209: }
210: PetscViewerASCIIPrintf(viewer, "\n");
211: PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
212: PetscViewerASCIIPopTab(viewer);
213: }
214: }
215: PetscDSGetConstants(prob, &numConstants, &constants);
216: if (numConstants) {
217: PetscViewerASCIIPrintf(viewer, "%D constants\n", numConstants);
218: PetscViewerASCIIPushTab(viewer);
219: for (f = 0; f < numConstants; ++f) {PetscViewerASCIIPrintf(viewer, "%g\n", (double) PetscRealPart(constants[f]));}
220: PetscViewerASCIIPopTab(viewer);
221: }
222: PetscViewerASCIIPopTab(viewer);
223: return(0);
224: }
226: /*@C
227: PetscDSView - Views a PetscDS
229: Collective on prob
231: Input Parameter:
232: + prob - the PetscDS object to view
233: - v - the viewer
235: Level: developer
237: .seealso PetscDSDestroy()
238: @*/
239: PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v)
240: {
241: PetscBool iascii;
246: if (!v) {PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject) prob), &v);}
248: PetscObjectTypeCompare((PetscObject) v, PETSCVIEWERASCII, &iascii);
249: if (iascii) {PetscDSView_Ascii(prob, v);}
250: if (prob->ops->view) {(*prob->ops->view)(prob, v);}
251: return(0);
252: }
254: /*@
255: PetscDSSetFromOptions - sets parameters in a PetscDS from the options database
257: Collective on prob
259: Input Parameter:
260: . prob - the PetscDS object to set options for
262: Options Database:
263: + -petscds_type <type> : Set the DS type
264: . -petscds_view <view opt> : View the DS
265: . -petscds_jac_pre : Turn formation of a separate Jacobian preconditioner on and off
266: . -bc_<name> <ids> : Specify a list of label ids for a boundary condition
267: - -bc_<name>_comp <comps> : Specify a list of field components to constrain for a boundary condition
269: Level: developer
271: .seealso PetscDSView()
272: @*/
273: PetscErrorCode PetscDSSetFromOptions(PetscDS prob)
274: {
275: DSBoundary b;
276: const char *defaultType;
277: char name[256];
278: PetscBool flg;
283: if (!((PetscObject) prob)->type_name) {
284: defaultType = PETSCDSBASIC;
285: } else {
286: defaultType = ((PetscObject) prob)->type_name;
287: }
288: PetscDSRegisterAll();
290: PetscObjectOptionsBegin((PetscObject) prob);
291: for (b = prob->boundary; b; b = b->next) {
292: char optname[1024];
293: PetscInt ids[1024], len = 1024;
294: PetscBool flg;
296: PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name);
297: PetscMemzero(ids, sizeof(ids));
298: PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg);
299: if (flg) {
300: b->numids = len;
301: PetscFree(b->ids);
302: PetscMalloc1(len, &b->ids);
303: PetscArraycpy(b->ids, ids, len);
304: }
305: len = 1024;
306: PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name);
307: PetscMemzero(ids, sizeof(ids));
308: PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg);
309: if (flg) {
310: b->numcomps = len;
311: PetscFree(b->comps);
312: PetscMalloc1(len, &b->comps);
313: PetscArraycpy(b->comps, ids, len);
314: }
315: }
316: PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg);
317: if (flg) {
318: PetscDSSetType(prob, name);
319: } else if (!((PetscObject) prob)->type_name) {
320: PetscDSSetType(prob, defaultType);
321: }
322: PetscOptionsBool("-petscds_jac_pre", "Discrete System", "PetscDSUseJacobianPreconditioner", prob->useJacPre, &prob->useJacPre, &flg);
323: if (prob->ops->setfromoptions) {(*prob->ops->setfromoptions)(prob);}
324: /* process any options handlers added with PetscObjectAddOptionsHandler() */
325: PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject) prob);
326: PetscOptionsEnd();
327: if (prob->Nf) {PetscDSViewFromOptions(prob, NULL, "-petscds_view");}
328: return(0);
329: }
331: /*@C
332: PetscDSSetUp - Construct data structures for the PetscDS
334: Collective on prob
336: Input Parameter:
337: . prob - the PetscDS object to setup
339: Level: developer
341: .seealso PetscDSView(), PetscDSDestroy()
342: @*/
343: PetscErrorCode PetscDSSetUp(PetscDS prob)
344: {
345: const PetscInt Nf = prob->Nf;
346: PetscInt dim, dimEmbed, NbMax = 0, NcMax = 0, NqMax = 0, NsMax = 1, f;
351: if (prob->setup) return(0);
352: /* Calculate sizes */
353: PetscDSGetSpatialDimension(prob, &dim);
354: PetscDSGetCoordinateDimension(prob, &dimEmbed);
355: prob->totDim = prob->totComp = 0;
356: PetscMalloc2(Nf,&prob->Nc,Nf,&prob->Nb);
357: PetscCalloc2(Nf+1,&prob->off,Nf+1,&prob->offDer);
358: PetscMalloc4(Nf,&prob->basis,Nf,&prob->basisDer,Nf,&prob->basisFace,Nf,&prob->basisDerFace);
359: for (f = 0; f < Nf; ++f) {
360: PetscObject obj;
361: PetscClassId id;
362: PetscQuadrature q;
363: PetscInt Nq = 0, Nb, Nc;
365: PetscDSGetDiscretization(prob, f, &obj);
366: PetscObjectGetClassId(obj, &id);
367: if (id == PETSCFE_CLASSID) {
368: PetscFE fe = (PetscFE) obj;
370: PetscFEGetQuadrature(fe, &q);
371: PetscFEGetDimension(fe, &Nb);
372: PetscFEGetNumComponents(fe, &Nc);
373: PetscFEGetDefaultTabulation(fe, &prob->basis[f], &prob->basisDer[f], NULL);
374: PetscFEGetFaceTabulation(fe, &prob->basisFace[f], &prob->basisDerFace[f], NULL);
375: } else if (id == PETSCFV_CLASSID) {
376: PetscFV fv = (PetscFV) obj;
378: PetscFVGetQuadrature(fv, &q);
379: PetscFVGetNumComponents(fv, &Nc);
380: Nb = Nc;
381: PetscFVGetDefaultTabulation(fv, &prob->basis[f], &prob->basisDer[f], NULL);
382: /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */
383: } else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", f);
384: prob->Nc[f] = Nc;
385: prob->Nb[f] = Nb;
386: prob->off[f+1] = Nc + prob->off[f];
387: prob->offDer[f+1] = Nc*dim + prob->offDer[f];
388: if (q) {PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL);}
389: NqMax = PetscMax(NqMax, Nq);
390: NbMax = PetscMax(NbMax, Nb);
391: NcMax = PetscMax(NcMax, Nc);
392: prob->totDim += Nb;
393: prob->totComp += Nc;
394: /* There are two faces for all fields but the cohesive field on a hybrid cell */
395: if (prob->isHybrid && (f < Nf-1)) prob->totDim += Nb;
396: }
397: /* Allocate works space */
398: if (prob->isHybrid) NsMax = 2;
399: PetscMalloc3(NsMax*prob->totComp,&prob->u,NsMax*prob->totComp,&prob->u_t,NsMax*prob->totComp*dimEmbed,&prob->u_x);
400: PetscMalloc5(dimEmbed,&prob->x,NbMax*NcMax,&prob->basisReal,NbMax*NcMax*dimEmbed,&prob->basisDerReal,NbMax*NcMax,&prob->testReal,NbMax*NcMax*dimEmbed,&prob->testDerReal);
401: PetscMalloc6(NsMax*NqMax*NcMax,&prob->f0,NsMax*NqMax*NcMax*dim,&prob->f1,
402: NsMax*NsMax*NqMax*NcMax*NcMax,&prob->g0,NsMax*NsMax*NqMax*NcMax*NcMax*dim,&prob->g1,
403: NsMax*NsMax*NqMax*NcMax*NcMax*dim,&prob->g2,NsMax*NsMax*NqMax*NcMax*NcMax*dim*dim,&prob->g3);
404: if (prob->ops->setup) {(*prob->ops->setup)(prob);}
405: prob->setup = PETSC_TRUE;
406: return(0);
407: }
409: static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob)
410: {
414: PetscFree2(prob->Nc,prob->Nb);
415: PetscFree2(prob->off,prob->offDer);
416: PetscFree4(prob->basis,prob->basisDer,prob->basisFace,prob->basisDerFace);
417: PetscFree3(prob->u,prob->u_t,prob->u_x);
418: PetscFree5(prob->x,prob->basisReal, prob->basisDerReal,prob->testReal,prob->testDerReal);
419: PetscFree6(prob->f0,prob->f1,prob->g0,prob->g1,prob->g2,prob->g3);
420: return(0);
421: }
423: static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew)
424: {
425: PetscObject *tmpd;
426: PetscBool *tmpi;
427: PetscPointFunc *tmpobj, *tmpf, *tmpup;
428: PetscPointJac *tmpg, *tmpgp, *tmpgt;
429: PetscBdPointFunc *tmpfbd;
430: PetscBdPointJac *tmpgbd;
431: PetscRiemannFunc *tmpr;
432: PetscSimplePointFunc *tmpexactSol;
433: void **tmpexactCtx;
434: void **tmpctx;
435: PetscInt Nf = prob->Nf, f;
436: PetscErrorCode ierr;
439: if (Nf >= NfNew) return(0);
440: prob->setup = PETSC_FALSE;
441: PetscDSDestroyStructs_Static(prob);
442: PetscMalloc2(NfNew, &tmpd, NfNew, &tmpi);
443: for (f = 0; f < Nf; ++f) {tmpd[f] = prob->disc[f]; tmpi[f] = prob->implicit[f];}
444: for (f = Nf; f < NfNew; ++f) {tmpd[f] = NULL; tmpi[f] = PETSC_TRUE;}
445: PetscFree2(prob->disc, prob->implicit);
446: prob->Nf = NfNew;
447: prob->disc = tmpd;
448: prob->implicit = tmpi;
449: PetscCalloc7(NfNew, &tmpobj, NfNew*2, &tmpf, NfNew*NfNew*4, &tmpg, NfNew*NfNew*4, &tmpgp, NfNew*NfNew*4, &tmpgt, NfNew, &tmpr, NfNew, &tmpctx);
450: PetscCalloc1(NfNew, &tmpup);
451: for (f = 0; f < Nf; ++f) tmpobj[f] = prob->obj[f];
452: for (f = 0; f < Nf*2; ++f) tmpf[f] = prob->f[f];
453: for (f = 0; f < Nf*Nf*4; ++f) tmpg[f] = prob->g[f];
454: for (f = 0; f < Nf*Nf*4; ++f) tmpgp[f] = prob->gp[f];
455: for (f = 0; f < Nf; ++f) tmpr[f] = prob->r[f];
456: for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f];
457: for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f];
458: for (f = Nf; f < NfNew; ++f) tmpobj[f] = NULL;
459: for (f = Nf*2; f < NfNew*2; ++f) tmpf[f] = NULL;
460: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpg[f] = NULL;
461: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgp[f] = NULL;
462: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgt[f] = NULL;
463: for (f = Nf; f < NfNew; ++f) tmpr[f] = NULL;
464: for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL;
465: for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL;
466: PetscFree7(prob->obj, prob->f, prob->g, prob->gp, prob->gt, prob->r, prob->ctx);
467: PetscFree(prob->update);
468: prob->obj = tmpobj;
469: prob->f = tmpf;
470: prob->g = tmpg;
471: prob->gp = tmpgp;
472: prob->gt = tmpgt;
473: prob->r = tmpr;
474: prob->update = tmpup;
475: prob->ctx = tmpctx;
476: PetscCalloc4(NfNew*2, &tmpfbd, NfNew*NfNew*4, &tmpgbd, NfNew, &tmpexactSol, NfNew, &tmpexactCtx);
477: for (f = 0; f < Nf*2; ++f) tmpfbd[f] = prob->fBd[f];
478: for (f = 0; f < Nf*Nf*4; ++f) tmpgbd[f] = prob->gBd[f];
479: for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f];
480: for (f = 0; f < Nf; ++f) tmpexactCtx[f] = prob->exactCtx[f];
481: for (f = Nf*2; f < NfNew*2; ++f) tmpfbd[f] = NULL;
482: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgbd[f] = NULL;
483: for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL;
484: for (f = Nf; f < NfNew; ++f) tmpexactCtx[f] = NULL;
485: PetscFree4(prob->fBd, prob->gBd, prob->exactSol, prob->exactCtx);
486: prob->fBd = tmpfbd;
487: prob->gBd = tmpgbd;
488: prob->exactSol = tmpexactSol;
489: prob->exactCtx = tmpexactCtx;
490: return(0);
491: }
493: /*@
494: PetscDSDestroy - Destroys a PetscDS object
496: Collective on prob
498: Input Parameter:
499: . prob - the PetscDS object to destroy
501: Level: developer
503: .seealso PetscDSView()
504: @*/
505: PetscErrorCode PetscDSDestroy(PetscDS *prob)
506: {
507: PetscInt f;
508: DSBoundary next;
512: if (!*prob) return(0);
515: if (--((PetscObject)(*prob))->refct > 0) {*prob = 0; return(0);}
516: ((PetscObject) (*prob))->refct = 0;
517: if ((*prob)->subprobs) {
518: PetscInt dim, d;
520: PetscDSGetSpatialDimension(*prob, &dim);
521: for (d = 0; d < dim; ++d) {PetscDSDestroy(&(*prob)->subprobs[d]);}
522: }
523: PetscFree((*prob)->subprobs);
524: PetscDSDestroyStructs_Static(*prob);
525: for (f = 0; f < (*prob)->Nf; ++f) {
526: PetscObjectDereference((*prob)->disc[f]);
527: }
528: PetscFree2((*prob)->disc, (*prob)->implicit);
529: PetscFree7((*prob)->obj,(*prob)->f,(*prob)->g,(*prob)->gp,(*prob)->gt,(*prob)->r,(*prob)->ctx);
530: PetscFree((*prob)->update);
531: PetscFree4((*prob)->fBd,(*prob)->gBd,(*prob)->exactSol,(*prob)->exactCtx);
532: if ((*prob)->ops->destroy) {(*(*prob)->ops->destroy)(*prob);}
533: next = (*prob)->boundary;
534: while (next) {
535: DSBoundary b = next;
537: next = b->next;
538: PetscFree(b->comps);
539: PetscFree(b->ids);
540: PetscFree(b->name);
541: PetscFree(b->labelname);
542: PetscFree(b);
543: }
544: PetscFree((*prob)->constants);
545: PetscHeaderDestroy(prob);
546: return(0);
547: }
549: /*@
550: PetscDSCreate - Creates an empty PetscDS object. The type can then be set with PetscDSSetType().
552: Collective
554: Input Parameter:
555: . comm - The communicator for the PetscDS object
557: Output Parameter:
558: . prob - The PetscDS object
560: Level: beginner
562: .seealso: PetscDSSetType(), PETSCDSBASIC
563: @*/
564: PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *prob)
565: {
566: PetscDS p;
571: *prob = NULL;
572: PetscDSInitializePackage();
574: PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView);
576: p->Nf = 0;
577: p->setup = PETSC_FALSE;
578: p->numConstants = 0;
579: p->constants = NULL;
580: p->dimEmbed = -1;
581: p->useJacPre = PETSC_TRUE;
583: *prob = p;
584: return(0);
585: }
587: /*@
588: PetscDSGetNumFields - Returns the number of fields in the DS
590: Not collective
592: Input Parameter:
593: . prob - The PetscDS object
595: Output Parameter:
596: . Nf - The number of fields
598: Level: beginner
600: .seealso: PetscDSGetSpatialDimension(), PetscDSCreate()
601: @*/
602: PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf)
603: {
607: *Nf = prob->Nf;
608: return(0);
609: }
611: /*@
612: PetscDSGetSpatialDimension - Returns the spatial dimension of the DS, meaning the topological dimension of the discretizations
614: Not collective
616: Input Parameter:
617: . prob - The PetscDS object
619: Output Parameter:
620: . dim - The spatial dimension
622: Level: beginner
624: .seealso: PetscDSGetCoordinateDimension(), PetscDSGetNumFields(), PetscDSCreate()
625: @*/
626: PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim)
627: {
633: *dim = 0;
634: if (prob->Nf) {
635: PetscObject obj;
636: PetscClassId id;
638: PetscDSGetDiscretization(prob, 0, &obj);
639: PetscObjectGetClassId(obj, &id);
640: if (id == PETSCFE_CLASSID) {PetscFEGetSpatialDimension((PetscFE) obj, dim);}
641: else if (id == PETSCFV_CLASSID) {PetscFVGetSpatialDimension((PetscFV) obj, dim);}
642: else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
643: }
644: return(0);
645: }
647: /*@
648: PetscDSGetCoordinateDimension - Returns the coordinate dimension of the DS, meaning the dimension of the space into which the discretiaztions are embedded
650: Not collective
652: Input Parameter:
653: . prob - The PetscDS object
655: Output Parameter:
656: . dimEmbed - The coordinate dimension
658: Level: beginner
660: .seealso: PetscDSSetCoordinateDimension(), PetscDSGetSpatialDimension(), PetscDSGetNumFields(), PetscDSCreate()
661: @*/
662: PetscErrorCode PetscDSGetCoordinateDimension(PetscDS prob, PetscInt *dimEmbed)
663: {
667: if (prob->dimEmbed < 0) SETERRQ(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONGSTATE, "No coordinate dimension set for this DS");
668: *dimEmbed = prob->dimEmbed;
669: return(0);
670: }
672: /*@
673: PetscDSSetCoordinateDimension - Set the coordinate dimension of the DS, meaning the dimension of the space into which the discretiaztions are embedded
675: Logically collective on prob
677: Input Parameters:
678: + prob - The PetscDS object
679: - dimEmbed - The coordinate dimension
681: Level: beginner
683: .seealso: PetscDSGetCoordinateDimension(), PetscDSGetSpatialDimension(), PetscDSGetNumFields(), PetscDSCreate()
684: @*/
685: PetscErrorCode PetscDSSetCoordinateDimension(PetscDS prob, PetscInt dimEmbed)
686: {
689: if (dimEmbed < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Coordinate dimension must be non-negative, not %D", dimEmbed);
690: prob->dimEmbed = dimEmbed;
691: return(0);
692: }
694: /*@
695: PetscDSGetHybrid - Returns the flag for a hybrid (cohesive) cell
697: Not collective
699: Input Parameter:
700: . prob - The PetscDS object
702: Output Parameter:
703: . isHybrid - The flag
705: Level: developer
707: .seealso: PetscDSSetHybrid(), PetscDSCreate()
708: @*/
709: PetscErrorCode PetscDSGetHybrid(PetscDS prob, PetscBool *isHybrid)
710: {
714: *isHybrid = prob->isHybrid;
715: return(0);
716: }
718: /*@
719: PetscDSSetHybrid - Set the flag for a hybrid (cohesive) cell
721: Not collective
723: Input Parameters:
724: + prob - The PetscDS object
725: - isHybrid - The flag
727: Level: developer
729: .seealso: PetscDSGetHybrid(), PetscDSCreate()
730: @*/
731: PetscErrorCode PetscDSSetHybrid(PetscDS prob, PetscBool isHybrid)
732: {
735: prob->isHybrid = isHybrid;
736: return(0);
737: }
739: /*@
740: PetscDSGetTotalDimension - Returns the total size of the approximation space for this system
742: Not collective
744: Input Parameter:
745: . prob - The PetscDS object
747: Output Parameter:
748: . dim - The total problem dimension
750: Level: beginner
752: .seealso: PetscDSGetNumFields(), PetscDSCreate()
753: @*/
754: PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim)
755: {
760: PetscDSSetUp(prob);
762: *dim = prob->totDim;
763: return(0);
764: }
766: /*@
767: PetscDSGetTotalComponents - Returns the total number of components in this system
769: Not collective
771: Input Parameter:
772: . prob - The PetscDS object
774: Output Parameter:
775: . dim - The total number of components
777: Level: beginner
779: .seealso: PetscDSGetNumFields(), PetscDSCreate()
780: @*/
781: PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc)
782: {
787: PetscDSSetUp(prob);
789: *Nc = prob->totComp;
790: return(0);
791: }
793: /*@
794: PetscDSGetDiscretization - Returns the discretization object for the given field
796: Not collective
798: Input Parameters:
799: + prob - The PetscDS object
800: - f - The field number
802: Output Parameter:
803: . disc - The discretization object
805: Level: beginner
807: .seealso: PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
808: @*/
809: PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc)
810: {
814: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
815: *disc = prob->disc[f];
816: return(0);
817: }
819: /*@
820: PetscDSSetDiscretization - Sets the discretization object for the given field
822: Not collective
824: Input Parameters:
825: + prob - The PetscDS object
826: . f - The field number
827: - disc - The discretization object
829: Level: beginner
831: .seealso: PetscDSGetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
832: @*/
833: PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc)
834: {
840: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
841: PetscDSEnlarge_Static(prob, f+1);
842: PetscObjectDereference(prob->disc[f]);
843: prob->disc[f] = disc;
844: PetscObjectReference(disc);
845: {
846: PetscClassId id;
848: PetscObjectGetClassId(disc, &id);
849: if (id == PETSCFE_CLASSID) {
850: PetscDSSetImplicit(prob, f, PETSC_TRUE);
851: } else if (id == PETSCFV_CLASSID) {
852: PetscDSSetImplicit(prob, f, PETSC_FALSE);
853: }
854: }
855: return(0);
856: }
858: /*@
859: PetscDSAddDiscretization - Adds a discretization object
861: Not collective
863: Input Parameters:
864: + prob - The PetscDS object
865: - disc - The boundary discretization object
867: Level: beginner
869: .seealso: PetscDSGetDiscretization(), PetscDSSetDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
870: @*/
871: PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc)
872: {
876: PetscDSSetDiscretization(prob, prob->Nf, disc);
877: return(0);
878: }
880: /*@
881: PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for IMEX
883: Not collective
885: Input Parameters:
886: + prob - The PetscDS object
887: - f - The field number
889: Output Parameter:
890: . implicit - The flag indicating what kind of solve to use for this field
892: Level: developer
894: .seealso: PetscDSSetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
895: @*/
896: PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit)
897: {
901: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
902: *implicit = prob->implicit[f];
903: return(0);
904: }
906: /*@
907: PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for IMEX
909: Not collective
911: Input Parameters:
912: + prob - The PetscDS object
913: . f - The field number
914: - implicit - The flag indicating what kind of solve to use for this field
916: Level: developer
918: .seealso: PetscDSGetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
919: @*/
920: PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit)
921: {
924: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
925: prob->implicit[f] = implicit;
926: return(0);
927: }
929: PetscErrorCode PetscDSGetObjective(PetscDS prob, PetscInt f,
930: void (**obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
931: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
932: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
933: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
934: {
938: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
939: *obj = prob->obj[f];
940: return(0);
941: }
943: PetscErrorCode PetscDSSetObjective(PetscDS prob, PetscInt f,
944: void (*obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
945: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
946: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
947: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
948: {
954: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
955: PetscDSEnlarge_Static(prob, f+1);
956: prob->obj[f] = obj;
957: return(0);
958: }
960: /*@C
961: PetscDSGetResidual - Get the pointwise residual function for a given test field
963: Not collective
965: Input Parameters:
966: + prob - The PetscDS
967: - f - The test field number
969: Output Parameters:
970: + f0 - integrand for the test function term
971: - f1 - integrand for the test function gradient term
973: Note: We are using a first order FEM model for the weak form:
975: \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)
977: The calling sequence for the callbacks f0 and f1 is given by:
979: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
980: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
981: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
982: $ PetscReal t, const PetscReal x[], PetscScalar f0[])
984: + dim - the spatial dimension
985: . Nf - the number of fields
986: . uOff - the offset into u[] and u_t[] for each field
987: . uOff_x - the offset into u_x[] for each field
988: . u - each field evaluated at the current point
989: . u_t - the time derivative of each field evaluated at the current point
990: . u_x - the gradient of each field evaluated at the current point
991: . aOff - the offset into a[] and a_t[] for each auxiliary field
992: . aOff_x - the offset into a_x[] for each auxiliary field
993: . a - each auxiliary field evaluated at the current point
994: . a_t - the time derivative of each auxiliary field evaluated at the current point
995: . a_x - the gradient of auxiliary each field evaluated at the current point
996: . t - current time
997: . x - coordinates of the current point
998: . numConstants - number of constant parameters
999: . constants - constant parameters
1000: - f0 - output values at the current point
1002: Level: intermediate
1004: .seealso: PetscDSSetResidual()
1005: @*/
1006: PetscErrorCode PetscDSGetResidual(PetscDS prob, PetscInt f,
1007: void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1008: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1009: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1010: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1011: void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1012: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1013: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1014: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1015: {
1018: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1021: return(0);
1022: }
1024: /*@C
1025: PetscDSSetResidual - Set the pointwise residual function for a given test field
1027: Not collective
1029: Input Parameters:
1030: + prob - The PetscDS
1031: . f - The test field number
1032: . f0 - integrand for the test function term
1033: - f1 - integrand for the test function gradient term
1035: Note: We are using a first order FEM model for the weak form:
1037: \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)
1039: The calling sequence for the callbacks f0 and f1 is given by:
1041: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1042: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1043: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1044: $ PetscReal t, const PetscReal x[], PetscScalar f0[])
1046: + dim - the spatial dimension
1047: . Nf - the number of fields
1048: . uOff - the offset into u[] and u_t[] for each field
1049: . uOff_x - the offset into u_x[] for each field
1050: . u - each field evaluated at the current point
1051: . u_t - the time derivative of each field evaluated at the current point
1052: . u_x - the gradient of each field evaluated at the current point
1053: . aOff - the offset into a[] and a_t[] for each auxiliary field
1054: . aOff_x - the offset into a_x[] for each auxiliary field
1055: . a - each auxiliary field evaluated at the current point
1056: . a_t - the time derivative of each auxiliary field evaluated at the current point
1057: . a_x - the gradient of auxiliary each field evaluated at the current point
1058: . t - current time
1059: . x - coordinates of the current point
1060: . numConstants - number of constant parameters
1061: . constants - constant parameters
1062: - f0 - output values at the current point
1064: Level: intermediate
1066: .seealso: PetscDSGetResidual()
1067: @*/
1068: PetscErrorCode PetscDSSetResidual(PetscDS prob, PetscInt f,
1069: void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1070: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1071: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1072: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1073: void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1074: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1075: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1076: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1077: {
1084: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1085: PetscDSEnlarge_Static(prob, f+1);
1086: prob->f[f*2+0] = f0;
1087: prob->f[f*2+1] = f1;
1088: return(0);
1089: }
1091: /*@C
1092: PetscDSHasJacobian - Signals that Jacobian functions have been set
1094: Not collective
1096: Input Parameter:
1097: . prob - The PetscDS
1099: Output Parameter:
1100: . hasJac - flag that pointwise function for the Jacobian has been set
1102: Level: intermediate
1104: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1105: @*/
1106: PetscErrorCode PetscDSHasJacobian(PetscDS prob, PetscBool *hasJac)
1107: {
1108: PetscInt f, g, h;
1112: *hasJac = PETSC_FALSE;
1113: for (f = 0; f < prob->Nf; ++f) {
1114: for (g = 0; g < prob->Nf; ++g) {
1115: for (h = 0; h < 4; ++h) {
1116: if (prob->g[(f*prob->Nf + g)*4+h]) *hasJac = PETSC_TRUE;
1117: }
1118: }
1119: }
1120: return(0);
1121: }
1123: /*@C
1124: PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field
1126: Not collective
1128: Input Parameters:
1129: + prob - The PetscDS
1130: . f - The test field number
1131: - g - The field number
1133: Output Parameters:
1134: + g0 - integrand for the test and basis function term
1135: . g1 - integrand for the test function and basis function gradient term
1136: . g2 - integrand for the test function gradient and basis function term
1137: - g3 - integrand for the test function gradient and basis function gradient term
1139: Note: We are using a first order FEM model for the weak form:
1141: \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1143: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1145: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1146: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1147: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1148: $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1150: + dim - the spatial dimension
1151: . Nf - the number of fields
1152: . uOff - the offset into u[] and u_t[] for each field
1153: . uOff_x - the offset into u_x[] for each field
1154: . u - each field evaluated at the current point
1155: . u_t - the time derivative of each field evaluated at the current point
1156: . u_x - the gradient of each field evaluated at the current point
1157: . aOff - the offset into a[] and a_t[] for each auxiliary field
1158: . aOff_x - the offset into a_x[] for each auxiliary field
1159: . a - each auxiliary field evaluated at the current point
1160: . a_t - the time derivative of each auxiliary field evaluated at the current point
1161: . a_x - the gradient of auxiliary each field evaluated at the current point
1162: . t - current time
1163: . u_tShift - the multiplier a for dF/dU_t
1164: . x - coordinates of the current point
1165: . numConstants - number of constant parameters
1166: . constants - constant parameters
1167: - g0 - output values at the current point
1169: Level: intermediate
1171: .seealso: PetscDSSetJacobian()
1172: @*/
1173: PetscErrorCode PetscDSGetJacobian(PetscDS prob, PetscInt f, PetscInt g,
1174: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1175: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1176: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1177: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1178: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1179: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1180: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1181: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1182: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1183: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1184: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1185: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1186: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1187: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1188: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1189: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1190: {
1193: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1194: if ((g < 0) || (g >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", g, prob->Nf);
1199: return(0);
1200: }
1202: /*@C
1203: PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields
1205: Not collective
1207: Input Parameters:
1208: + prob - The PetscDS
1209: . f - The test field number
1210: . g - The field number
1211: . g0 - integrand for the test and basis function term
1212: . g1 - integrand for the test function and basis function gradient term
1213: . g2 - integrand for the test function gradient and basis function term
1214: - g3 - integrand for the test function gradient and basis function gradient term
1216: Note: We are using a first order FEM model for the weak form:
1218: \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1220: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1222: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1223: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1224: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1225: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1227: + dim - the spatial dimension
1228: . Nf - the number of fields
1229: . uOff - the offset into u[] and u_t[] for each field
1230: . uOff_x - the offset into u_x[] for each field
1231: . u - each field evaluated at the current point
1232: . u_t - the time derivative of each field evaluated at the current point
1233: . u_x - the gradient of each field evaluated at the current point
1234: . aOff - the offset into a[] and a_t[] for each auxiliary field
1235: . aOff_x - the offset into a_x[] for each auxiliary field
1236: . a - each auxiliary field evaluated at the current point
1237: . a_t - the time derivative of each auxiliary field evaluated at the current point
1238: . a_x - the gradient of auxiliary each field evaluated at the current point
1239: . t - current time
1240: . u_tShift - the multiplier a for dF/dU_t
1241: . x - coordinates of the current point
1242: . numConstants - number of constant parameters
1243: . constants - constant parameters
1244: - g0 - output values at the current point
1246: Level: intermediate
1248: .seealso: PetscDSGetJacobian()
1249: @*/
1250: PetscErrorCode PetscDSSetJacobian(PetscDS prob, PetscInt f, PetscInt g,
1251: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1252: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1253: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1254: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1255: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1256: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1257: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1258: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1259: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1260: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1261: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1262: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1263: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1264: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1265: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1266: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1267: {
1276: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1277: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1278: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1279: prob->g[(f*prob->Nf + g)*4+0] = g0;
1280: prob->g[(f*prob->Nf + g)*4+1] = g1;
1281: prob->g[(f*prob->Nf + g)*4+2] = g2;
1282: prob->g[(f*prob->Nf + g)*4+3] = g3;
1283: return(0);
1284: }
1286: /*@C
1287: PetscDSUseJacobianPreconditioner - Whether to construct a Jacobian preconditioner
1289: Not collective
1291: Input Parameters:
1292: + prob - The PetscDS
1293: - useJacPre - flag that enables construction of a Jacobian preconditioner
1295: Level: intermediate
1297: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1298: @*/
1299: PetscErrorCode PetscDSUseJacobianPreconditioner(PetscDS prob, PetscBool useJacPre)
1300: {
1303: prob->useJacPre = useJacPre;
1304: return(0);
1305: }
1307: /*@C
1308: PetscDSHasJacobianPreconditioner - Signals that a Jacobian preconditioner matrix has been set
1310: Not collective
1312: Input Parameter:
1313: . prob - The PetscDS
1315: Output Parameter:
1316: . hasJacPre - flag that pointwise function for Jacobian preconditioner matrix has been set
1318: Level: intermediate
1320: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1321: @*/
1322: PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS prob, PetscBool *hasJacPre)
1323: {
1324: PetscInt f, g, h;
1328: *hasJacPre = PETSC_FALSE;
1329: if (!prob->useJacPre) return(0);
1330: for (f = 0; f < prob->Nf; ++f) {
1331: for (g = 0; g < prob->Nf; ++g) {
1332: for (h = 0; h < 4; ++h) {
1333: if (prob->gp[(f*prob->Nf + g)*4+h]) *hasJacPre = PETSC_TRUE;
1334: }
1335: }
1336: }
1337: return(0);
1338: }
1340: /*@C
1341: PetscDSGetJacobianPreconditioner - Get the pointwise Jacobian preconditioner function for given test and basis field. If this is missing, the system matrix is used to build the preconditioner.
1343: Not collective
1345: Input Parameters:
1346: + prob - The PetscDS
1347: . f - The test field number
1348: - g - The field number
1350: Output Parameters:
1351: + g0 - integrand for the test and basis function term
1352: . g1 - integrand for the test function and basis function gradient term
1353: . g2 - integrand for the test function gradient and basis function term
1354: - g3 - integrand for the test function gradient and basis function gradient term
1356: Note: We are using a first order FEM model for the weak form:
1358: \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1360: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1362: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1363: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1364: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1365: $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1367: + dim - the spatial dimension
1368: . Nf - the number of fields
1369: . uOff - the offset into u[] and u_t[] for each field
1370: . uOff_x - the offset into u_x[] for each field
1371: . u - each field evaluated at the current point
1372: . u_t - the time derivative of each field evaluated at the current point
1373: . u_x - the gradient of each field evaluated at the current point
1374: . aOff - the offset into a[] and a_t[] for each auxiliary field
1375: . aOff_x - the offset into a_x[] for each auxiliary field
1376: . a - each auxiliary field evaluated at the current point
1377: . a_t - the time derivative of each auxiliary field evaluated at the current point
1378: . a_x - the gradient of auxiliary each field evaluated at the current point
1379: . t - current time
1380: . u_tShift - the multiplier a for dF/dU_t
1381: . x - coordinates of the current point
1382: . numConstants - number of constant parameters
1383: . constants - constant parameters
1384: - g0 - output values at the current point
1386: Level: intermediate
1388: .seealso: PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1389: @*/
1390: PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g,
1391: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1392: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1393: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1394: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1395: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1396: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1397: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1398: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1399: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1400: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1401: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1402: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1403: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1404: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1405: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1406: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1407: {
1410: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1411: if ((g < 0) || (g >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", g, prob->Nf);
1416: return(0);
1417: }
1419: /*@C
1420: PetscDSSetJacobianPreconditioner - Set the pointwise Jacobian preconditioner function for given test and basis fields. If this is missing, the system matrix is used to build the preconditioner.
1422: Not collective
1424: Input Parameters:
1425: + prob - The PetscDS
1426: . f - The test field number
1427: . g - The field number
1428: . g0 - integrand for the test and basis function term
1429: . g1 - integrand for the test function and basis function gradient term
1430: . g2 - integrand for the test function gradient and basis function term
1431: - g3 - integrand for the test function gradient and basis function gradient term
1433: Note: We are using a first order FEM model for the weak form:
1435: \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1437: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1439: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1440: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1441: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1442: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1444: + dim - the spatial dimension
1445: . Nf - the number of fields
1446: . uOff - the offset into u[] and u_t[] for each field
1447: . uOff_x - the offset into u_x[] for each field
1448: . u - each field evaluated at the current point
1449: . u_t - the time derivative of each field evaluated at the current point
1450: . u_x - the gradient of each field evaluated at the current point
1451: . aOff - the offset into a[] and a_t[] for each auxiliary field
1452: . aOff_x - the offset into a_x[] for each auxiliary field
1453: . a - each auxiliary field evaluated at the current point
1454: . a_t - the time derivative of each auxiliary field evaluated at the current point
1455: . a_x - the gradient of auxiliary each field evaluated at the current point
1456: . t - current time
1457: . u_tShift - the multiplier a for dF/dU_t
1458: . x - coordinates of the current point
1459: . numConstants - number of constant parameters
1460: . constants - constant parameters
1461: - g0 - output values at the current point
1463: Level: intermediate
1465: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobian()
1466: @*/
1467: PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g,
1468: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1469: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1470: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1471: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1472: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1473: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1474: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1475: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1476: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1477: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1478: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1479: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1480: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1481: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1482: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1483: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1484: {
1493: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1494: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1495: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1496: prob->gp[(f*prob->Nf + g)*4+0] = g0;
1497: prob->gp[(f*prob->Nf + g)*4+1] = g1;
1498: prob->gp[(f*prob->Nf + g)*4+2] = g2;
1499: prob->gp[(f*prob->Nf + g)*4+3] = g3;
1500: return(0);
1501: }
1503: /*@C
1504: PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, dF/du_t, has been set
1506: Not collective
1508: Input Parameter:
1509: . prob - The PetscDS
1511: Output Parameter:
1512: . hasDynJac - flag that pointwise function for dynamic Jacobian has been set
1514: Level: intermediate
1516: .seealso: PetscDSGetDynamicJacobian(), PetscDSSetDynamicJacobian(), PetscDSGetJacobian()
1517: @*/
1518: PetscErrorCode PetscDSHasDynamicJacobian(PetscDS prob, PetscBool *hasDynJac)
1519: {
1520: PetscInt f, g, h;
1524: *hasDynJac = PETSC_FALSE;
1525: for (f = 0; f < prob->Nf; ++f) {
1526: for (g = 0; g < prob->Nf; ++g) {
1527: for (h = 0; h < 4; ++h) {
1528: if (prob->gt[(f*prob->Nf + g)*4+h]) *hasDynJac = PETSC_TRUE;
1529: }
1530: }
1531: }
1532: return(0);
1533: }
1535: /*@C
1536: PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, dF/du_t, function for given test and basis field
1538: Not collective
1540: Input Parameters:
1541: + prob - The PetscDS
1542: . f - The test field number
1543: - g - The field number
1545: Output Parameters:
1546: + g0 - integrand for the test and basis function term
1547: . g1 - integrand for the test function and basis function gradient term
1548: . g2 - integrand for the test function gradient and basis function term
1549: - g3 - integrand for the test function gradient and basis function gradient term
1551: Note: We are using a first order FEM model for the weak form:
1553: \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1555: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1557: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1558: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1559: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1560: $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1562: + dim - the spatial dimension
1563: . Nf - the number of fields
1564: . uOff - the offset into u[] and u_t[] for each field
1565: . uOff_x - the offset into u_x[] for each field
1566: . u - each field evaluated at the current point
1567: . u_t - the time derivative of each field evaluated at the current point
1568: . u_x - the gradient of each field evaluated at the current point
1569: . aOff - the offset into a[] and a_t[] for each auxiliary field
1570: . aOff_x - the offset into a_x[] for each auxiliary field
1571: . a - each auxiliary field evaluated at the current point
1572: . a_t - the time derivative of each auxiliary field evaluated at the current point
1573: . a_x - the gradient of auxiliary each field evaluated at the current point
1574: . t - current time
1575: . u_tShift - the multiplier a for dF/dU_t
1576: . x - coordinates of the current point
1577: . numConstants - number of constant parameters
1578: . constants - constant parameters
1579: - g0 - output values at the current point
1581: Level: intermediate
1583: .seealso: PetscDSSetJacobian()
1584: @*/
1585: PetscErrorCode PetscDSGetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g,
1586: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1587: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1588: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1589: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1590: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1591: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1592: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1593: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1594: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1595: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1596: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1597: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1598: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1599: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1600: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1601: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1602: {
1605: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1606: if ((g < 0) || (g >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", g, prob->Nf);
1611: return(0);
1612: }
1614: /*@C
1615: PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, dF/du_t, function for given test and basis fields
1617: Not collective
1619: Input Parameters:
1620: + prob - The PetscDS
1621: . f - The test field number
1622: . g - The field number
1623: . g0 - integrand for the test and basis function term
1624: . g1 - integrand for the test function and basis function gradient term
1625: . g2 - integrand for the test function gradient and basis function term
1626: - g3 - integrand for the test function gradient and basis function gradient term
1628: Note: We are using a first order FEM model for the weak form:
1630: \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1632: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1634: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1635: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1636: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1637: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1639: + dim - the spatial dimension
1640: . Nf - the number of fields
1641: . uOff - the offset into u[] and u_t[] for each field
1642: . uOff_x - the offset into u_x[] for each field
1643: . u - each field evaluated at the current point
1644: . u_t - the time derivative of each field evaluated at the current point
1645: . u_x - the gradient of each field evaluated at the current point
1646: . aOff - the offset into a[] and a_t[] for each auxiliary field
1647: . aOff_x - the offset into a_x[] for each auxiliary field
1648: . a - each auxiliary field evaluated at the current point
1649: . a_t - the time derivative of each auxiliary field evaluated at the current point
1650: . a_x - the gradient of auxiliary each field evaluated at the current point
1651: . t - current time
1652: . u_tShift - the multiplier a for dF/dU_t
1653: . x - coordinates of the current point
1654: . numConstants - number of constant parameters
1655: . constants - constant parameters
1656: - g0 - output values at the current point
1658: Level: intermediate
1660: .seealso: PetscDSGetJacobian()
1661: @*/
1662: PetscErrorCode PetscDSSetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g,
1663: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1664: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1665: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1666: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1667: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1668: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1669: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1670: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1671: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1672: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1673: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1674: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1675: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1676: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1677: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1678: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1679: {
1688: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1689: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1690: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1691: prob->gt[(f*prob->Nf + g)*4+0] = g0;
1692: prob->gt[(f*prob->Nf + g)*4+1] = g1;
1693: prob->gt[(f*prob->Nf + g)*4+2] = g2;
1694: prob->gt[(f*prob->Nf + g)*4+3] = g3;
1695: return(0);
1696: }
1698: /*@C
1699: PetscDSGetRiemannSolver - Returns the Riemann solver for the given field
1701: Not collective
1703: Input Arguments:
1704: + prob - The PetscDS object
1705: - f - The field number
1707: Output Argument:
1708: . r - Riemann solver
1710: Calling sequence for r:
1712: $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)
1714: + dim - The spatial dimension
1715: . Nf - The number of fields
1716: . x - The coordinates at a point on the interface
1717: . n - The normal vector to the interface
1718: . uL - The state vector to the left of the interface
1719: . uR - The state vector to the right of the interface
1720: . flux - output array of flux through the interface
1721: . numConstants - number of constant parameters
1722: . constants - constant parameters
1723: - ctx - optional user context
1725: Level: intermediate
1727: .seealso: PetscDSSetRiemannSolver()
1728: @*/
1729: PetscErrorCode PetscDSGetRiemannSolver(PetscDS prob, PetscInt f,
1730: void (**r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscInt numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx))
1731: {
1734: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1736: *r = prob->r[f];
1737: return(0);
1738: }
1740: /*@C
1741: PetscDSSetRiemannSolver - Sets the Riemann solver for the given field
1743: Not collective
1745: Input Arguments:
1746: + prob - The PetscDS object
1747: . f - The field number
1748: - r - Riemann solver
1750: Calling sequence for r:
1752: $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)
1754: + dim - The spatial dimension
1755: . Nf - The number of fields
1756: . x - The coordinates at a point on the interface
1757: . n - The normal vector to the interface
1758: . uL - The state vector to the left of the interface
1759: . uR - The state vector to the right of the interface
1760: . flux - output array of flux through the interface
1761: . numConstants - number of constant parameters
1762: . constants - constant parameters
1763: - ctx - optional user context
1765: Level: intermediate
1767: .seealso: PetscDSGetRiemannSolver()
1768: @*/
1769: PetscErrorCode PetscDSSetRiemannSolver(PetscDS prob, PetscInt f,
1770: void (*r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscInt numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx))
1771: {
1777: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1778: PetscDSEnlarge_Static(prob, f+1);
1779: prob->r[f] = r;
1780: return(0);
1781: }
1783: /*@C
1784: PetscDSGetUpdate - Get the pointwise update function for a given field
1786: Not collective
1788: Input Parameters:
1789: + prob - The PetscDS
1790: - f - The field number
1792: Output Parameters:
1793: . update - update function
1795: Note: The calling sequence for the callback update is given by:
1797: $ update(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1798: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1799: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1800: $ PetscReal t, const PetscReal x[], PetscScalar uNew[])
1802: + dim - the spatial dimension
1803: . Nf - the number of fields
1804: . uOff - the offset into u[] and u_t[] for each field
1805: . uOff_x - the offset into u_x[] for each field
1806: . u - each field evaluated at the current point
1807: . u_t - the time derivative of each field evaluated at the current point
1808: . u_x - the gradient of each field evaluated at the current point
1809: . aOff - the offset into a[] and a_t[] for each auxiliary field
1810: . aOff_x - the offset into a_x[] for each auxiliary field
1811: . a - each auxiliary field evaluated at the current point
1812: . a_t - the time derivative of each auxiliary field evaluated at the current point
1813: . a_x - the gradient of auxiliary each field evaluated at the current point
1814: . t - current time
1815: . x - coordinates of the current point
1816: - uNew - new value for field at the current point
1818: Level: intermediate
1820: .seealso: PetscDSSetUpdate(), PetscDSSetResidual()
1821: @*/
1822: PetscErrorCode PetscDSGetUpdate(PetscDS prob, PetscInt f,
1823: void (**update)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1824: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1825: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1826: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
1827: {
1830: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1832: return(0);
1833: }
1835: /*@C
1836: PetscDSSetUpdate - Set the pointwise update function for a given field
1838: Not collective
1840: Input Parameters:
1841: + prob - The PetscDS
1842: . f - The field number
1843: - update - update function
1845: Note: The calling sequence for the callback update is given by:
1847: $ update(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1848: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1849: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1850: $ PetscReal t, const PetscReal x[], PetscScalar uNew[])
1852: + dim - the spatial dimension
1853: . Nf - the number of fields
1854: . uOff - the offset into u[] and u_t[] for each field
1855: . uOff_x - the offset into u_x[] for each field
1856: . u - each field evaluated at the current point
1857: . u_t - the time derivative of each field evaluated at the current point
1858: . u_x - the gradient of each field evaluated at the current point
1859: . aOff - the offset into a[] and a_t[] for each auxiliary field
1860: . aOff_x - the offset into a_x[] for each auxiliary field
1861: . a - each auxiliary field evaluated at the current point
1862: . a_t - the time derivative of each auxiliary field evaluated at the current point
1863: . a_x - the gradient of auxiliary each field evaluated at the current point
1864: . t - current time
1865: . x - coordinates of the current point
1866: - uNew - new field values at the current point
1868: Level: intermediate
1870: .seealso: PetscDSGetResidual()
1871: @*/
1872: PetscErrorCode PetscDSSetUpdate(PetscDS prob, PetscInt f,
1873: void (*update)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1874: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1875: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1876: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
1877: {
1883: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1884: PetscDSEnlarge_Static(prob, f+1);
1885: prob->update[f] = update;
1886: return(0);
1887: }
1889: PetscErrorCode PetscDSGetContext(PetscDS prob, PetscInt f, void **ctx)
1890: {
1893: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1895: *ctx = prob->ctx[f];
1896: return(0);
1897: }
1899: PetscErrorCode PetscDSSetContext(PetscDS prob, PetscInt f, void *ctx)
1900: {
1905: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1906: PetscDSEnlarge_Static(prob, f+1);
1907: prob->ctx[f] = ctx;
1908: return(0);
1909: }
1911: /*@C
1912: PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field
1914: Not collective
1916: Input Parameters:
1917: + prob - The PetscDS
1918: - f - The test field number
1920: Output Parameters:
1921: + f0 - boundary integrand for the test function term
1922: - f1 - boundary integrand for the test function gradient term
1924: Note: We are using a first order FEM model for the weak form:
1926: \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
1928: The calling sequence for the callbacks f0 and f1 is given by:
1930: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1931: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1932: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1933: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[])
1935: + dim - the spatial dimension
1936: . Nf - the number of fields
1937: . uOff - the offset into u[] and u_t[] for each field
1938: . uOff_x - the offset into u_x[] for each field
1939: . u - each field evaluated at the current point
1940: . u_t - the time derivative of each field evaluated at the current point
1941: . u_x - the gradient of each field evaluated at the current point
1942: . aOff - the offset into a[] and a_t[] for each auxiliary field
1943: . aOff_x - the offset into a_x[] for each auxiliary field
1944: . a - each auxiliary field evaluated at the current point
1945: . a_t - the time derivative of each auxiliary field evaluated at the current point
1946: . a_x - the gradient of auxiliary each field evaluated at the current point
1947: . t - current time
1948: . x - coordinates of the current point
1949: . n - unit normal at the current point
1950: . numConstants - number of constant parameters
1951: . constants - constant parameters
1952: - f0 - output values at the current point
1954: Level: intermediate
1956: .seealso: PetscDSSetBdResidual()
1957: @*/
1958: PetscErrorCode PetscDSGetBdResidual(PetscDS prob, PetscInt f,
1959: void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1960: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1961: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1962: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1963: void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1964: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1965: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1966: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1967: {
1970: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1973: return(0);
1974: }
1976: /*@C
1977: PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field
1979: Not collective
1981: Input Parameters:
1982: + prob - The PetscDS
1983: . f - The test field number
1984: . f0 - boundary integrand for the test function term
1985: - f1 - boundary integrand for the test function gradient term
1987: Note: We are using a first order FEM model for the weak form:
1989: \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
1991: The calling sequence for the callbacks f0 and f1 is given by:
1993: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1994: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1995: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1996: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[])
1998: + dim - the spatial dimension
1999: . Nf - the number of fields
2000: . uOff - the offset into u[] and u_t[] for each field
2001: . uOff_x - the offset into u_x[] for each field
2002: . u - each field evaluated at the current point
2003: . u_t - the time derivative of each field evaluated at the current point
2004: . u_x - the gradient of each field evaluated at the current point
2005: . aOff - the offset into a[] and a_t[] for each auxiliary field
2006: . aOff_x - the offset into a_x[] for each auxiliary field
2007: . a - each auxiliary field evaluated at the current point
2008: . a_t - the time derivative of each auxiliary field evaluated at the current point
2009: . a_x - the gradient of auxiliary each field evaluated at the current point
2010: . t - current time
2011: . x - coordinates of the current point
2012: . n - unit normal at the current point
2013: . numConstants - number of constant parameters
2014: . constants - constant parameters
2015: - f0 - output values at the current point
2017: Level: intermediate
2019: .seealso: PetscDSGetBdResidual()
2020: @*/
2021: PetscErrorCode PetscDSSetBdResidual(PetscDS prob, PetscInt f,
2022: void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2023: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2024: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2025: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
2026: void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2027: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2028: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2029: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
2030: {
2035: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2036: PetscDSEnlarge_Static(prob, f+1);
2039: return(0);
2040: }
2042: /*@C
2043: PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field
2045: Not collective
2047: Input Parameters:
2048: + prob - The PetscDS
2049: . f - The test field number
2050: - g - The field number
2052: Output Parameters:
2053: + g0 - integrand for the test and basis function term
2054: . g1 - integrand for the test function and basis function gradient term
2055: . g2 - integrand for the test function gradient and basis function term
2056: - g3 - integrand for the test function gradient and basis function gradient term
2058: Note: We are using a first order FEM model for the weak form:
2060: \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi
2062: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
2064: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2065: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2066: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2067: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[])
2069: + dim - the spatial dimension
2070: . Nf - the number of fields
2071: . uOff - the offset into u[] and u_t[] for each field
2072: . uOff_x - the offset into u_x[] for each field
2073: . u - each field evaluated at the current point
2074: . u_t - the time derivative of each field evaluated at the current point
2075: . u_x - the gradient of each field evaluated at the current point
2076: . aOff - the offset into a[] and a_t[] for each auxiliary field
2077: . aOff_x - the offset into a_x[] for each auxiliary field
2078: . a - each auxiliary field evaluated at the current point
2079: . a_t - the time derivative of each auxiliary field evaluated at the current point
2080: . a_x - the gradient of auxiliary each field evaluated at the current point
2081: . t - current time
2082: . u_tShift - the multiplier a for dF/dU_t
2083: . x - coordinates of the current point
2084: . n - normal at the current point
2085: . numConstants - number of constant parameters
2086: . constants - constant parameters
2087: - g0 - output values at the current point
2089: Level: intermediate
2091: .seealso: PetscDSSetBdJacobian()
2092: @*/
2093: PetscErrorCode PetscDSGetBdJacobian(PetscDS prob, PetscInt f, PetscInt g,
2094: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2095: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2096: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2097: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
2098: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2099: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2100: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2101: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
2102: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2103: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2104: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2105: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
2106: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2107: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2108: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2109: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2110: {
2113: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
2114: if ((g < 0) || (g >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", g, prob->Nf);
2119: return(0);
2120: }
2122: /*@C
2123: PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field
2125: Not collective
2127: Input Parameters:
2128: + prob - The PetscDS
2129: . f - The test field number
2130: . g - The field number
2131: . g0 - integrand for the test and basis function term
2132: . g1 - integrand for the test function and basis function gradient term
2133: . g2 - integrand for the test function gradient and basis function term
2134: - g3 - integrand for the test function gradient and basis function gradient term
2136: Note: We are using a first order FEM model for the weak form:
2138: \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi
2140: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
2142: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2143: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2144: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2145: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[])
2147: + dim - the spatial dimension
2148: . Nf - the number of fields
2149: . uOff - the offset into u[] and u_t[] for each field
2150: . uOff_x - the offset into u_x[] for each field
2151: . u - each field evaluated at the current point
2152: . u_t - the time derivative of each field evaluated at the current point
2153: . u_x - the gradient of each field evaluated at the current point
2154: . aOff - the offset into a[] and a_t[] for each auxiliary field
2155: . aOff_x - the offset into a_x[] for each auxiliary field
2156: . a - each auxiliary field evaluated at the current point
2157: . a_t - the time derivative of each auxiliary field evaluated at the current point
2158: . a_x - the gradient of auxiliary each field evaluated at the current point
2159: . t - current time
2160: . u_tShift - the multiplier a for dF/dU_t
2161: . x - coordinates of the current point
2162: . n - normal at the current point
2163: . numConstants - number of constant parameters
2164: . constants - constant parameters
2165: - g0 - output values at the current point
2167: Level: intermediate
2169: .seealso: PetscDSGetBdJacobian()
2170: @*/
2171: PetscErrorCode PetscDSSetBdJacobian(PetscDS prob, PetscInt f, PetscInt g,
2172: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2173: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2174: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2175: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
2176: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2177: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2178: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2179: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
2180: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2181: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2182: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2183: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
2184: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2185: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2186: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2187: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2188: {
2197: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2198: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
2199: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
2200: prob->gBd[(f*prob->Nf + g)*4+0] = g0;
2201: prob->gBd[(f*prob->Nf + g)*4+1] = g1;
2202: prob->gBd[(f*prob->Nf + g)*4+2] = g2;
2203: prob->gBd[(f*prob->Nf + g)*4+3] = g3;
2204: return(0);
2205: }
2207: /*@C
2208: PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field
2210: Not collective
2212: Input Parameters:
2213: + prob - The PetscDS
2214: - f - The test field number
2216: Output Parameter:
2217: + exactSol - exact solution for the test field
2218: - exactCtx - exact solution context
2220: Note: The calling sequence for the solution functions is given by:
2222: $ sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2224: + dim - the spatial dimension
2225: . t - current time
2226: . x - coordinates of the current point
2227: . Nc - the number of field components
2228: . u - the solution field evaluated at the current point
2229: - ctx - a user context
2231: Level: intermediate
2233: .seealso: PetscDSSetExactSolution()
2234: @*/
2235: PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx)
2236: {
2239: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
2242: return(0);
2243: }
2245: /*@C
2246: PetscDSSetExactSolution - Set the pointwise exact solution function for a given test field
2248: Not collective
2250: Input Parameters:
2251: + prob - The PetscDS
2252: . f - The test field number
2253: . sol - solution function for the test fields
2254: - ctx - solution context or NULL
2256: Note: The calling sequence for solution functions is given by:
2258: $ sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2260: + dim - the spatial dimension
2261: . t - current time
2262: . x - coordinates of the current point
2263: . Nc - the number of field components
2264: . u - the solution field evaluated at the current point
2265: - ctx - a user context
2267: Level: intermediate
2269: .seealso: PetscDSGetExactSolution()
2270: @*/
2271: PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx)
2272: {
2277: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2278: PetscDSEnlarge_Static(prob, f+1);
2281: return(0);
2282: }
2284: /*@C
2285: PetscDSGetConstants - Returns the array of constants passed to point functions
2287: Not collective
2289: Input Parameter:
2290: . prob - The PetscDS object
2292: Output Parameters:
2293: + numConstants - The number of constants
2294: - constants - The array of constants, NULL if there are none
2296: Level: intermediate
2298: .seealso: PetscDSSetConstants(), PetscDSCreate()
2299: @*/
2300: PetscErrorCode PetscDSGetConstants(PetscDS prob, PetscInt *numConstants, const PetscScalar *constants[])
2301: {
2306: return(0);
2307: }
2309: /*@C
2310: PetscDSSetConstants - Set the array of constants passed to point functions
2312: Not collective
2314: Input Parameters:
2315: + prob - The PetscDS object
2316: . numConstants - The number of constants
2317: - constants - The array of constants, NULL if there are none
2319: Level: intermediate
2321: .seealso: PetscDSGetConstants(), PetscDSCreate()
2322: @*/
2323: PetscErrorCode PetscDSSetConstants(PetscDS prob, PetscInt numConstants, PetscScalar constants[])
2324: {
2329: if (numConstants != prob->numConstants) {
2330: PetscFree(prob->constants);
2331: prob->numConstants = numConstants;
2332: if (prob->numConstants) {
2333: PetscMalloc1(prob->numConstants, &prob->constants);
2334: } else {
2335: prob->constants = NULL;
2336: }
2337: }
2338: if (prob->numConstants) {
2340: PetscArraycpy(prob->constants, constants, prob->numConstants);
2341: }
2342: return(0);
2343: }
2345: /*@
2346: PetscDSGetFieldIndex - Returns the index of the given field
2348: Not collective
2350: Input Parameters:
2351: + prob - The PetscDS object
2352: - disc - The discretization object
2354: Output Parameter:
2355: . f - The field number
2357: Level: beginner
2359: .seealso: PetscGetDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
2360: @*/
2361: PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f)
2362: {
2363: PetscInt g;
2368: *f = -1;
2369: for (g = 0; g < prob->Nf; ++g) {if (disc == prob->disc[g]) break;}
2370: if (g == prob->Nf) SETERRQ(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS.");
2371: *f = g;
2372: return(0);
2373: }
2375: /*@
2376: PetscDSGetFieldSize - Returns the size of the given field in the full space basis
2378: Not collective
2380: Input Parameters:
2381: + prob - The PetscDS object
2382: - f - The field number
2384: Output Parameter:
2385: . size - The size
2387: Level: beginner
2389: .seealso: PetscDSGetFieldOffset(), PetscDSGetNumFields(), PetscDSCreate()
2390: @*/
2391: PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size)
2392: {
2398: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
2399: PetscDSSetUp(prob);
2400: *size = prob->Nb[f];
2401: return(0);
2402: }
2404: /*@
2405: PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis
2407: Not collective
2409: Input Parameters:
2410: + prob - The PetscDS object
2411: - f - The field number
2413: Output Parameter:
2414: . off - The offset
2416: Level: beginner
2418: .seealso: PetscDSGetFieldSize(), PetscDSGetNumFields(), PetscDSCreate()
2419: @*/
2420: PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off)
2421: {
2422: PetscInt size, g;
2428: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
2429: *off = 0;
2430: for (g = 0; g < f; ++g) {
2431: PetscDSGetFieldSize(prob, g, &size);
2432: *off += size;
2433: }
2434: return(0);
2435: }
2437: /*@
2438: PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point
2440: Not collective
2442: Input Parameter:
2443: . prob - The PetscDS object
2445: Output Parameter:
2446: . dimensions - The number of dimensions
2448: Level: beginner
2450: .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate()
2451: @*/
2452: PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[])
2453: {
2458: PetscDSSetUp(prob);
2460: *dimensions = prob->Nb;
2461: return(0);
2462: }
2464: /*@
2465: PetscDSGetComponents - Returns the number of components for each field on an evaluation point
2467: Not collective
2469: Input Parameter:
2470: . prob - The PetscDS object
2472: Output Parameter:
2473: . components - The number of components
2475: Level: beginner
2477: .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate()
2478: @*/
2479: PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[])
2480: {
2485: PetscDSSetUp(prob);
2487: *components = prob->Nc;
2488: return(0);
2489: }
2491: /*@
2492: PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point
2494: Not collective
2496: Input Parameters:
2497: + prob - The PetscDS object
2498: - f - The field number
2500: Output Parameter:
2501: . off - The offset
2503: Level: beginner
2505: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2506: @*/
2507: PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off)
2508: {
2512: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
2513: *off = prob->off[f];
2514: return(0);
2515: }
2517: /*@
2518: PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point
2520: Not collective
2522: Input Parameter:
2523: . prob - The PetscDS object
2525: Output Parameter:
2526: . offsets - The offsets
2528: Level: beginner
2530: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2531: @*/
2532: PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[])
2533: {
2537: *offsets = prob->off;
2538: return(0);
2539: }
2541: /*@
2542: PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point
2544: Not collective
2546: Input Parameter:
2547: . prob - The PetscDS object
2549: Output Parameter:
2550: . offsets - The offsets
2552: Level: beginner
2554: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2555: @*/
2556: PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[])
2557: {
2561: *offsets = prob->offDer;
2562: return(0);
2563: }
2565: /*@C
2566: PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization
2568: Not collective
2570: Input Parameter:
2571: . prob - The PetscDS object
2573: Output Parameters:
2574: + basis - The basis function tabulation at quadrature points
2575: - basisDer - The basis function derivative tabulation at quadrature points
2577: Level: intermediate
2579: .seealso: PetscDSCreate()
2580: @*/
2581: PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer)
2582: {
2587: PetscDSSetUp(prob);
2590: return(0);
2591: }
2593: /*@C
2594: PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces
2596: Not collective
2598: Input Parameter:
2599: . prob - The PetscDS object
2601: Output Parameters:
2602: + basisFace - The basis function tabulation at quadrature points
2603: - basisDerFace - The basis function derivative tabulation at quadrature points
2605: Level: intermediate
2607: .seealso: PetscDSGetTabulation(), PetscDSCreate()
2608: @*/
2609: PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer)
2610: {
2615: PetscDSSetUp(prob);
2618: return(0);
2619: }
2621: PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x)
2622: {
2627: PetscDSSetUp(prob);
2631: return(0);
2632: }
2634: PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3)
2635: {
2640: PetscDSSetUp(prob);
2647: return(0);
2648: }
2650: PetscErrorCode PetscDSGetWorkspace(PetscDS prob, PetscReal **x, PetscScalar **basisReal, PetscScalar **basisDerReal, PetscScalar **testReal, PetscScalar **testDerReal)
2651: {
2656: PetscDSSetUp(prob);
2662: return(0);
2663: }
2665: /*@C
2666: PetscDSAddBoundary - Add a boundary condition to the model
2668: Input Parameters:
2669: + ds - The PetscDS object
2670: . type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2671: . name - The BC name
2672: . labelname - The label defining constrained points
2673: . field - The field to constrain
2674: . numcomps - The number of constrained field components (0 will constrain all fields)
2675: . comps - An array of constrained component numbers
2676: . bcFunc - A pointwise function giving boundary values
2677: . numids - The number of DMLabel ids for constrained points
2678: . ids - An array of ids for constrained points
2679: - ctx - An optional user context for bcFunc
2681: Options Database Keys:
2682: + -bc_<boundary name> <num> - Overrides the boundary ids
2683: - -bc_<boundary name>_comp <num> - Overrides the boundary components
2685: Level: developer
2687: .seealso: PetscDSGetBoundary()
2688: @*/
2689: PetscErrorCode PetscDSAddBoundary(PetscDS ds, DMBoundaryConditionType type, const char name[], const char labelname[], PetscInt field, PetscInt numcomps, const PetscInt *comps, void (*bcFunc)(void), PetscInt numids, const PetscInt *ids, void *ctx)
2690: {
2691: DSBoundary b;
2696: PetscNew(&b);
2697: PetscStrallocpy(name, (char **) &b->name);
2698: PetscStrallocpy(labelname, (char **) &b->labelname);
2699: PetscMalloc1(numcomps, &b->comps);
2700: if (numcomps) {PetscArraycpy(b->comps, comps, numcomps);}
2701: PetscMalloc1(numids, &b->ids);
2702: if (numids) {PetscArraycpy(b->ids, ids, numids);}
2703: b->type = type;
2704: b->field = field;
2705: b->numcomps = numcomps;
2706: b->func = bcFunc;
2707: b->numids = numids;
2708: b->ctx = ctx;
2709: b->next = ds->boundary;
2710: ds->boundary = b;
2711: return(0);
2712: }
2714: /*@C
2715: PetscDSUpdateBoundary - Change a boundary condition for the model
2717: Input Parameters:
2718: + ds - The PetscDS object
2719: . bd - The boundary condition number
2720: . type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2721: . name - The BC name
2722: . labelname - The label defining constrained points
2723: . field - The field to constrain
2724: . numcomps - The number of constrained field components
2725: . comps - An array of constrained component numbers
2726: . bcFunc - A pointwise function giving boundary values
2727: . numids - The number of DMLabel ids for constrained points
2728: . ids - An array of ids for constrained points
2729: - ctx - An optional user context for bcFunc
2731: Note: The boundary condition number is the order in which it was registered. The user can get the number of boundary conditions from PetscDSGetNumBoundary().
2733: Level: developer
2735: .seealso: PetscDSAddBoundary(), PetscDSGetBoundary(), PetscDSGetNumBoundary()
2736: @*/
2737: PetscErrorCode PetscDSUpdateBoundary(PetscDS ds, PetscInt bd, DMBoundaryConditionType type, const char name[], const char labelname[], PetscInt field, PetscInt numcomps, const PetscInt *comps, void (*bcFunc)(void), PetscInt numids, const PetscInt *ids, void *ctx)
2738: {
2739: DSBoundary b = ds->boundary;
2740: PetscInt n = 0;
2745: while (b) {
2746: if (n == bd) break;
2747: b = b->next;
2748: ++n;
2749: }
2750: if (!b) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %d is not in [0, %d)", bd, n);
2751: if (name) {
2752: PetscFree(b->name);
2753: PetscStrallocpy(name, (char **) &b->name);
2754: }
2755: if (labelname) {
2756: PetscFree(b->labelname);
2757: PetscStrallocpy(labelname, (char **) &b->labelname);
2758: }
2759: if (numcomps >= 0 && numcomps != b->numcomps) {
2760: b->numcomps = numcomps;
2761: PetscFree(b->comps);
2762: PetscMalloc1(numcomps, &b->comps);
2763: if (numcomps) {PetscArraycpy(b->comps, comps, numcomps);}
2764: }
2765: if (numids >= 0 && numids != b->numids) {
2766: b->numids = numids;
2767: PetscFree(b->ids);
2768: PetscMalloc1(numids, &b->ids);
2769: if (numids) {PetscArraycpy(b->ids, ids, numids);}
2770: }
2771: b->type = type;
2772: if (field >= 0) {b->field = field;}
2773: if (bcFunc) {b->func = bcFunc;}
2774: if (ctx) {b->ctx = ctx;}
2775: return(0);
2776: }
2778: /*@
2779: PetscDSGetNumBoundary - Get the number of registered BC
2781: Input Parameters:
2782: . ds - The PetscDS object
2784: Output Parameters:
2785: . numBd - The number of BC
2787: Level: intermediate
2789: .seealso: PetscDSAddBoundary(), PetscDSGetBoundary()
2790: @*/
2791: PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd)
2792: {
2793: DSBoundary b = ds->boundary;
2798: *numBd = 0;
2799: while (b) {++(*numBd); b = b->next;}
2800: return(0);
2801: }
2803: /*@C
2804: PetscDSGetBoundary - Gets a boundary condition to the model
2806: Input Parameters:
2807: + ds - The PetscDS object
2808: - bd - The BC number
2810: Output Parameters:
2811: + type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2812: . name - The BC name
2813: . labelname - The label defining constrained points
2814: . field - The field to constrain
2815: . numcomps - The number of constrained field components
2816: . comps - An array of constrained component numbers
2817: . bcFunc - A pointwise function giving boundary values
2818: . numids - The number of DMLabel ids for constrained points
2819: . ids - An array of ids for constrained points
2820: - ctx - An optional user context for bcFunc
2822: Options Database Keys:
2823: + -bc_<boundary name> <num> - Overrides the boundary ids
2824: - -bc_<boundary name>_comp <num> - Overrides the boundary components
2826: Level: developer
2828: .seealso: PetscDSAddBoundary()
2829: @*/
2830: PetscErrorCode PetscDSGetBoundary(PetscDS ds, PetscInt bd, DMBoundaryConditionType *type, const char **name, const char **labelname, PetscInt *field, PetscInt *numcomps, const PetscInt **comps, void (**func)(void), PetscInt *numids, const PetscInt **ids, void **ctx)
2831: {
2832: DSBoundary b = ds->boundary;
2833: PetscInt n = 0;
2837: while (b) {
2838: if (n == bd) break;
2839: b = b->next;
2840: ++n;
2841: }
2842: if (!b) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %d is not in [0, %d)", bd, n);
2843: if (type) {
2845: *type = b->type;
2846: }
2847: if (name) {
2849: *name = b->name;
2850: }
2851: if (labelname) {
2853: *labelname = b->labelname;
2854: }
2855: if (field) {
2857: *field = b->field;
2858: }
2859: if (numcomps) {
2861: *numcomps = b->numcomps;
2862: }
2863: if (comps) {
2865: *comps = b->comps;
2866: }
2867: if (func) {
2869: *func = b->func;
2870: }
2871: if (numids) {
2873: *numids = b->numids;
2874: }
2875: if (ids) {
2877: *ids = b->ids;
2878: }
2879: if (ctx) {
2881: *ctx = b->ctx;
2882: }
2883: return(0);
2884: }
2886: /*@
2887: PetscDSCopyBoundary - Copy all boundary condition objects to the new problem
2889: Not collective
2891: Input Parameter:
2892: . prob - The PetscDS object
2894: Output Parameter:
2895: . newprob - The PetscDS copy
2897: Level: intermediate
2899: .seealso: PetscDSCopyEquations(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
2900: @*/
2901: PetscErrorCode PetscDSCopyBoundary(PetscDS probA, PetscDS probB)
2902: {
2903: DSBoundary b, next, *lastnext;
2909: if (probA == probB) return(0);
2910: next = probB->boundary;
2911: while (next) {
2912: DSBoundary b = next;
2914: next = b->next;
2915: PetscFree(b->comps);
2916: PetscFree(b->ids);
2917: PetscFree(b->name);
2918: PetscFree(b->labelname);
2919: PetscFree(b);
2920: }
2921: lastnext = &(probB->boundary);
2922: for (b = probA->boundary; b; b = b->next) {
2923: DSBoundary bNew;
2925: PetscNew(&bNew);
2926: bNew->numcomps = b->numcomps;
2927: PetscMalloc1(bNew->numcomps, &bNew->comps);
2928: PetscArraycpy(bNew->comps, b->comps, bNew->numcomps);
2929: bNew->numids = b->numids;
2930: PetscMalloc1(bNew->numids, &bNew->ids);
2931: PetscArraycpy(bNew->ids, b->ids, bNew->numids);
2932: PetscStrallocpy(b->labelname,(char **) &(bNew->labelname));
2933: PetscStrallocpy(b->name,(char **) &(bNew->name));
2934: bNew->ctx = b->ctx;
2935: bNew->type = b->type;
2936: bNew->field = b->field;
2937: bNew->func = b->func;
2939: *lastnext = bNew;
2940: lastnext = &(bNew->next);
2941: }
2942: return(0);
2943: }
2945: /*@C
2946: PetscDSSelectEquations - Copy pointwise function pointers to the new problem with different field layout
2948: Not collective
2950: Input Parameter:
2951: + prob - The PetscDS object
2952: . numFields - Number of new fields
2953: - fields - Old field number for each new field
2955: Output Parameter:
2956: . newprob - The PetscDS copy
2958: Level: intermediate
2960: .seealso: PetscDSCopyBoundary(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
2961: @*/
2962: PetscErrorCode PetscDSSelectEquations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
2963: {
2964: PetscInt Nf, Nfn, fn, gn;
2971: PetscDSGetNumFields(prob, &Nf);
2972: PetscDSGetNumFields(newprob, &Nfn);
2973: if (numFields > Nfn) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Number of fields %D to transfer must not be greater then the total number of fields %D", numFields, Nfn);
2974: for (fn = 0; fn < numFields; ++fn) {
2975: const PetscInt f = fields ? fields[fn] : fn;
2976: PetscPointFunc obj;
2977: PetscPointFunc f0, f1;
2978: PetscBdPointFunc f0Bd, f1Bd;
2979: PetscRiemannFunc r;
2981: if (f >= Nf) continue;
2982: PetscDSGetObjective(prob, f, &obj);
2983: PetscDSGetResidual(prob, f, &f0, &f1);
2984: PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd);
2985: PetscDSGetRiemannSolver(prob, f, &r);
2986: PetscDSSetObjective(newprob, fn, obj);
2987: PetscDSSetResidual(newprob, fn, f0, f1);
2988: PetscDSSetBdResidual(newprob, fn, f0Bd, f1Bd);
2989: PetscDSSetRiemannSolver(newprob, fn, r);
2990: for (gn = 0; gn < numFields; ++gn) {
2991: const PetscInt g = fields ? fields[gn] : gn;
2992: PetscPointJac g0, g1, g2, g3;
2993: PetscPointJac g0p, g1p, g2p, g3p;
2994: PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd;
2996: if (g >= Nf) continue;
2997: PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3);
2998: PetscDSGetJacobianPreconditioner(prob, f, g, &g0p, &g1p, &g2p, &g3p);
2999: PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd);
3000: PetscDSSetJacobian(newprob, fn, gn, g0, g1, g2, g3);
3001: PetscDSSetJacobianPreconditioner(prob, fn, gn, g0p, g1p, g2p, g3p);
3002: PetscDSSetBdJacobian(newprob, fn, gn, g0Bd, g1Bd, g2Bd, g3Bd);
3003: }
3004: }
3005: return(0);
3006: }
3008: /*@
3009: PetscDSCopyEquations - Copy all pointwise function pointers to the new problem
3011: Not collective
3013: Input Parameter:
3014: . prob - The PetscDS object
3016: Output Parameter:
3017: . newprob - The PetscDS copy
3019: Level: intermediate
3021: .seealso: PetscDSCopyBoundary(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
3022: @*/
3023: PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob)
3024: {
3025: PetscInt Nf, Ng;
3031: PetscDSGetNumFields(prob, &Nf);
3032: PetscDSGetNumFields(newprob, &Ng);
3033: if (Nf != Ng) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %D != %D", Nf, Ng);
3034: PetscDSSelectEquations(prob, Nf, NULL, newprob);
3035: return(0);
3036: }
3037: /*@
3038: PetscDSCopyConstants - Copy all constants to the new problem
3040: Not collective
3042: Input Parameter:
3043: . prob - The PetscDS object
3045: Output Parameter:
3046: . newprob - The PetscDS copy
3048: Level: intermediate
3050: .seealso: PetscDSCopyBoundary(), PetscDSCopyEquations(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
3051: @*/
3052: PetscErrorCode PetscDSCopyConstants(PetscDS prob, PetscDS newprob)
3053: {
3054: PetscInt Nc;
3055: const PetscScalar *constants;
3056: PetscErrorCode ierr;
3061: PetscDSGetConstants(prob, &Nc, &constants);
3062: PetscDSSetConstants(newprob, Nc, (PetscScalar *) constants);
3063: return(0);
3064: }
3066: PetscErrorCode PetscDSGetHeightSubspace(PetscDS prob, PetscInt height, PetscDS *subprob)
3067: {
3068: PetscInt dim, Nf, f;
3074: if (height == 0) {*subprob = prob; return(0);}
3075: PetscDSGetNumFields(prob, &Nf);
3076: PetscDSGetSpatialDimension(prob, &dim);
3077: if (height > dim) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_OUTOFRANGE, "DS can only handle height in [0, %D], not %D", dim, height);
3078: if (!prob->subprobs) {PetscCalloc1(dim, &prob->subprobs);}
3079: if (!prob->subprobs[height-1]) {
3080: PetscInt cdim;
3082: PetscDSCreate(PetscObjectComm((PetscObject) prob), &prob->subprobs[height-1]);
3083: PetscDSGetCoordinateDimension(prob, &cdim);
3084: PetscDSSetCoordinateDimension(prob->subprobs[height-1], cdim);
3085: for (f = 0; f < Nf; ++f) {
3086: PetscFE subfe;
3087: PetscObject obj;
3088: PetscClassId id;
3090: PetscDSGetDiscretization(prob, f, &obj);
3091: PetscObjectGetClassId(obj, &id);
3092: if (id == PETSCFE_CLASSID) {PetscFEGetHeightSubspace((PetscFE) obj, height, &subfe);}
3093: else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unsupported discretization type for field %d", f);
3094: PetscDSSetDiscretization(prob->subprobs[height-1], f, (PetscObject) subfe);
3095: }
3096: }
3097: *subprob = prob->subprobs[height-1];
3098: return(0);
3099: }
3101: static PetscErrorCode PetscDSDestroy_Basic(PetscDS prob)
3102: {
3103: PetscErrorCode ierr;
3106: PetscFree(prob->data);
3107: return(0);
3108: }
3110: static PetscErrorCode PetscDSInitialize_Basic(PetscDS prob)
3111: {
3113: prob->ops->setfromoptions = NULL;
3114: prob->ops->setup = NULL;
3115: prob->ops->view = NULL;
3116: prob->ops->destroy = PetscDSDestroy_Basic;
3117: return(0);
3118: }
3120: /*MC
3121: PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions
3123: Level: intermediate
3125: .seealso: PetscDSType, PetscDSCreate(), PetscDSSetType()
3126: M*/
3128: PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS prob)
3129: {
3130: PetscDS_Basic *b;
3131: PetscErrorCode ierr;
3135: PetscNewLog(prob, &b);
3136: prob->data = b;
3138: PetscDSInitialize_Basic(prob);
3139: return(0);
3140: }