Actual source code: dtds.c
petsc-3.13.6 2020-09-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: PetscDSViewFromOptions - View from Options
229: Collective on PetscDS
231: Input Parameters:
232: + A - the PetscDS object
233: . obj - Optional object
234: - name - command line option
236: Level: intermediate
237: .seealso: PetscDS, PetscDSView, PetscObjectViewFromOptions(), PetscDSCreate()
238: @*/
239: PetscErrorCode PetscDSViewFromOptions(PetscDS A,PetscObject obj,const char name[])
240: {
245: PetscObjectViewFromOptions((PetscObject)A,obj,name);
246: return(0);
247: }
249: /*@C
250: PetscDSView - Views a PetscDS
252: Collective on prob
254: Input Parameter:
255: + prob - the PetscDS object to view
256: - v - the viewer
258: Level: developer
260: .seealso PetscDSDestroy()
261: @*/
262: PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v)
263: {
264: PetscBool iascii;
269: if (!v) {PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject) prob), &v);}
271: PetscObjectTypeCompare((PetscObject) v, PETSCVIEWERASCII, &iascii);
272: if (iascii) {PetscDSView_Ascii(prob, v);}
273: if (prob->ops->view) {(*prob->ops->view)(prob, v);}
274: return(0);
275: }
277: /*@
278: PetscDSSetFromOptions - sets parameters in a PetscDS from the options database
280: Collective on prob
282: Input Parameter:
283: . prob - the PetscDS object to set options for
285: Options Database:
286: + -petscds_type <type> : Set the DS type
287: . -petscds_view <view opt> : View the DS
288: . -petscds_jac_pre : Turn formation of a separate Jacobian preconditioner on and off
289: . -bc_<name> <ids> : Specify a list of label ids for a boundary condition
290: - -bc_<name>_comp <comps> : Specify a list of field components to constrain for a boundary condition
292: Level: developer
294: .seealso PetscDSView()
295: @*/
296: PetscErrorCode PetscDSSetFromOptions(PetscDS prob)
297: {
298: DSBoundary b;
299: const char *defaultType;
300: char name[256];
301: PetscBool flg;
306: if (!((PetscObject) prob)->type_name) {
307: defaultType = PETSCDSBASIC;
308: } else {
309: defaultType = ((PetscObject) prob)->type_name;
310: }
311: PetscDSRegisterAll();
313: PetscObjectOptionsBegin((PetscObject) prob);
314: for (b = prob->boundary; b; b = b->next) {
315: char optname[1024];
316: PetscInt ids[1024], len = 1024;
317: PetscBool flg;
319: PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name);
320: PetscMemzero(ids, sizeof(ids));
321: PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg);
322: if (flg) {
323: b->numids = len;
324: PetscFree(b->ids);
325: PetscMalloc1(len, &b->ids);
326: PetscArraycpy(b->ids, ids, len);
327: }
328: len = 1024;
329: PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name);
330: PetscMemzero(ids, sizeof(ids));
331: PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg);
332: if (flg) {
333: b->numcomps = len;
334: PetscFree(b->comps);
335: PetscMalloc1(len, &b->comps);
336: PetscArraycpy(b->comps, ids, len);
337: }
338: }
339: PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg);
340: if (flg) {
341: PetscDSSetType(prob, name);
342: } else if (!((PetscObject) prob)->type_name) {
343: PetscDSSetType(prob, defaultType);
344: }
345: PetscOptionsBool("-petscds_jac_pre", "Discrete System", "PetscDSUseJacobianPreconditioner", prob->useJacPre, &prob->useJacPre, &flg);
346: if (prob->ops->setfromoptions) {(*prob->ops->setfromoptions)(prob);}
347: /* process any options handlers added with PetscObjectAddOptionsHandler() */
348: PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject) prob);
349: PetscOptionsEnd();
350: if (prob->Nf) {PetscDSViewFromOptions(prob, NULL, "-petscds_view");}
351: return(0);
352: }
354: /*@C
355: PetscDSSetUp - Construct data structures for the PetscDS
357: Collective on prob
359: Input Parameter:
360: . prob - the PetscDS object to setup
362: Level: developer
364: .seealso PetscDSView(), PetscDSDestroy()
365: @*/
366: PetscErrorCode PetscDSSetUp(PetscDS prob)
367: {
368: const PetscInt Nf = prob->Nf;
369: PetscInt dim, dimEmbed, NbMax = 0, NcMax = 0, NqMax = 0, NsMax = 1, f;
374: if (prob->setup) return(0);
375: /* Calculate sizes */
376: PetscDSGetSpatialDimension(prob, &dim);
377: PetscDSGetCoordinateDimension(prob, &dimEmbed);
378: prob->totDim = prob->totComp = 0;
379: PetscMalloc2(Nf,&prob->Nc,Nf,&prob->Nb);
380: PetscCalloc2(Nf+1,&prob->off,Nf+1,&prob->offDer);
381: PetscMalloc2(Nf,&prob->T,Nf,&prob->Tf);
382: for (f = 0; f < Nf; ++f) {
383: PetscObject obj;
384: PetscClassId id;
385: PetscQuadrature q;
386: PetscInt Nq = 0, Nb, Nc;
388: PetscDSGetDiscretization(prob, f, &obj);
389: PetscObjectGetClassId(obj, &id);
390: if (id == PETSCFE_CLASSID) {
391: PetscFE fe = (PetscFE) obj;
393: PetscFEGetQuadrature(fe, &q);
394: PetscFEGetDimension(fe, &Nb);
395: PetscFEGetNumComponents(fe, &Nc);
396: PetscFEGetCellTabulation(fe, &prob->T[f]);
397: PetscFEGetFaceTabulation(fe, &prob->Tf[f]);
398: } else if (id == PETSCFV_CLASSID) {
399: PetscFV fv = (PetscFV) obj;
401: PetscFVGetQuadrature(fv, &q);
402: PetscFVGetNumComponents(fv, &Nc);
403: Nb = Nc;
404: PetscFVGetCellTabulation(fv, &prob->T[f]);
405: /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */
406: } else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", f);
407: prob->Nc[f] = Nc;
408: prob->Nb[f] = Nb;
409: prob->off[f+1] = Nc + prob->off[f];
410: prob->offDer[f+1] = Nc*dim + prob->offDer[f];
411: if (q) {PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL);}
412: NqMax = PetscMax(NqMax, Nq);
413: NbMax = PetscMax(NbMax, Nb);
414: NcMax = PetscMax(NcMax, Nc);
415: prob->totDim += Nb;
416: prob->totComp += Nc;
417: /* There are two faces for all fields but the cohesive field on a hybrid cell */
418: if (prob->isHybrid && (f < Nf-1)) prob->totDim += Nb;
419: }
420: /* Allocate works space */
421: if (prob->isHybrid) NsMax = 2;
422: PetscMalloc3(NsMax*prob->totComp,&prob->u,NsMax*prob->totComp,&prob->u_t,NsMax*prob->totComp*dimEmbed,&prob->u_x);
423: PetscMalloc5(dimEmbed,&prob->x,NbMax*NcMax,&prob->basisReal,NbMax*NcMax*dimEmbed,&prob->basisDerReal,NbMax*NcMax,&prob->testReal,NbMax*NcMax*dimEmbed,&prob->testDerReal);
424: PetscMalloc6(NsMax*NqMax*NcMax,&prob->f0,NsMax*NqMax*NcMax*dim,&prob->f1,
425: NsMax*NsMax*NqMax*NcMax*NcMax,&prob->g0,NsMax*NsMax*NqMax*NcMax*NcMax*dim,&prob->g1,
426: NsMax*NsMax*NqMax*NcMax*NcMax*dim,&prob->g2,NsMax*NsMax*NqMax*NcMax*NcMax*dim*dim,&prob->g3);
427: if (prob->ops->setup) {(*prob->ops->setup)(prob);}
428: prob->setup = PETSC_TRUE;
429: return(0);
430: }
432: static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob)
433: {
437: PetscFree2(prob->Nc,prob->Nb);
438: PetscFree2(prob->off,prob->offDer);
439: PetscFree2(prob->T,prob->Tf);
440: PetscFree3(prob->u,prob->u_t,prob->u_x);
441: PetscFree5(prob->x,prob->basisReal, prob->basisDerReal,prob->testReal,prob->testDerReal);
442: PetscFree6(prob->f0,prob->f1,prob->g0,prob->g1,prob->g2,prob->g3);
443: return(0);
444: }
446: static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew)
447: {
448: PetscObject *tmpd;
449: PetscBool *tmpi;
450: PetscPointFunc *tmpobj, *tmpf, *tmpup;
451: PetscPointJac *tmpg, *tmpgp, *tmpgt;
452: PetscBdPointFunc *tmpfbd;
453: PetscBdPointJac *tmpgbd;
454: PetscRiemannFunc *tmpr;
455: PetscSimplePointFunc *tmpexactSol;
456: void **tmpexactCtx;
457: void **tmpctx;
458: PetscInt Nf = prob->Nf, f;
459: PetscErrorCode ierr;
462: if (Nf >= NfNew) return(0);
463: prob->setup = PETSC_FALSE;
464: PetscDSDestroyStructs_Static(prob);
465: PetscMalloc2(NfNew, &tmpd, NfNew, &tmpi);
466: for (f = 0; f < Nf; ++f) {tmpd[f] = prob->disc[f]; tmpi[f] = prob->implicit[f];}
467: for (f = Nf; f < NfNew; ++f) {tmpd[f] = NULL; tmpi[f] = PETSC_TRUE;}
468: PetscFree2(prob->disc, prob->implicit);
469: prob->Nf = NfNew;
470: prob->disc = tmpd;
471: prob->implicit = tmpi;
472: PetscCalloc7(NfNew, &tmpobj, NfNew*2, &tmpf, NfNew*NfNew*4, &tmpg, NfNew*NfNew*4, &tmpgp, NfNew*NfNew*4, &tmpgt, NfNew, &tmpr, NfNew, &tmpctx);
473: PetscCalloc1(NfNew, &tmpup);
474: for (f = 0; f < Nf; ++f) tmpobj[f] = prob->obj[f];
475: for (f = 0; f < Nf*2; ++f) tmpf[f] = prob->f[f];
476: for (f = 0; f < Nf*Nf*4; ++f) tmpg[f] = prob->g[f];
477: for (f = 0; f < Nf*Nf*4; ++f) tmpgp[f] = prob->gp[f];
478: for (f = 0; f < Nf; ++f) tmpr[f] = prob->r[f];
479: for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f];
480: for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f];
481: for (f = Nf; f < NfNew; ++f) tmpobj[f] = NULL;
482: for (f = Nf*2; f < NfNew*2; ++f) tmpf[f] = NULL;
483: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpg[f] = NULL;
484: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgp[f] = NULL;
485: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgt[f] = NULL;
486: for (f = Nf; f < NfNew; ++f) tmpr[f] = NULL;
487: for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL;
488: for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL;
489: PetscFree7(prob->obj, prob->f, prob->g, prob->gp, prob->gt, prob->r, prob->ctx);
490: PetscFree(prob->update);
491: prob->obj = tmpobj;
492: prob->f = tmpf;
493: prob->g = tmpg;
494: prob->gp = tmpgp;
495: prob->gt = tmpgt;
496: prob->r = tmpr;
497: prob->update = tmpup;
498: prob->ctx = tmpctx;
499: PetscCalloc4(NfNew*2, &tmpfbd, NfNew*NfNew*4, &tmpgbd, NfNew, &tmpexactSol, NfNew, &tmpexactCtx);
500: for (f = 0; f < Nf*2; ++f) tmpfbd[f] = prob->fBd[f];
501: for (f = 0; f < Nf*Nf*4; ++f) tmpgbd[f] = prob->gBd[f];
502: for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f];
503: for (f = 0; f < Nf; ++f) tmpexactCtx[f] = prob->exactCtx[f];
504: for (f = Nf*2; f < NfNew*2; ++f) tmpfbd[f] = NULL;
505: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgbd[f] = NULL;
506: for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL;
507: for (f = Nf; f < NfNew; ++f) tmpexactCtx[f] = NULL;
508: PetscFree4(prob->fBd, prob->gBd, prob->exactSol, prob->exactCtx);
509: prob->fBd = tmpfbd;
510: prob->gBd = tmpgbd;
511: prob->exactSol = tmpexactSol;
512: prob->exactCtx = tmpexactCtx;
513: return(0);
514: }
516: /*@
517: PetscDSDestroy - Destroys a PetscDS object
519: Collective on prob
521: Input Parameter:
522: . prob - the PetscDS object to destroy
524: Level: developer
526: .seealso PetscDSView()
527: @*/
528: PetscErrorCode PetscDSDestroy(PetscDS *prob)
529: {
530: PetscInt f;
531: DSBoundary next;
535: if (!*prob) return(0);
538: if (--((PetscObject)(*prob))->refct > 0) {*prob = 0; return(0);}
539: ((PetscObject) (*prob))->refct = 0;
540: if ((*prob)->subprobs) {
541: PetscInt dim, d;
543: PetscDSGetSpatialDimension(*prob, &dim);
544: for (d = 0; d < dim; ++d) {PetscDSDestroy(&(*prob)->subprobs[d]);}
545: }
546: PetscFree((*prob)->subprobs);
547: PetscDSDestroyStructs_Static(*prob);
548: for (f = 0; f < (*prob)->Nf; ++f) {
549: PetscObjectDereference((*prob)->disc[f]);
550: }
551: PetscFree2((*prob)->disc, (*prob)->implicit);
552: PetscFree7((*prob)->obj,(*prob)->f,(*prob)->g,(*prob)->gp,(*prob)->gt,(*prob)->r,(*prob)->ctx);
553: PetscFree((*prob)->update);
554: PetscFree4((*prob)->fBd,(*prob)->gBd,(*prob)->exactSol,(*prob)->exactCtx);
555: if ((*prob)->ops->destroy) {(*(*prob)->ops->destroy)(*prob);}
556: next = (*prob)->boundary;
557: while (next) {
558: DSBoundary b = next;
560: next = b->next;
561: PetscFree(b->comps);
562: PetscFree(b->ids);
563: PetscFree(b->name);
564: PetscFree(b->labelname);
565: PetscFree(b);
566: }
567: PetscFree((*prob)->constants);
568: PetscHeaderDestroy(prob);
569: return(0);
570: }
572: /*@
573: PetscDSCreate - Creates an empty PetscDS object. The type can then be set with PetscDSSetType().
575: Collective
577: Input Parameter:
578: . comm - The communicator for the PetscDS object
580: Output Parameter:
581: . prob - The PetscDS object
583: Level: beginner
585: .seealso: PetscDSSetType(), PETSCDSBASIC
586: @*/
587: PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *prob)
588: {
589: PetscDS p;
594: *prob = NULL;
595: PetscDSInitializePackage();
597: PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView);
599: p->Nf = 0;
600: p->setup = PETSC_FALSE;
601: p->numConstants = 0;
602: p->constants = NULL;
603: p->dimEmbed = -1;
604: p->useJacPre = PETSC_TRUE;
606: *prob = p;
607: return(0);
608: }
610: /*@
611: PetscDSGetNumFields - Returns the number of fields in the DS
613: Not collective
615: Input Parameter:
616: . prob - The PetscDS object
618: Output Parameter:
619: . Nf - The number of fields
621: Level: beginner
623: .seealso: PetscDSGetSpatialDimension(), PetscDSCreate()
624: @*/
625: PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf)
626: {
630: *Nf = prob->Nf;
631: return(0);
632: }
634: /*@
635: PetscDSGetSpatialDimension - Returns the spatial dimension of the DS, meaning the topological dimension of the discretizations
637: Not collective
639: Input Parameter:
640: . prob - The PetscDS object
642: Output Parameter:
643: . dim - The spatial dimension
645: Level: beginner
647: .seealso: PetscDSGetCoordinateDimension(), PetscDSGetNumFields(), PetscDSCreate()
648: @*/
649: PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim)
650: {
656: *dim = 0;
657: if (prob->Nf) {
658: PetscObject obj;
659: PetscClassId id;
661: PetscDSGetDiscretization(prob, 0, &obj);
662: PetscObjectGetClassId(obj, &id);
663: if (id == PETSCFE_CLASSID) {PetscFEGetSpatialDimension((PetscFE) obj, dim);}
664: else if (id == PETSCFV_CLASSID) {PetscFVGetSpatialDimension((PetscFV) obj, dim);}
665: else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
666: }
667: return(0);
668: }
670: /*@
671: PetscDSGetCoordinateDimension - Returns the coordinate dimension of the DS, meaning the dimension of the space into which the discretiaztions are embedded
673: Not collective
675: Input Parameter:
676: . prob - The PetscDS object
678: Output Parameter:
679: . dimEmbed - The coordinate dimension
681: Level: beginner
683: .seealso: PetscDSSetCoordinateDimension(), PetscDSGetSpatialDimension(), PetscDSGetNumFields(), PetscDSCreate()
684: @*/
685: PetscErrorCode PetscDSGetCoordinateDimension(PetscDS prob, PetscInt *dimEmbed)
686: {
690: if (prob->dimEmbed < 0) SETERRQ(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONGSTATE, "No coordinate dimension set for this DS");
691: *dimEmbed = prob->dimEmbed;
692: return(0);
693: }
695: /*@
696: PetscDSSetCoordinateDimension - Set the coordinate dimension of the DS, meaning the dimension of the space into which the discretiaztions are embedded
698: Logically collective on prob
700: Input Parameters:
701: + prob - The PetscDS object
702: - dimEmbed - The coordinate dimension
704: Level: beginner
706: .seealso: PetscDSGetCoordinateDimension(), PetscDSGetSpatialDimension(), PetscDSGetNumFields(), PetscDSCreate()
707: @*/
708: PetscErrorCode PetscDSSetCoordinateDimension(PetscDS prob, PetscInt dimEmbed)
709: {
712: if (dimEmbed < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Coordinate dimension must be non-negative, not %D", dimEmbed);
713: prob->dimEmbed = dimEmbed;
714: return(0);
715: }
717: /*@
718: PetscDSGetHybrid - Returns the flag for a hybrid (cohesive) cell
720: Not collective
722: Input Parameter:
723: . prob - The PetscDS object
725: Output Parameter:
726: . isHybrid - The flag
728: Level: developer
730: .seealso: PetscDSSetHybrid(), PetscDSCreate()
731: @*/
732: PetscErrorCode PetscDSGetHybrid(PetscDS prob, PetscBool *isHybrid)
733: {
737: *isHybrid = prob->isHybrid;
738: return(0);
739: }
741: /*@
742: PetscDSSetHybrid - Set the flag for a hybrid (cohesive) cell
744: Not collective
746: Input Parameters:
747: + prob - The PetscDS object
748: - isHybrid - The flag
750: Level: developer
752: .seealso: PetscDSGetHybrid(), PetscDSCreate()
753: @*/
754: PetscErrorCode PetscDSSetHybrid(PetscDS prob, PetscBool isHybrid)
755: {
758: prob->isHybrid = isHybrid;
759: return(0);
760: }
762: /*@
763: PetscDSGetTotalDimension - Returns the total size of the approximation space for this system
765: Not collective
767: Input Parameter:
768: . prob - The PetscDS object
770: Output Parameter:
771: . dim - The total problem dimension
773: Level: beginner
775: .seealso: PetscDSGetNumFields(), PetscDSCreate()
776: @*/
777: PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim)
778: {
783: PetscDSSetUp(prob);
785: *dim = prob->totDim;
786: return(0);
787: }
789: /*@
790: PetscDSGetTotalComponents - Returns the total number of components in this system
792: Not collective
794: Input Parameter:
795: . prob - The PetscDS object
797: Output Parameter:
798: . dim - The total number of components
800: Level: beginner
802: .seealso: PetscDSGetNumFields(), PetscDSCreate()
803: @*/
804: PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc)
805: {
810: PetscDSSetUp(prob);
812: *Nc = prob->totComp;
813: return(0);
814: }
816: /*@
817: PetscDSGetDiscretization - Returns the discretization object for the given field
819: Not collective
821: Input Parameters:
822: + prob - The PetscDS object
823: - f - The field number
825: Output Parameter:
826: . disc - The discretization object
828: Level: beginner
830: .seealso: PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
831: @*/
832: PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc)
833: {
837: 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);
838: *disc = prob->disc[f];
839: return(0);
840: }
842: /*@
843: PetscDSSetDiscretization - Sets the discretization object for the given field
845: Not collective
847: Input Parameters:
848: + prob - The PetscDS object
849: . f - The field number
850: - disc - The discretization object
852: Level: beginner
854: .seealso: PetscDSGetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
855: @*/
856: PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc)
857: {
863: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
864: PetscDSEnlarge_Static(prob, f+1);
865: PetscObjectDereference(prob->disc[f]);
866: prob->disc[f] = disc;
867: PetscObjectReference(disc);
868: {
869: PetscClassId id;
871: PetscObjectGetClassId(disc, &id);
872: if (id == PETSCFE_CLASSID) {
873: PetscDSSetImplicit(prob, f, PETSC_TRUE);
874: } else if (id == PETSCFV_CLASSID) {
875: PetscDSSetImplicit(prob, f, PETSC_FALSE);
876: }
877: }
878: return(0);
879: }
881: /*@
882: PetscDSAddDiscretization - Adds a discretization object
884: Not collective
886: Input Parameters:
887: + prob - The PetscDS object
888: - disc - The boundary discretization object
890: Level: beginner
892: .seealso: PetscDSGetDiscretization(), PetscDSSetDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
893: @*/
894: PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc)
895: {
899: PetscDSSetDiscretization(prob, prob->Nf, disc);
900: return(0);
901: }
903: /*@
904: PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for IMEX
906: Not collective
908: Input Parameters:
909: + prob - The PetscDS object
910: - f - The field number
912: Output Parameter:
913: . implicit - The flag indicating what kind of solve to use for this field
915: Level: developer
917: .seealso: PetscDSSetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
918: @*/
919: PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit)
920: {
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: *implicit = prob->implicit[f];
926: return(0);
927: }
929: /*@
930: PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for IMEX
932: Not collective
934: Input Parameters:
935: + prob - The PetscDS object
936: . f - The field number
937: - implicit - The flag indicating what kind of solve to use for this field
939: Level: developer
941: .seealso: PetscDSGetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
942: @*/
943: PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit)
944: {
947: 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);
948: prob->implicit[f] = implicit;
949: return(0);
950: }
952: PetscErrorCode PetscDSGetObjective(PetscDS prob, PetscInt f,
953: void (**obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
954: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
955: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
956: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
957: {
961: 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);
962: *obj = prob->obj[f];
963: return(0);
964: }
966: PetscErrorCode PetscDSSetObjective(PetscDS prob, PetscInt f,
967: void (*obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
968: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
969: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
970: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
971: {
977: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
978: PetscDSEnlarge_Static(prob, f+1);
979: prob->obj[f] = obj;
980: return(0);
981: }
983: /*@C
984: PetscDSGetResidual - Get the pointwise residual function for a given test field
986: Not collective
988: Input Parameters:
989: + prob - The PetscDS
990: - f - The test field number
992: Output Parameters:
993: + f0 - integrand for the test function term
994: - f1 - integrand for the test function gradient term
996: Note: We are using a first order FEM model for the weak form:
998: \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)
1000: The calling sequence for the callbacks f0 and f1 is given by:
1002: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1003: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1004: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1005: $ PetscReal t, const PetscReal x[], PetscScalar f0[])
1007: + dim - the spatial dimension
1008: . Nf - the number of fields
1009: . uOff - the offset into u[] and u_t[] for each field
1010: . uOff_x - the offset into u_x[] for each field
1011: . u - each field evaluated at the current point
1012: . u_t - the time derivative of each field evaluated at the current point
1013: . u_x - the gradient of each field evaluated at the current point
1014: . aOff - the offset into a[] and a_t[] for each auxiliary field
1015: . aOff_x - the offset into a_x[] for each auxiliary field
1016: . a - each auxiliary field evaluated at the current point
1017: . a_t - the time derivative of each auxiliary field evaluated at the current point
1018: . a_x - the gradient of auxiliary each field evaluated at the current point
1019: . t - current time
1020: . x - coordinates of the current point
1021: . numConstants - number of constant parameters
1022: . constants - constant parameters
1023: - f0 - output values at the current point
1025: Level: intermediate
1027: .seealso: PetscDSSetResidual()
1028: @*/
1029: PetscErrorCode PetscDSGetResidual(PetscDS prob, PetscInt f,
1030: void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1031: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1032: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1033: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1034: void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1035: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1036: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1037: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1038: {
1041: 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);
1044: return(0);
1045: }
1047: /*@C
1048: PetscDSSetResidual - Set the pointwise residual function for a given test field
1050: Not collective
1052: Input Parameters:
1053: + prob - The PetscDS
1054: . f - The test field number
1055: . f0 - integrand for the test function term
1056: - f1 - integrand for the test function gradient term
1058: Note: We are using a first order FEM model for the weak form:
1060: \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)
1062: The calling sequence for the callbacks f0 and f1 is given by:
1064: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1065: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1066: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1067: $ PetscReal t, const PetscReal x[], PetscScalar f0[])
1069: + dim - the spatial dimension
1070: . Nf - the number of fields
1071: . uOff - the offset into u[] and u_t[] for each field
1072: . uOff_x - the offset into u_x[] for each field
1073: . u - each field evaluated at the current point
1074: . u_t - the time derivative of each field evaluated at the current point
1075: . u_x - the gradient of each field evaluated at the current point
1076: . aOff - the offset into a[] and a_t[] for each auxiliary field
1077: . aOff_x - the offset into a_x[] for each auxiliary field
1078: . a - each auxiliary field evaluated at the current point
1079: . a_t - the time derivative of each auxiliary field evaluated at the current point
1080: . a_x - the gradient of auxiliary each field evaluated at the current point
1081: . t - current time
1082: . x - coordinates of the current point
1083: . numConstants - number of constant parameters
1084: . constants - constant parameters
1085: - f0 - output values at the current point
1087: Level: intermediate
1089: .seealso: PetscDSGetResidual()
1090: @*/
1091: PetscErrorCode PetscDSSetResidual(PetscDS prob, PetscInt f,
1092: void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1093: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1094: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1095: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1096: void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1097: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1098: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1099: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1100: {
1107: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1108: PetscDSEnlarge_Static(prob, f+1);
1109: prob->f[f*2+0] = f0;
1110: prob->f[f*2+1] = f1;
1111: return(0);
1112: }
1114: /*@C
1115: PetscDSHasJacobian - Signals that Jacobian functions have been set
1117: Not collective
1119: Input Parameter:
1120: . prob - The PetscDS
1122: Output Parameter:
1123: . hasJac - flag that pointwise function for the Jacobian has been set
1125: Level: intermediate
1127: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1128: @*/
1129: PetscErrorCode PetscDSHasJacobian(PetscDS prob, PetscBool *hasJac)
1130: {
1131: PetscInt f, g, h;
1135: *hasJac = PETSC_FALSE;
1136: for (f = 0; f < prob->Nf; ++f) {
1137: for (g = 0; g < prob->Nf; ++g) {
1138: for (h = 0; h < 4; ++h) {
1139: if (prob->g[(f*prob->Nf + g)*4+h]) *hasJac = PETSC_TRUE;
1140: }
1141: }
1142: }
1143: return(0);
1144: }
1146: /*@C
1147: PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field
1149: Not collective
1151: Input Parameters:
1152: + prob - The PetscDS
1153: . f - The test field number
1154: - g - The field number
1156: Output Parameters:
1157: + g0 - integrand for the test and basis function term
1158: . g1 - integrand for the test function and basis function gradient term
1159: . g2 - integrand for the test function gradient and basis function term
1160: - g3 - integrand for the test function gradient and basis function gradient term
1162: Note: We are using a first order FEM model for the weak form:
1164: \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
1166: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1168: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1169: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1170: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1171: $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1173: + dim - the spatial dimension
1174: . Nf - the number of fields
1175: . uOff - the offset into u[] and u_t[] for each field
1176: . uOff_x - the offset into u_x[] for each field
1177: . u - each field evaluated at the current point
1178: . u_t - the time derivative of each field evaluated at the current point
1179: . u_x - the gradient of each field evaluated at the current point
1180: . aOff - the offset into a[] and a_t[] for each auxiliary field
1181: . aOff_x - the offset into a_x[] for each auxiliary field
1182: . a - each auxiliary field evaluated at the current point
1183: . a_t - the time derivative of each auxiliary field evaluated at the current point
1184: . a_x - the gradient of auxiliary each field evaluated at the current point
1185: . t - current time
1186: . u_tShift - the multiplier a for dF/dU_t
1187: . x - coordinates of the current point
1188: . numConstants - number of constant parameters
1189: . constants - constant parameters
1190: - g0 - output values at the current point
1192: Level: intermediate
1194: .seealso: PetscDSSetJacobian()
1195: @*/
1196: PetscErrorCode PetscDSGetJacobian(PetscDS prob, PetscInt f, PetscInt g,
1197: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1198: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1199: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1200: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1201: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1202: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1203: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1204: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1205: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1206: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1207: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1208: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1209: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1210: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1211: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1212: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1213: {
1216: 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);
1217: 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);
1222: return(0);
1223: }
1225: /*@C
1226: PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields
1228: Not collective
1230: Input Parameters:
1231: + prob - The PetscDS
1232: . f - The test field number
1233: . g - The field number
1234: . g0 - integrand for the test and basis function term
1235: . g1 - integrand for the test function and basis function gradient term
1236: . g2 - integrand for the test function gradient and basis function term
1237: - g3 - integrand for the test function gradient and basis function gradient term
1239: Note: We are using a first order FEM model for the weak form:
1241: \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
1243: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1245: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1246: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1247: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1248: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1250: + dim - the spatial dimension
1251: . Nf - the number of fields
1252: . uOff - the offset into u[] and u_t[] for each field
1253: . uOff_x - the offset into u_x[] for each field
1254: . u - each field evaluated at the current point
1255: . u_t - the time derivative of each field evaluated at the current point
1256: . u_x - the gradient of each field evaluated at the current point
1257: . aOff - the offset into a[] and a_t[] for each auxiliary field
1258: . aOff_x - the offset into a_x[] for each auxiliary field
1259: . a - each auxiliary field evaluated at the current point
1260: . a_t - the time derivative of each auxiliary field evaluated at the current point
1261: . a_x - the gradient of auxiliary each field evaluated at the current point
1262: . t - current time
1263: . u_tShift - the multiplier a for dF/dU_t
1264: . x - coordinates of the current point
1265: . numConstants - number of constant parameters
1266: . constants - constant parameters
1267: - g0 - output values at the current point
1269: Level: intermediate
1271: .seealso: PetscDSGetJacobian()
1272: @*/
1273: PetscErrorCode PetscDSSetJacobian(PetscDS prob, PetscInt f, PetscInt g,
1274: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1275: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1276: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1277: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1278: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1279: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1280: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1281: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1282: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1283: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1284: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1285: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1286: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1287: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1288: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1289: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1290: {
1299: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1300: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1301: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1302: prob->g[(f*prob->Nf + g)*4+0] = g0;
1303: prob->g[(f*prob->Nf + g)*4+1] = g1;
1304: prob->g[(f*prob->Nf + g)*4+2] = g2;
1305: prob->g[(f*prob->Nf + g)*4+3] = g3;
1306: return(0);
1307: }
1309: /*@C
1310: PetscDSUseJacobianPreconditioner - Whether to construct a Jacobian preconditioner
1312: Not collective
1314: Input Parameters:
1315: + prob - The PetscDS
1316: - useJacPre - flag that enables construction of a Jacobian preconditioner
1318: Level: intermediate
1320: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1321: @*/
1322: PetscErrorCode PetscDSUseJacobianPreconditioner(PetscDS prob, PetscBool useJacPre)
1323: {
1326: prob->useJacPre = useJacPre;
1327: return(0);
1328: }
1330: /*@C
1331: PetscDSHasJacobianPreconditioner - Signals that a Jacobian preconditioner matrix has been set
1333: Not collective
1335: Input Parameter:
1336: . prob - The PetscDS
1338: Output Parameter:
1339: . hasJacPre - flag that pointwise function for Jacobian preconditioner matrix has been set
1341: Level: intermediate
1343: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1344: @*/
1345: PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS prob, PetscBool *hasJacPre)
1346: {
1347: PetscInt f, g, h;
1351: *hasJacPre = PETSC_FALSE;
1352: if (!prob->useJacPre) return(0);
1353: for (f = 0; f < prob->Nf; ++f) {
1354: for (g = 0; g < prob->Nf; ++g) {
1355: for (h = 0; h < 4; ++h) {
1356: if (prob->gp[(f*prob->Nf + g)*4+h]) *hasJacPre = PETSC_TRUE;
1357: }
1358: }
1359: }
1360: return(0);
1361: }
1363: /*@C
1364: 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.
1366: Not collective
1368: Input Parameters:
1369: + prob - The PetscDS
1370: . f - The test field number
1371: - g - The field number
1373: Output Parameters:
1374: + g0 - integrand for the test and basis function term
1375: . g1 - integrand for the test function and basis function gradient term
1376: . g2 - integrand for the test function gradient and basis function term
1377: - g3 - integrand for the test function gradient and basis function gradient term
1379: Note: We are using a first order FEM model for the weak form:
1381: \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
1383: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1385: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1386: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1387: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1388: $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1390: + dim - the spatial dimension
1391: . Nf - the number of fields
1392: . uOff - the offset into u[] and u_t[] for each field
1393: . uOff_x - the offset into u_x[] for each field
1394: . u - each field evaluated at the current point
1395: . u_t - the time derivative of each field evaluated at the current point
1396: . u_x - the gradient of each field evaluated at the current point
1397: . aOff - the offset into a[] and a_t[] for each auxiliary field
1398: . aOff_x - the offset into a_x[] for each auxiliary field
1399: . a - each auxiliary field evaluated at the current point
1400: . a_t - the time derivative of each auxiliary field evaluated at the current point
1401: . a_x - the gradient of auxiliary each field evaluated at the current point
1402: . t - current time
1403: . u_tShift - the multiplier a for dF/dU_t
1404: . x - coordinates of the current point
1405: . numConstants - number of constant parameters
1406: . constants - constant parameters
1407: - g0 - output values at the current point
1409: Level: intermediate
1411: .seealso: PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1412: @*/
1413: PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g,
1414: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1415: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1416: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1417: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1418: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1419: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1420: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1421: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1422: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1423: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1424: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1425: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1426: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1427: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1428: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1429: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1430: {
1433: 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);
1434: 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);
1439: return(0);
1440: }
1442: /*@C
1443: 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.
1445: Not collective
1447: Input Parameters:
1448: + prob - The PetscDS
1449: . f - The test field number
1450: . g - The field number
1451: . g0 - integrand for the test and basis function term
1452: . g1 - integrand for the test function and basis function gradient term
1453: . g2 - integrand for the test function gradient and basis function term
1454: - g3 - integrand for the test function gradient and basis function gradient term
1456: Note: We are using a first order FEM model for the weak form:
1458: \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
1460: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1462: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1463: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1464: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1465: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1467: + dim - the spatial dimension
1468: . Nf - the number of fields
1469: . uOff - the offset into u[] and u_t[] for each field
1470: . uOff_x - the offset into u_x[] for each field
1471: . u - each field evaluated at the current point
1472: . u_t - the time derivative of each field evaluated at the current point
1473: . u_x - the gradient of each field evaluated at the current point
1474: . aOff - the offset into a[] and a_t[] for each auxiliary field
1475: . aOff_x - the offset into a_x[] for each auxiliary field
1476: . a - each auxiliary field evaluated at the current point
1477: . a_t - the time derivative of each auxiliary field evaluated at the current point
1478: . a_x - the gradient of auxiliary each field evaluated at the current point
1479: . t - current time
1480: . u_tShift - the multiplier a for dF/dU_t
1481: . x - coordinates of the current point
1482: . numConstants - number of constant parameters
1483: . constants - constant parameters
1484: - g0 - output values at the current point
1486: Level: intermediate
1488: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobian()
1489: @*/
1490: PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g,
1491: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1492: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1493: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1494: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1495: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1496: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1497: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1498: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1499: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1500: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1501: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1502: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1503: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1504: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1505: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1506: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1507: {
1516: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1517: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1518: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1519: prob->gp[(f*prob->Nf + g)*4+0] = g0;
1520: prob->gp[(f*prob->Nf + g)*4+1] = g1;
1521: prob->gp[(f*prob->Nf + g)*4+2] = g2;
1522: prob->gp[(f*prob->Nf + g)*4+3] = g3;
1523: return(0);
1524: }
1526: /*@C
1527: PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, dF/du_t, has been set
1529: Not collective
1531: Input Parameter:
1532: . prob - The PetscDS
1534: Output Parameter:
1535: . hasDynJac - flag that pointwise function for dynamic Jacobian has been set
1537: Level: intermediate
1539: .seealso: PetscDSGetDynamicJacobian(), PetscDSSetDynamicJacobian(), PetscDSGetJacobian()
1540: @*/
1541: PetscErrorCode PetscDSHasDynamicJacobian(PetscDS prob, PetscBool *hasDynJac)
1542: {
1543: PetscInt f, g, h;
1547: *hasDynJac = PETSC_FALSE;
1548: for (f = 0; f < prob->Nf; ++f) {
1549: for (g = 0; g < prob->Nf; ++g) {
1550: for (h = 0; h < 4; ++h) {
1551: if (prob->gt[(f*prob->Nf + g)*4+h]) *hasDynJac = PETSC_TRUE;
1552: }
1553: }
1554: }
1555: return(0);
1556: }
1558: /*@C
1559: PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, dF/du_t, function for given test and basis field
1561: Not collective
1563: Input Parameters:
1564: + prob - The PetscDS
1565: . f - The test field number
1566: - g - The field number
1568: Output Parameters:
1569: + g0 - integrand for the test and basis function term
1570: . g1 - integrand for the test function and basis function gradient term
1571: . g2 - integrand for the test function gradient and basis function term
1572: - g3 - integrand for the test function gradient and basis function gradient term
1574: Note: We are using a first order FEM model for the weak form:
1576: \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
1578: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1580: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1581: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1582: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1583: $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1585: + dim - the spatial dimension
1586: . Nf - the number of fields
1587: . uOff - the offset into u[] and u_t[] for each field
1588: . uOff_x - the offset into u_x[] for each field
1589: . u - each field evaluated at the current point
1590: . u_t - the time derivative of each field evaluated at the current point
1591: . u_x - the gradient of each field evaluated at the current point
1592: . aOff - the offset into a[] and a_t[] for each auxiliary field
1593: . aOff_x - the offset into a_x[] for each auxiliary field
1594: . a - each auxiliary field evaluated at the current point
1595: . a_t - the time derivative of each auxiliary field evaluated at the current point
1596: . a_x - the gradient of auxiliary each field evaluated at the current point
1597: . t - current time
1598: . u_tShift - the multiplier a for dF/dU_t
1599: . x - coordinates of the current point
1600: . numConstants - number of constant parameters
1601: . constants - constant parameters
1602: - g0 - output values at the current point
1604: Level: intermediate
1606: .seealso: PetscDSSetJacobian()
1607: @*/
1608: PetscErrorCode PetscDSGetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g,
1609: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1610: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1611: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1612: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1613: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1614: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1615: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1616: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1617: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1618: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1619: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1620: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1621: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1622: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1623: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1624: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1625: {
1628: 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);
1629: 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);
1634: return(0);
1635: }
1637: /*@C
1638: PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, dF/du_t, function for given test and basis fields
1640: Not collective
1642: Input Parameters:
1643: + prob - The PetscDS
1644: . f - The test field number
1645: . g - The field number
1646: . g0 - integrand for the test and basis function term
1647: . g1 - integrand for the test function and basis function gradient term
1648: . g2 - integrand for the test function gradient and basis function term
1649: - g3 - integrand for the test function gradient and basis function gradient term
1651: Note: We are using a first order FEM model for the weak form:
1653: \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi
1655: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1657: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1658: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1659: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1660: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1662: + dim - the spatial dimension
1663: . Nf - the number of fields
1664: . uOff - the offset into u[] and u_t[] for each field
1665: . uOff_x - the offset into u_x[] for each field
1666: . u - each field evaluated at the current point
1667: . u_t - the time derivative of each field evaluated at the current point
1668: . u_x - the gradient of each field evaluated at the current point
1669: . aOff - the offset into a[] and a_t[] for each auxiliary field
1670: . aOff_x - the offset into a_x[] for each auxiliary field
1671: . a - each auxiliary field evaluated at the current point
1672: . a_t - the time derivative of each auxiliary field evaluated at the current point
1673: . a_x - the gradient of auxiliary each field evaluated at the current point
1674: . t - current time
1675: . u_tShift - the multiplier a for dF/dU_t
1676: . x - coordinates of the current point
1677: . numConstants - number of constant parameters
1678: . constants - constant parameters
1679: - g0 - output values at the current point
1681: Level: intermediate
1683: .seealso: PetscDSGetJacobian()
1684: @*/
1685: PetscErrorCode PetscDSSetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g,
1686: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1687: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1688: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1689: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1690: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1691: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1692: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1693: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1694: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1695: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1696: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1697: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1698: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1699: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1700: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1701: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1702: {
1711: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1712: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1713: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1714: prob->gt[(f*prob->Nf + g)*4+0] = g0;
1715: prob->gt[(f*prob->Nf + g)*4+1] = g1;
1716: prob->gt[(f*prob->Nf + g)*4+2] = g2;
1717: prob->gt[(f*prob->Nf + g)*4+3] = g3;
1718: return(0);
1719: }
1721: /*@C
1722: PetscDSGetRiemannSolver - Returns the Riemann solver for the given field
1724: Not collective
1726: Input Arguments:
1727: + prob - The PetscDS object
1728: - f - The field number
1730: Output Argument:
1731: . r - Riemann solver
1733: Calling sequence for r:
1735: $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)
1737: + dim - The spatial dimension
1738: . Nf - The number of fields
1739: . x - The coordinates at a point on the interface
1740: . n - The normal vector to the interface
1741: . uL - The state vector to the left of the interface
1742: . uR - The state vector to the right of the interface
1743: . flux - output array of flux through the interface
1744: . numConstants - number of constant parameters
1745: . constants - constant parameters
1746: - ctx - optional user context
1748: Level: intermediate
1750: .seealso: PetscDSSetRiemannSolver()
1751: @*/
1752: PetscErrorCode PetscDSGetRiemannSolver(PetscDS prob, PetscInt f,
1753: 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))
1754: {
1757: 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);
1759: *r = prob->r[f];
1760: return(0);
1761: }
1763: /*@C
1764: PetscDSSetRiemannSolver - Sets the Riemann solver for the given field
1766: Not collective
1768: Input Arguments:
1769: + prob - The PetscDS object
1770: . f - The field number
1771: - r - Riemann solver
1773: Calling sequence for r:
1775: $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)
1777: + dim - The spatial dimension
1778: . Nf - The number of fields
1779: . x - The coordinates at a point on the interface
1780: . n - The normal vector to the interface
1781: . uL - The state vector to the left of the interface
1782: . uR - The state vector to the right of the interface
1783: . flux - output array of flux through the interface
1784: . numConstants - number of constant parameters
1785: . constants - constant parameters
1786: - ctx - optional user context
1788: Level: intermediate
1790: .seealso: PetscDSGetRiemannSolver()
1791: @*/
1792: PetscErrorCode PetscDSSetRiemannSolver(PetscDS prob, PetscInt f,
1793: 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))
1794: {
1800: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1801: PetscDSEnlarge_Static(prob, f+1);
1802: prob->r[f] = r;
1803: return(0);
1804: }
1806: /*@C
1807: PetscDSGetUpdate - Get the pointwise update function for a given field
1809: Not collective
1811: Input Parameters:
1812: + prob - The PetscDS
1813: - f - The field number
1815: Output Parameters:
1816: . update - update function
1818: Note: The calling sequence for the callback update is given by:
1820: $ update(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1821: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1822: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1823: $ PetscReal t, const PetscReal x[], PetscScalar uNew[])
1825: + dim - the spatial dimension
1826: . Nf - the number of fields
1827: . uOff - the offset into u[] and u_t[] for each field
1828: . uOff_x - the offset into u_x[] for each field
1829: . u - each field evaluated at the current point
1830: . u_t - the time derivative of each field evaluated at the current point
1831: . u_x - the gradient of each field evaluated at the current point
1832: . aOff - the offset into a[] and a_t[] for each auxiliary field
1833: . aOff_x - the offset into a_x[] for each auxiliary field
1834: . a - each auxiliary field evaluated at the current point
1835: . a_t - the time derivative of each auxiliary field evaluated at the current point
1836: . a_x - the gradient of auxiliary each field evaluated at the current point
1837: . t - current time
1838: . x - coordinates of the current point
1839: - uNew - new value for field at the current point
1841: Level: intermediate
1843: .seealso: PetscDSSetUpdate(), PetscDSSetResidual()
1844: @*/
1845: PetscErrorCode PetscDSGetUpdate(PetscDS prob, PetscInt f,
1846: void (**update)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1847: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1848: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1849: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
1850: {
1853: 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);
1855: return(0);
1856: }
1858: /*@C
1859: PetscDSSetUpdate - Set the pointwise update function for a given field
1861: Not collective
1863: Input Parameters:
1864: + prob - The PetscDS
1865: . f - The field number
1866: - update - update function
1868: Note: The calling sequence for the callback update is given by:
1870: $ update(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1871: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1872: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1873: $ PetscReal t, const PetscReal x[], PetscScalar uNew[])
1875: + dim - the spatial dimension
1876: . Nf - the number of fields
1877: . uOff - the offset into u[] and u_t[] for each field
1878: . uOff_x - the offset into u_x[] for each field
1879: . u - each field evaluated at the current point
1880: . u_t - the time derivative of each field evaluated at the current point
1881: . u_x - the gradient of each field evaluated at the current point
1882: . aOff - the offset into a[] and a_t[] for each auxiliary field
1883: . aOff_x - the offset into a_x[] for each auxiliary field
1884: . a - each auxiliary field evaluated at the current point
1885: . a_t - the time derivative of each auxiliary field evaluated at the current point
1886: . a_x - the gradient of auxiliary each field evaluated at the current point
1887: . t - current time
1888: . x - coordinates of the current point
1889: - uNew - new field values at the current point
1891: Level: intermediate
1893: .seealso: PetscDSGetResidual()
1894: @*/
1895: PetscErrorCode PetscDSSetUpdate(PetscDS prob, PetscInt f,
1896: void (*update)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1897: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1898: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1899: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
1900: {
1906: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1907: PetscDSEnlarge_Static(prob, f+1);
1908: prob->update[f] = update;
1909: return(0);
1910: }
1912: PetscErrorCode PetscDSGetContext(PetscDS prob, PetscInt f, void **ctx)
1913: {
1916: 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);
1918: *ctx = prob->ctx[f];
1919: return(0);
1920: }
1922: PetscErrorCode PetscDSSetContext(PetscDS prob, PetscInt f, void *ctx)
1923: {
1928: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1929: PetscDSEnlarge_Static(prob, f+1);
1930: prob->ctx[f] = ctx;
1931: return(0);
1932: }
1934: /*@C
1935: PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field
1937: Not collective
1939: Input Parameters:
1940: + prob - The PetscDS
1941: - f - The test field number
1943: Output Parameters:
1944: + f0 - boundary integrand for the test function term
1945: - f1 - boundary integrand for the test function gradient term
1947: Note: We are using a first order FEM model for the weak form:
1949: \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
1951: The calling sequence for the callbacks f0 and f1 is given by:
1953: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1954: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1955: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1956: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[])
1958: + dim - the spatial dimension
1959: . Nf - the number of fields
1960: . uOff - the offset into u[] and u_t[] for each field
1961: . uOff_x - the offset into u_x[] for each field
1962: . u - each field evaluated at the current point
1963: . u_t - the time derivative of each field evaluated at the current point
1964: . u_x - the gradient of each field evaluated at the current point
1965: . aOff - the offset into a[] and a_t[] for each auxiliary field
1966: . aOff_x - the offset into a_x[] for each auxiliary field
1967: . a - each auxiliary field evaluated at the current point
1968: . a_t - the time derivative of each auxiliary field evaluated at the current point
1969: . a_x - the gradient of auxiliary each field evaluated at the current point
1970: . t - current time
1971: . x - coordinates of the current point
1972: . n - unit normal at the current point
1973: . numConstants - number of constant parameters
1974: . constants - constant parameters
1975: - f0 - output values at the current point
1977: Level: intermediate
1979: .seealso: PetscDSSetBdResidual()
1980: @*/
1981: PetscErrorCode PetscDSGetBdResidual(PetscDS prob, PetscInt f,
1982: void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1983: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1984: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1985: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1986: void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1987: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1988: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1989: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1990: {
1993: 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);
1996: return(0);
1997: }
1999: /*@C
2000: PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field
2002: Not collective
2004: Input Parameters:
2005: + prob - The PetscDS
2006: . f - The test field number
2007: . f0 - boundary integrand for the test function term
2008: - f1 - boundary integrand for the test function gradient term
2010: Note: We are using a first order FEM model for the weak form:
2012: \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
2014: The calling sequence for the callbacks f0 and f1 is given by:
2016: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2017: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2018: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2019: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[])
2021: + dim - the spatial dimension
2022: . Nf - the number of fields
2023: . uOff - the offset into u[] and u_t[] for each field
2024: . uOff_x - the offset into u_x[] for each field
2025: . u - each field evaluated at the current point
2026: . u_t - the time derivative of each field evaluated at the current point
2027: . u_x - the gradient of each field evaluated at the current point
2028: . aOff - the offset into a[] and a_t[] for each auxiliary field
2029: . aOff_x - the offset into a_x[] for each auxiliary field
2030: . a - each auxiliary field evaluated at the current point
2031: . a_t - the time derivative of each auxiliary field evaluated at the current point
2032: . a_x - the gradient of auxiliary each field evaluated at the current point
2033: . t - current time
2034: . x - coordinates of the current point
2035: . n - unit normal at the current point
2036: . numConstants - number of constant parameters
2037: . constants - constant parameters
2038: - f0 - output values at the current point
2040: Level: intermediate
2042: .seealso: PetscDSGetBdResidual()
2043: @*/
2044: PetscErrorCode PetscDSSetBdResidual(PetscDS prob, PetscInt f,
2045: void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2046: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2047: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2048: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
2049: void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2050: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2051: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2052: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
2053: {
2058: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2059: PetscDSEnlarge_Static(prob, f+1);
2062: return(0);
2063: }
2065: /*@C
2066: PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field
2068: Not collective
2070: Input Parameters:
2071: + prob - The PetscDS
2072: . f - The test field number
2073: - g - The field number
2075: Output Parameters:
2076: + g0 - integrand for the test and basis function term
2077: . g1 - integrand for the test function and basis function gradient term
2078: . g2 - integrand for the test function gradient and basis function term
2079: - g3 - integrand for the test function gradient and basis function gradient term
2081: Note: We are using a first order FEM model for the weak form:
2083: \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
2085: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
2087: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2088: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2089: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2090: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[])
2092: + dim - the spatial dimension
2093: . Nf - the number of fields
2094: . uOff - the offset into u[] and u_t[] for each field
2095: . uOff_x - the offset into u_x[] for each field
2096: . u - each field evaluated at the current point
2097: . u_t - the time derivative of each field evaluated at the current point
2098: . u_x - the gradient of each field evaluated at the current point
2099: . aOff - the offset into a[] and a_t[] for each auxiliary field
2100: . aOff_x - the offset into a_x[] for each auxiliary field
2101: . a - each auxiliary field evaluated at the current point
2102: . a_t - the time derivative of each auxiliary field evaluated at the current point
2103: . a_x - the gradient of auxiliary each field evaluated at the current point
2104: . t - current time
2105: . u_tShift - the multiplier a for dF/dU_t
2106: . x - coordinates of the current point
2107: . n - normal at the current point
2108: . numConstants - number of constant parameters
2109: . constants - constant parameters
2110: - g0 - output values at the current point
2112: Level: intermediate
2114: .seealso: PetscDSSetBdJacobian()
2115: @*/
2116: PetscErrorCode PetscDSGetBdJacobian(PetscDS prob, PetscInt f, PetscInt g,
2117: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2118: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2119: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2120: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
2121: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2122: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2123: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2124: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
2125: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2126: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2127: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2128: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
2129: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2130: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2131: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2132: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2133: {
2136: 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);
2137: 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);
2142: return(0);
2143: }
2145: /*@C
2146: PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field
2148: Not collective
2150: Input Parameters:
2151: + prob - The PetscDS
2152: . f - The test field number
2153: . g - The field number
2154: . g0 - integrand for the test and basis function term
2155: . g1 - integrand for the test function and basis function gradient term
2156: . g2 - integrand for the test function gradient and basis function term
2157: - g3 - integrand for the test function gradient and basis function gradient term
2159: Note: We are using a first order FEM model for the weak form:
2161: \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
2163: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
2165: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2166: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2167: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2168: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[])
2170: + dim - the spatial dimension
2171: . Nf - the number of fields
2172: . uOff - the offset into u[] and u_t[] for each field
2173: . uOff_x - the offset into u_x[] for each field
2174: . u - each field evaluated at the current point
2175: . u_t - the time derivative of each field evaluated at the current point
2176: . u_x - the gradient of each field evaluated at the current point
2177: . aOff - the offset into a[] and a_t[] for each auxiliary field
2178: . aOff_x - the offset into a_x[] for each auxiliary field
2179: . a - each auxiliary field evaluated at the current point
2180: . a_t - the time derivative of each auxiliary field evaluated at the current point
2181: . a_x - the gradient of auxiliary each field evaluated at the current point
2182: . t - current time
2183: . u_tShift - the multiplier a for dF/dU_t
2184: . x - coordinates of the current point
2185: . n - normal at the current point
2186: . numConstants - number of constant parameters
2187: . constants - constant parameters
2188: - g0 - output values at the current point
2190: Level: intermediate
2192: .seealso: PetscDSGetBdJacobian()
2193: @*/
2194: PetscErrorCode PetscDSSetBdJacobian(PetscDS prob, PetscInt f, PetscInt g,
2195: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2196: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2197: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2198: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
2199: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2200: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2201: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2202: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
2203: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2204: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2205: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2206: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
2207: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2208: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2209: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2210: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2211: {
2220: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2221: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
2222: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
2223: prob->gBd[(f*prob->Nf + g)*4+0] = g0;
2224: prob->gBd[(f*prob->Nf + g)*4+1] = g1;
2225: prob->gBd[(f*prob->Nf + g)*4+2] = g2;
2226: prob->gBd[(f*prob->Nf + g)*4+3] = g3;
2227: return(0);
2228: }
2230: /*@C
2231: PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field
2233: Not collective
2235: Input Parameters:
2236: + prob - The PetscDS
2237: - f - The test field number
2239: Output Parameter:
2240: + exactSol - exact solution for the test field
2241: - exactCtx - exact solution context
2243: Note: The calling sequence for the solution functions is given by:
2245: $ sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2247: + dim - the spatial dimension
2248: . t - current time
2249: . x - coordinates of the current point
2250: . Nc - the number of field components
2251: . u - the solution field evaluated at the current point
2252: - ctx - a user context
2254: Level: intermediate
2256: .seealso: PetscDSSetExactSolution()
2257: @*/
2258: PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx)
2259: {
2262: 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);
2265: return(0);
2266: }
2268: /*@C
2269: PetscDSSetExactSolution - Set the pointwise exact solution function for a given test field
2271: Not collective
2273: Input Parameters:
2274: + prob - The PetscDS
2275: . f - The test field number
2276: . sol - solution function for the test fields
2277: - ctx - solution context or NULL
2279: Note: The calling sequence for solution functions is given by:
2281: $ sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2283: + dim - the spatial dimension
2284: . t - current time
2285: . x - coordinates of the current point
2286: . Nc - the number of field components
2287: . u - the solution field evaluated at the current point
2288: - ctx - a user context
2290: Level: intermediate
2292: .seealso: PetscDSGetExactSolution()
2293: @*/
2294: PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx)
2295: {
2300: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2301: PetscDSEnlarge_Static(prob, f+1);
2304: return(0);
2305: }
2307: /*@C
2308: PetscDSGetConstants - Returns the array of constants passed to point functions
2310: Not collective
2312: Input Parameter:
2313: . prob - The PetscDS object
2315: Output Parameters:
2316: + numConstants - The number of constants
2317: - constants - The array of constants, NULL if there are none
2319: Level: intermediate
2321: .seealso: PetscDSSetConstants(), PetscDSCreate()
2322: @*/
2323: PetscErrorCode PetscDSGetConstants(PetscDS prob, PetscInt *numConstants, const PetscScalar *constants[])
2324: {
2329: return(0);
2330: }
2332: /*@C
2333: PetscDSSetConstants - Set the array of constants passed to point functions
2335: Not collective
2337: Input Parameters:
2338: + prob - The PetscDS object
2339: . numConstants - The number of constants
2340: - constants - The array of constants, NULL if there are none
2342: Level: intermediate
2344: .seealso: PetscDSGetConstants(), PetscDSCreate()
2345: @*/
2346: PetscErrorCode PetscDSSetConstants(PetscDS prob, PetscInt numConstants, PetscScalar constants[])
2347: {
2352: if (numConstants != prob->numConstants) {
2353: PetscFree(prob->constants);
2354: prob->numConstants = numConstants;
2355: if (prob->numConstants) {
2356: PetscMalloc1(prob->numConstants, &prob->constants);
2357: } else {
2358: prob->constants = NULL;
2359: }
2360: }
2361: if (prob->numConstants) {
2363: PetscArraycpy(prob->constants, constants, prob->numConstants);
2364: }
2365: return(0);
2366: }
2368: /*@
2369: PetscDSGetFieldIndex - Returns the index of the given field
2371: Not collective
2373: Input Parameters:
2374: + prob - The PetscDS object
2375: - disc - The discretization object
2377: Output Parameter:
2378: . f - The field number
2380: Level: beginner
2382: .seealso: PetscGetDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
2383: @*/
2384: PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f)
2385: {
2386: PetscInt g;
2391: *f = -1;
2392: for (g = 0; g < prob->Nf; ++g) {if (disc == prob->disc[g]) break;}
2393: if (g == prob->Nf) SETERRQ(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS.");
2394: *f = g;
2395: return(0);
2396: }
2398: /*@
2399: PetscDSGetFieldSize - Returns the size of the given field in the full space basis
2401: Not collective
2403: Input Parameters:
2404: + prob - The PetscDS object
2405: - f - The field number
2407: Output Parameter:
2408: . size - The size
2410: Level: beginner
2412: .seealso: PetscDSGetFieldOffset(), PetscDSGetNumFields(), PetscDSCreate()
2413: @*/
2414: PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size)
2415: {
2421: 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);
2422: PetscDSSetUp(prob);
2423: *size = prob->Nb[f];
2424: return(0);
2425: }
2427: /*@
2428: PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis
2430: Not collective
2432: Input Parameters:
2433: + prob - The PetscDS object
2434: - f - The field number
2436: Output Parameter:
2437: . off - The offset
2439: Level: beginner
2441: .seealso: PetscDSGetFieldSize(), PetscDSGetNumFields(), PetscDSCreate()
2442: @*/
2443: PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off)
2444: {
2445: PetscInt size, g;
2451: 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);
2452: *off = 0;
2453: for (g = 0; g < f; ++g) {
2454: PetscDSGetFieldSize(prob, g, &size);
2455: *off += size;
2456: }
2457: return(0);
2458: }
2460: /*@
2461: PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point
2463: Not collective
2465: Input Parameter:
2466: . prob - The PetscDS object
2468: Output Parameter:
2469: . dimensions - The number of dimensions
2471: Level: beginner
2473: .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate()
2474: @*/
2475: PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[])
2476: {
2481: PetscDSSetUp(prob);
2483: *dimensions = prob->Nb;
2484: return(0);
2485: }
2487: /*@
2488: PetscDSGetComponents - Returns the number of components for each field on an evaluation point
2490: Not collective
2492: Input Parameter:
2493: . prob - The PetscDS object
2495: Output Parameter:
2496: . components - The number of components
2498: Level: beginner
2500: .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate()
2501: @*/
2502: PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[])
2503: {
2508: PetscDSSetUp(prob);
2510: *components = prob->Nc;
2511: return(0);
2512: }
2514: /*@
2515: PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point
2517: Not collective
2519: Input Parameters:
2520: + prob - The PetscDS object
2521: - f - The field number
2523: Output Parameter:
2524: . off - The offset
2526: Level: beginner
2528: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2529: @*/
2530: PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off)
2531: {
2535: 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);
2536: *off = prob->off[f];
2537: return(0);
2538: }
2540: /*@
2541: PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point
2543: Not collective
2545: Input Parameter:
2546: . prob - The PetscDS object
2548: Output Parameter:
2549: . offsets - The offsets
2551: Level: beginner
2553: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2554: @*/
2555: PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[])
2556: {
2560: *offsets = prob->off;
2561: return(0);
2562: }
2564: /*@
2565: PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point
2567: Not collective
2569: Input Parameter:
2570: . prob - The PetscDS object
2572: Output Parameter:
2573: . offsets - The offsets
2575: Level: beginner
2577: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2578: @*/
2579: PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[])
2580: {
2584: *offsets = prob->offDer;
2585: return(0);
2586: }
2588: /*@C
2589: PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization
2591: Not collective
2593: Input Parameter:
2594: . prob - The PetscDS object
2596: Output Parameter:
2597: . T - The basis function and derivatives tabulation at quadrature points for each field
2599: Level: intermediate
2601: .seealso: PetscDSCreate()
2602: @*/
2603: PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscTabulation *T[])
2604: {
2610: PetscDSSetUp(prob);
2611: *T = prob->T;
2612: return(0);
2613: }
2615: /*@C
2616: PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces
2618: Not collective
2620: Input Parameter:
2621: . prob - The PetscDS object
2623: Output Parameter:
2624: . Tf - The basis function and derviative tabulation on each local face at quadrature points for each and field
2626: Level: intermediate
2628: .seealso: PetscDSGetTabulation(), PetscDSCreate()
2629: @*/
2630: PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscTabulation *Tf[])
2631: {
2637: PetscDSSetUp(prob);
2638: *Tf = prob->Tf;
2639: return(0);
2640: }
2642: PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x)
2643: {
2648: PetscDSSetUp(prob);
2652: return(0);
2653: }
2655: PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3)
2656: {
2661: PetscDSSetUp(prob);
2668: return(0);
2669: }
2671: PetscErrorCode PetscDSGetWorkspace(PetscDS prob, PetscReal **x, PetscScalar **basisReal, PetscScalar **basisDerReal, PetscScalar **testReal, PetscScalar **testDerReal)
2672: {
2677: PetscDSSetUp(prob);
2683: return(0);
2684: }
2686: /*@C
2687: PetscDSAddBoundary - Add a boundary condition to the model
2689: Input Parameters:
2690: + ds - The PetscDS object
2691: . type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2692: . name - The BC name
2693: . labelname - The label defining constrained points
2694: . field - The field to constrain
2695: . numcomps - The number of constrained field components (0 will constrain all fields)
2696: . comps - An array of constrained component numbers
2697: . bcFunc - A pointwise function giving boundary values
2698: . numids - The number of DMLabel ids for constrained points
2699: . ids - An array of ids for constrained points
2700: - ctx - An optional user context for bcFunc
2702: Options Database Keys:
2703: + -bc_<boundary name> <num> - Overrides the boundary ids
2704: - -bc_<boundary name>_comp <num> - Overrides the boundary components
2706: Level: developer
2708: .seealso: PetscDSGetBoundary()
2709: @*/
2710: 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)
2711: {
2712: DSBoundary b;
2717: PetscNew(&b);
2718: PetscStrallocpy(name, (char **) &b->name);
2719: PetscStrallocpy(labelname, (char **) &b->labelname);
2720: PetscMalloc1(numcomps, &b->comps);
2721: if (numcomps) {PetscArraycpy(b->comps, comps, numcomps);}
2722: PetscMalloc1(numids, &b->ids);
2723: if (numids) {PetscArraycpy(b->ids, ids, numids);}
2724: b->type = type;
2725: b->field = field;
2726: b->numcomps = numcomps;
2727: b->func = bcFunc;
2728: b->numids = numids;
2729: b->ctx = ctx;
2730: b->next = ds->boundary;
2731: ds->boundary = b;
2732: return(0);
2733: }
2735: /*@C
2736: PetscDSUpdateBoundary - Change a boundary condition for the model
2738: Input Parameters:
2739: + ds - The PetscDS object
2740: . bd - The boundary condition number
2741: . type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2742: . name - The BC name
2743: . labelname - The label defining constrained points
2744: . field - The field to constrain
2745: . numcomps - The number of constrained field components
2746: . comps - An array of constrained component numbers
2747: . bcFunc - A pointwise function giving boundary values
2748: . numids - The number of DMLabel ids for constrained points
2749: . ids - An array of ids for constrained points
2750: - ctx - An optional user context for bcFunc
2752: Note: The boundary condition number is the order in which it was registered. The user can get the number of boundary conditions from PetscDSGetNumBoundary().
2754: Level: developer
2756: .seealso: PetscDSAddBoundary(), PetscDSGetBoundary(), PetscDSGetNumBoundary()
2757: @*/
2758: 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)
2759: {
2760: DSBoundary b = ds->boundary;
2761: PetscInt n = 0;
2766: while (b) {
2767: if (n == bd) break;
2768: b = b->next;
2769: ++n;
2770: }
2771: if (!b) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %d is not in [0, %d)", bd, n);
2772: if (name) {
2773: PetscFree(b->name);
2774: PetscStrallocpy(name, (char **) &b->name);
2775: }
2776: if (labelname) {
2777: PetscFree(b->labelname);
2778: PetscStrallocpy(labelname, (char **) &b->labelname);
2779: }
2780: if (numcomps >= 0 && numcomps != b->numcomps) {
2781: b->numcomps = numcomps;
2782: PetscFree(b->comps);
2783: PetscMalloc1(numcomps, &b->comps);
2784: if (numcomps) {PetscArraycpy(b->comps, comps, numcomps);}
2785: }
2786: if (numids >= 0 && numids != b->numids) {
2787: b->numids = numids;
2788: PetscFree(b->ids);
2789: PetscMalloc1(numids, &b->ids);
2790: if (numids) {PetscArraycpy(b->ids, ids, numids);}
2791: }
2792: b->type = type;
2793: if (field >= 0) {b->field = field;}
2794: if (bcFunc) {b->func = bcFunc;}
2795: if (ctx) {b->ctx = ctx;}
2796: return(0);
2797: }
2799: /*@
2800: PetscDSGetNumBoundary - Get the number of registered BC
2802: Input Parameters:
2803: . ds - The PetscDS object
2805: Output Parameters:
2806: . numBd - The number of BC
2808: Level: intermediate
2810: .seealso: PetscDSAddBoundary(), PetscDSGetBoundary()
2811: @*/
2812: PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd)
2813: {
2814: DSBoundary b = ds->boundary;
2819: *numBd = 0;
2820: while (b) {++(*numBd); b = b->next;}
2821: return(0);
2822: }
2824: /*@C
2825: PetscDSGetBoundary - Gets a boundary condition to the model
2827: Input Parameters:
2828: + ds - The PetscDS object
2829: - bd - The BC number
2831: Output Parameters:
2832: + type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2833: . name - The BC name
2834: . labelname - The label defining constrained points
2835: . field - The field to constrain
2836: . numcomps - The number of constrained field components
2837: . comps - An array of constrained component numbers
2838: . bcFunc - A pointwise function giving boundary values
2839: . numids - The number of DMLabel ids for constrained points
2840: . ids - An array of ids for constrained points
2841: - ctx - An optional user context for bcFunc
2843: Options Database Keys:
2844: + -bc_<boundary name> <num> - Overrides the boundary ids
2845: - -bc_<boundary name>_comp <num> - Overrides the boundary components
2847: Level: developer
2849: .seealso: PetscDSAddBoundary()
2850: @*/
2851: 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)
2852: {
2853: DSBoundary b = ds->boundary;
2854: PetscInt n = 0;
2858: while (b) {
2859: if (n == bd) break;
2860: b = b->next;
2861: ++n;
2862: }
2863: if (!b) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %d is not in [0, %d)", bd, n);
2864: if (type) {
2866: *type = b->type;
2867: }
2868: if (name) {
2870: *name = b->name;
2871: }
2872: if (labelname) {
2874: *labelname = b->labelname;
2875: }
2876: if (field) {
2878: *field = b->field;
2879: }
2880: if (numcomps) {
2882: *numcomps = b->numcomps;
2883: }
2884: if (comps) {
2886: *comps = b->comps;
2887: }
2888: if (func) {
2890: *func = b->func;
2891: }
2892: if (numids) {
2894: *numids = b->numids;
2895: }
2896: if (ids) {
2898: *ids = b->ids;
2899: }
2900: if (ctx) {
2902: *ctx = b->ctx;
2903: }
2904: return(0);
2905: }
2907: /*@
2908: PetscDSCopyBoundary - Copy all boundary condition objects to the new problem
2910: Not collective
2912: Input Parameter:
2913: . prob - The PetscDS object
2915: Output Parameter:
2916: . newprob - The PetscDS copy
2918: Level: intermediate
2920: .seealso: PetscDSCopyEquations(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
2921: @*/
2922: PetscErrorCode PetscDSCopyBoundary(PetscDS probA, PetscDS probB)
2923: {
2924: DSBoundary b, next, *lastnext;
2930: if (probA == probB) return(0);
2931: next = probB->boundary;
2932: while (next) {
2933: DSBoundary b = next;
2935: next = b->next;
2936: PetscFree(b->comps);
2937: PetscFree(b->ids);
2938: PetscFree(b->name);
2939: PetscFree(b->labelname);
2940: PetscFree(b);
2941: }
2942: lastnext = &(probB->boundary);
2943: for (b = probA->boundary; b; b = b->next) {
2944: DSBoundary bNew;
2946: PetscNew(&bNew);
2947: bNew->numcomps = b->numcomps;
2948: PetscMalloc1(bNew->numcomps, &bNew->comps);
2949: PetscArraycpy(bNew->comps, b->comps, bNew->numcomps);
2950: bNew->numids = b->numids;
2951: PetscMalloc1(bNew->numids, &bNew->ids);
2952: PetscArraycpy(bNew->ids, b->ids, bNew->numids);
2953: PetscStrallocpy(b->labelname,(char **) &(bNew->labelname));
2954: PetscStrallocpy(b->name,(char **) &(bNew->name));
2955: bNew->ctx = b->ctx;
2956: bNew->type = b->type;
2957: bNew->field = b->field;
2958: bNew->func = b->func;
2960: *lastnext = bNew;
2961: lastnext = &(bNew->next);
2962: }
2963: return(0);
2964: }
2966: /*@C
2967: PetscDSSelectEquations - Copy pointwise function pointers to the new problem with different field layout
2969: Not collective
2971: Input Parameter:
2972: + prob - The PetscDS object
2973: . numFields - Number of new fields
2974: - fields - Old field number for each new field
2976: Output Parameter:
2977: . newprob - The PetscDS copy
2979: Level: intermediate
2981: .seealso: PetscDSCopyBoundary(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
2982: @*/
2983: PetscErrorCode PetscDSSelectEquations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
2984: {
2985: PetscInt Nf, Nfn, fn, gn;
2992: PetscDSGetNumFields(prob, &Nf);
2993: PetscDSGetNumFields(newprob, &Nfn);
2994: 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);
2995: for (fn = 0; fn < numFields; ++fn) {
2996: const PetscInt f = fields ? fields[fn] : fn;
2997: PetscPointFunc obj;
2998: PetscPointFunc f0, f1;
2999: PetscBdPointFunc f0Bd, f1Bd;
3000: PetscRiemannFunc r;
3002: if (f >= Nf) continue;
3003: PetscDSGetObjective(prob, f, &obj);
3004: PetscDSGetResidual(prob, f, &f0, &f1);
3005: PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd);
3006: PetscDSGetRiemannSolver(prob, f, &r);
3007: PetscDSSetObjective(newprob, fn, obj);
3008: PetscDSSetResidual(newprob, fn, f0, f1);
3009: PetscDSSetBdResidual(newprob, fn, f0Bd, f1Bd);
3010: PetscDSSetRiemannSolver(newprob, fn, r);
3011: for (gn = 0; gn < numFields; ++gn) {
3012: const PetscInt g = fields ? fields[gn] : gn;
3013: PetscPointJac g0, g1, g2, g3;
3014: PetscPointJac g0p, g1p, g2p, g3p;
3015: PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd;
3017: if (g >= Nf) continue;
3018: PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3);
3019: PetscDSGetJacobianPreconditioner(prob, f, g, &g0p, &g1p, &g2p, &g3p);
3020: PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd);
3021: PetscDSSetJacobian(newprob, fn, gn, g0, g1, g2, g3);
3022: PetscDSSetJacobianPreconditioner(prob, fn, gn, g0p, g1p, g2p, g3p);
3023: PetscDSSetBdJacobian(newprob, fn, gn, g0Bd, g1Bd, g2Bd, g3Bd);
3024: }
3025: }
3026: return(0);
3027: }
3029: /*@
3030: PetscDSCopyEquations - Copy all pointwise function pointers to the new problem
3032: Not collective
3034: Input Parameter:
3035: . prob - The PetscDS object
3037: Output Parameter:
3038: . newprob - The PetscDS copy
3040: Level: intermediate
3042: .seealso: PetscDSCopyBoundary(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
3043: @*/
3044: PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob)
3045: {
3046: PetscInt Nf, Ng;
3052: PetscDSGetNumFields(prob, &Nf);
3053: PetscDSGetNumFields(newprob, &Ng);
3054: if (Nf != Ng) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %D != %D", Nf, Ng);
3055: PetscDSSelectEquations(prob, Nf, NULL, newprob);
3056: return(0);
3057: }
3058: /*@
3059: PetscDSCopyConstants - Copy all constants to the new problem
3061: Not collective
3063: Input Parameter:
3064: . prob - The PetscDS object
3066: Output Parameter:
3067: . newprob - The PetscDS copy
3069: Level: intermediate
3071: .seealso: PetscDSCopyBoundary(), PetscDSCopyEquations(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
3072: @*/
3073: PetscErrorCode PetscDSCopyConstants(PetscDS prob, PetscDS newprob)
3074: {
3075: PetscInt Nc;
3076: const PetscScalar *constants;
3077: PetscErrorCode ierr;
3082: PetscDSGetConstants(prob, &Nc, &constants);
3083: PetscDSSetConstants(newprob, Nc, (PetscScalar *) constants);
3084: return(0);
3085: }
3087: PetscErrorCode PetscDSGetHeightSubspace(PetscDS prob, PetscInt height, PetscDS *subprob)
3088: {
3089: PetscInt dim, Nf, f;
3095: if (height == 0) {*subprob = prob; return(0);}
3096: PetscDSGetNumFields(prob, &Nf);
3097: PetscDSGetSpatialDimension(prob, &dim);
3098: if (height > dim) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_OUTOFRANGE, "DS can only handle height in [0, %D], not %D", dim, height);
3099: if (!prob->subprobs) {PetscCalloc1(dim, &prob->subprobs);}
3100: if (!prob->subprobs[height-1]) {
3101: PetscInt cdim;
3103: PetscDSCreate(PetscObjectComm((PetscObject) prob), &prob->subprobs[height-1]);
3104: PetscDSGetCoordinateDimension(prob, &cdim);
3105: PetscDSSetCoordinateDimension(prob->subprobs[height-1], cdim);
3106: for (f = 0; f < Nf; ++f) {
3107: PetscFE subfe;
3108: PetscObject obj;
3109: PetscClassId id;
3111: PetscDSGetDiscretization(prob, f, &obj);
3112: PetscObjectGetClassId(obj, &id);
3113: if (id == PETSCFE_CLASSID) {PetscFEGetHeightSubspace((PetscFE) obj, height, &subfe);}
3114: else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unsupported discretization type for field %d", f);
3115: PetscDSSetDiscretization(prob->subprobs[height-1], f, (PetscObject) subfe);
3116: }
3117: }
3118: *subprob = prob->subprobs[height-1];
3119: return(0);
3120: }
3122: PetscErrorCode PetscDSIsFE_Internal(PetscDS ds, PetscInt f, PetscBool *isFE)
3123: {
3124: PetscObject obj;
3125: PetscClassId id;
3126: PetscInt Nf;
3132: PetscDSGetNumFields(ds, &Nf);
3133: if (f >= Nf) SETERRQ2(PetscObjectComm((PetscObject) ds), PETSC_ERR_ARG_SIZ, "Field %D must be in [0, %D)", f, Nf);
3134: PetscDSGetDiscretization(ds, f, &obj);
3135: PetscObjectGetClassId(obj, &id);
3136: if (id == PETSCFE_CLASSID) *isFE = PETSC_TRUE;
3137: else *isFE = PETSC_FALSE;
3138: return(0);
3139: }
3141: static PetscErrorCode PetscDSDestroy_Basic(PetscDS prob)
3142: {
3143: PetscErrorCode ierr;
3146: PetscFree(prob->data);
3147: return(0);
3148: }
3150: static PetscErrorCode PetscDSInitialize_Basic(PetscDS prob)
3151: {
3153: prob->ops->setfromoptions = NULL;
3154: prob->ops->setup = NULL;
3155: prob->ops->view = NULL;
3156: prob->ops->destroy = PetscDSDestroy_Basic;
3157: return(0);
3158: }
3160: /*MC
3161: PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions
3163: Level: intermediate
3165: .seealso: PetscDSType, PetscDSCreate(), PetscDSSetType()
3166: M*/
3168: PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS prob)
3169: {
3170: PetscDS_Basic *b;
3171: PetscErrorCode ierr;
3175: PetscNewLog(prob, &b);
3176: prob->data = b;
3178: PetscDSInitialize_Basic(prob);
3179: return(0);
3180: }