Actual source code: dtds.c
petsc-3.11.4 2019-09-28
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: .keywords: PetscDS, register
53: .seealso: PetscDSRegisterAll(), PetscDSRegisterDestroy()
55: @*/
56: PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS))
57: {
61: PetscFunctionListAdd(&PetscDSList, sname, function);
62: return(0);
63: }
65: /*@C
66: PetscDSSetType - Builds a particular PetscDS
68: Collective on PetscDS
70: Input Parameters:
71: + prob - The PetscDS object
72: - name - The kind of system
74: Options Database Key:
75: . -petscds_type <type> - Sets the PetscDS type; use -help for a list of available types
77: Level: intermediate
79: Not available from Fortran
81: .keywords: PetscDS, set, type
82: .seealso: PetscDSGetType(), PetscDSCreate()
83: @*/
84: PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name)
85: {
86: PetscErrorCode (*r)(PetscDS);
87: PetscBool match;
92: PetscObjectTypeCompare((PetscObject) prob, name, &match);
93: if (match) return(0);
95: PetscDSRegisterAll();
96: PetscFunctionListFind(PetscDSList, name, &r);
97: if (!r) SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscDS type: %s", name);
99: if (prob->ops->destroy) {
100: (*prob->ops->destroy)(prob);
101: prob->ops->destroy = NULL;
102: }
103: (*r)(prob);
104: PetscObjectChangeTypeName((PetscObject) prob, name);
105: return(0);
106: }
108: /*@C
109: PetscDSGetType - Gets the PetscDS type name (as a string) from the object.
111: Not Collective
113: Input Parameter:
114: . prob - The PetscDS
116: Output Parameter:
117: . name - The PetscDS type name
119: Level: intermediate
121: Not available from Fortran
123: .keywords: PetscDS, get, type, name
124: .seealso: PetscDSSetType(), PetscDSCreate()
125: @*/
126: PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name)
127: {
133: PetscDSRegisterAll();
134: *name = ((PetscObject) prob)->type_name;
135: return(0);
136: }
138: static PetscErrorCode PetscDSView_Ascii(PetscDS prob, PetscViewer viewer)
139: {
140: PetscViewerFormat format;
141: const PetscScalar *constants;
142: PetscInt numConstants, f;
143: PetscErrorCode ierr;
146: PetscViewerGetFormat(viewer, &format);
147: PetscViewerASCIIPrintf(viewer, "Discrete System with %d fields\n", prob->Nf);
148: PetscViewerASCIIPushTab(viewer);
149: PetscViewerASCIIPrintf(viewer, " cell total dim %D total comp %D\n", prob->totDim, prob->totComp);
150: if (prob->isHybrid) {PetscViewerASCIIPrintf(viewer, " hybrid cell\n");}
151: for (f = 0; f < prob->Nf; ++f) {
152: DSBoundary b;
153: PetscObject obj;
154: PetscClassId id;
155: PetscQuadrature q;
156: const char *name;
157: PetscInt Nc, Nq, Nqc;
159: PetscDSGetDiscretization(prob, f, &obj);
160: PetscObjectGetClassId(obj, &id);
161: PetscObjectGetName(obj, &name);
162: PetscViewerASCIIPrintf(viewer, "Field %s", name ? name : "<unknown>");
163: PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
164: if (id == PETSCFE_CLASSID) {
165: PetscFEGetNumComponents((PetscFE) obj, &Nc);
166: PetscFEGetQuadrature((PetscFE) obj, &q);
167: PetscViewerASCIIPrintf(viewer, " FEM");
168: } else if (id == PETSCFV_CLASSID) {
169: PetscFVGetNumComponents((PetscFV) obj, &Nc);
170: PetscFVGetQuadrature((PetscFV) obj, &q);
171: PetscViewerASCIIPrintf(viewer, " FVM");
172: }
173: else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %D", f);
174: if (Nc > 1) {PetscViewerASCIIPrintf(viewer, "%D components", Nc);}
175: else {PetscViewerASCIIPrintf(viewer, "%D component ", Nc);}
176: if (prob->implicit[f]) {PetscViewerASCIIPrintf(viewer, " (implicit)");}
177: else {PetscViewerASCIIPrintf(viewer, " (explicit)");}
178: if (q) {
179: PetscQuadratureGetData(q, NULL, &Nqc, &Nq, NULL, NULL);
180: PetscViewerASCIIPrintf(viewer, " (Nq %D Nqc %D)", Nq, Nqc);
181: }
182: PetscViewerASCIIPrintf(viewer, "\n");
183: PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
184: PetscViewerASCIIPushTab(viewer);
185: if (id == PETSCFE_CLASSID) {PetscFEView((PetscFE) obj, viewer);}
186: else if (id == PETSCFV_CLASSID) {PetscFVView((PetscFV) obj, viewer);}
187: PetscViewerASCIIPopTab(viewer);
189: for (b = prob->boundary; b; b = b->next) {
190: PetscInt c, i;
192: if (b->field != f) continue;
193: PetscViewerASCIIPushTab(viewer);
194: PetscViewerASCIIPrintf(viewer, "Boundary %s (%s) %s\n", b->name, b->labelname, DMBoundaryConditionTypes[b->type]);
195: if (!b->numcomps) {
196: PetscViewerASCIIPrintf(viewer, " all components\n");
197: } else {
198: PetscViewerASCIIPrintf(viewer, " components: ");
199: PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
200: for (c = 0; c < b->numcomps; ++c) {
201: if (c > 0) {PetscViewerASCIIPrintf(viewer, ", ");}
202: PetscViewerASCIIPrintf(viewer, "%D", b->comps[c]);
203: }
204: PetscViewerASCIIPrintf(viewer, "\n");
205: PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
206: }
207: PetscViewerASCIIPrintf(viewer, " ids: ");
208: PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
209: for (i = 0; i < b->numids; ++i) {
210: if (i > 0) {PetscViewerASCIIPrintf(viewer, ", ");}
211: PetscViewerASCIIPrintf(viewer, "%D", b->ids[i]);
212: }
213: PetscViewerASCIIPrintf(viewer, "\n");
214: PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
215: PetscViewerASCIIPopTab(viewer);
216: }
217: }
218: PetscDSGetConstants(prob, &numConstants, &constants);
219: if (numConstants) {
220: PetscViewerASCIIPrintf(viewer, "%D constants\n", numConstants);
221: PetscViewerASCIIPushTab(viewer);
222: for (f = 0; f < numConstants; ++f) {PetscViewerASCIIPrintf(viewer, "%g\n", (double) PetscRealPart(constants[f]));}
223: PetscViewerASCIIPopTab(viewer);
224: }
225: PetscViewerASCIIPopTab(viewer);
226: return(0);
227: }
229: /*@C
230: PetscDSView - Views a PetscDS
232: Collective on PetscDS
234: Input Parameter:
235: + prob - the PetscDS object to view
236: - v - the viewer
238: Level: developer
240: .seealso PetscDSDestroy()
241: @*/
242: PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v)
243: {
244: PetscBool iascii;
249: if (!v) {PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject) prob), &v);}
251: PetscObjectTypeCompare((PetscObject) v, PETSCVIEWERASCII, &iascii);
252: if (iascii) {PetscDSView_Ascii(prob, v);}
253: if (prob->ops->view) {(*prob->ops->view)(prob, v);}
254: return(0);
255: }
257: /*@
258: PetscDSSetFromOptions - sets parameters in a PetscDS from the options database
260: Collective on PetscDS
262: Input Parameter:
263: . prob - the PetscDS object to set options for
265: Options Database:
266: + -petscds_type <type> : Set the DS type
267: . -petscds_view <view opt> : View the DS
268: . -petscds_jac_pre : Turn formation of a separate Jacobian preconditioner on and off
269: . -bc_<name> <ids> : Specify a list of label ids for a boundary condition
270: - -bc_<name>_comp <comps> : Specify a list of field components to constrain for a boundary condition
272: Level: developer
274: .seealso PetscDSView()
275: @*/
276: PetscErrorCode PetscDSSetFromOptions(PetscDS prob)
277: {
278: DSBoundary b;
279: const char *defaultType;
280: char name[256];
281: PetscBool flg;
286: if (!((PetscObject) prob)->type_name) {
287: defaultType = PETSCDSBASIC;
288: } else {
289: defaultType = ((PetscObject) prob)->type_name;
290: }
291: PetscDSRegisterAll();
293: PetscObjectOptionsBegin((PetscObject) prob);
294: for (b = prob->boundary; b; b = b->next) {
295: char optname[1024];
296: PetscInt ids[1024], len = 1024;
297: PetscBool flg;
299: PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name);
300: PetscMemzero(ids, sizeof(ids));
301: PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg);
302: if (flg) {
303: b->numids = len;
304: PetscFree(b->ids);
305: PetscMalloc1(len, &b->ids);
306: PetscMemcpy(b->ids, ids, len*sizeof(PetscInt));
307: }
308: len = 1024;
309: PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name);
310: PetscMemzero(ids, sizeof(ids));
311: PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg);
312: if (flg) {
313: b->numcomps = len;
314: PetscFree(b->comps);
315: PetscMalloc1(len, &b->comps);
316: PetscMemcpy(b->comps, ids, len*sizeof(PetscInt));
317: }
318: }
319: PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg);
320: if (flg) {
321: PetscDSSetType(prob, name);
322: } else if (!((PetscObject) prob)->type_name) {
323: PetscDSSetType(prob, defaultType);
324: }
325: PetscOptionsBool("-petscds_jac_pre", "Discrete System", "PetscDSUseJacobianPreconditioner", prob->useJacPre, &prob->useJacPre, &flg);
326: if (prob->ops->setfromoptions) {(*prob->ops->setfromoptions)(prob);}
327: /* process any options handlers added with PetscObjectAddOptionsHandler() */
328: PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject) prob);
329: PetscOptionsEnd();
330: if (prob->Nf) {PetscDSViewFromOptions(prob, NULL, "-petscds_view");}
331: return(0);
332: }
334: /*@C
335: PetscDSSetUp - Construct data structures for the PetscDS
337: Collective on PetscDS
339: Input Parameter:
340: . prob - the PetscDS object to setup
342: Level: developer
344: .seealso PetscDSView(), PetscDSDestroy()
345: @*/
346: PetscErrorCode PetscDSSetUp(PetscDS prob)
347: {
348: const PetscInt Nf = prob->Nf;
349: PetscInt dim, dimEmbed, work, NcMax = 0, NqMax = 0, NsMax = 1, f;
354: if (prob->setup) return(0);
355: /* Calculate sizes */
356: PetscDSGetSpatialDimension(prob, &dim);
357: PetscDSGetCoordinateDimension(prob, &dimEmbed);
358: prob->totDim = prob->totComp = 0;
359: PetscMalloc2(Nf,&prob->Nc,Nf,&prob->Nb);
360: PetscCalloc2(Nf+1,&prob->off,Nf+1,&prob->offDer);
361: PetscMalloc4(Nf,&prob->basis,Nf,&prob->basisDer,Nf,&prob->basisFace,Nf,&prob->basisDerFace);
362: for (f = 0; f < Nf; ++f) {
363: PetscObject obj;
364: PetscClassId id;
365: PetscQuadrature q;
366: PetscInt Nq = 0, Nb, Nc;
368: PetscDSGetDiscretization(prob, f, &obj);
369: PetscObjectGetClassId(obj, &id);
370: if (id == PETSCFE_CLASSID) {
371: PetscFE fe = (PetscFE) obj;
373: PetscFEGetQuadrature(fe, &q);
374: PetscFEGetDimension(fe, &Nb);
375: PetscFEGetNumComponents(fe, &Nc);
376: PetscFEGetDefaultTabulation(fe, &prob->basis[f], &prob->basisDer[f], NULL);
377: PetscFEGetFaceTabulation(fe, &prob->basisFace[f], &prob->basisDerFace[f], NULL);
378: } else if (id == PETSCFV_CLASSID) {
379: PetscFV fv = (PetscFV) obj;
381: PetscFVGetQuadrature(fv, &q);
382: PetscFVGetNumComponents(fv, &Nc);
383: Nb = Nc;
384: PetscFVGetDefaultTabulation(fv, &prob->basis[f], &prob->basisDer[f], NULL);
385: /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */
386: } else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", f);
387: prob->Nc[f] = Nc;
388: prob->Nb[f] = Nb;
389: prob->off[f+1] = Nc + prob->off[f];
390: prob->offDer[f+1] = Nc*dim + prob->offDer[f];
391: if (q) {PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL);}
392: NqMax = PetscMax(NqMax, Nq);
393: NcMax = PetscMax(NcMax, Nc);
394: prob->totDim += Nb;
395: prob->totComp += Nc;
396: /* There are two faces for all fields but the cohesive field on a hybrid cell */
397: if (prob->isHybrid && (f < Nf-1)) prob->totDim += Nb;
398: }
399: work = PetscMax(prob->totComp*dim, PetscSqr(NcMax*dim));
400: /* Allocate works space */
401: if (prob->isHybrid) NsMax = 2;
402: PetscMalloc5(NsMax*prob->totComp,&prob->u,NsMax*prob->totComp,&prob->u_t,NsMax*prob->totComp*dimEmbed,&prob->u_x,dimEmbed,&prob->x,work,&prob->refSpaceDer);
403: PetscMalloc6(NsMax*NqMax*NcMax,&prob->f0,NsMax*NqMax*NcMax*dim,&prob->f1,
404: NsMax*NsMax*NqMax*NcMax*NcMax,&prob->g0,NsMax*NsMax*NqMax*NcMax*NcMax*dim,&prob->g1,
405: NsMax*NsMax*NqMax*NcMax*NcMax*dim,&prob->g2,NsMax*NsMax*NqMax*NcMax*NcMax*dim*dim,&prob->g3);
406: if (prob->ops->setup) {(*prob->ops->setup)(prob);}
407: prob->setup = PETSC_TRUE;
408: return(0);
409: }
411: static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob)
412: {
416: PetscFree2(prob->Nc,prob->Nb);
417: PetscFree2(prob->off,prob->offDer);
418: PetscFree4(prob->basis,prob->basisDer,prob->basisFace,prob->basisDerFace);
419: PetscFree5(prob->u,prob->u_t,prob->u_x,prob->x,prob->refSpaceDer);
420: PetscFree6(prob->f0,prob->f1,prob->g0,prob->g1,prob->g2,prob->g3);
421: return(0);
422: }
424: static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew)
425: {
426: PetscObject *tmpd;
427: PetscBool *tmpi;
428: PetscPointFunc *tmpobj, *tmpf, *tmpup;
429: PetscPointJac *tmpg, *tmpgp, *tmpgt;
430: PetscBdPointFunc *tmpfbd;
431: PetscBdPointJac *tmpgbd;
432: PetscRiemannFunc *tmpr;
433: PetscSimplePointFunc *tmpexactSol;
434: void **tmpctx;
435: PetscInt Nf = prob->Nf, f;
436: PetscErrorCode ierr;
439: if (Nf >= NfNew) return(0);
440: prob->setup = PETSC_FALSE;
441: PetscDSDestroyStructs_Static(prob);
442: PetscMalloc2(NfNew, &tmpd, NfNew, &tmpi);
443: for (f = 0; f < Nf; ++f) {tmpd[f] = prob->disc[f]; tmpi[f] = prob->implicit[f];}
444: for (f = Nf; f < NfNew; ++f) {tmpd[f] = NULL; tmpi[f] = PETSC_TRUE;}
445: PetscFree2(prob->disc, prob->implicit);
446: prob->Nf = NfNew;
447: prob->disc = tmpd;
448: prob->implicit = tmpi;
449: PetscCalloc7(NfNew, &tmpobj, NfNew*2, &tmpf, NfNew*NfNew*4, &tmpg, NfNew*NfNew*4, &tmpgp, NfNew*NfNew*4, &tmpgt, NfNew, &tmpr, NfNew, &tmpctx);
450: PetscCalloc1(NfNew, &tmpup);
451: for (f = 0; f < Nf; ++f) tmpobj[f] = prob->obj[f];
452: for (f = 0; f < Nf*2; ++f) tmpf[f] = prob->f[f];
453: for (f = 0; f < Nf*Nf*4; ++f) tmpg[f] = prob->g[f];
454: for (f = 0; f < Nf*Nf*4; ++f) tmpgp[f] = prob->gp[f];
455: for (f = 0; f < Nf; ++f) tmpr[f] = prob->r[f];
456: for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f];
457: for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f];
458: for (f = Nf; f < NfNew; ++f) tmpobj[f] = NULL;
459: for (f = Nf*2; f < NfNew*2; ++f) tmpf[f] = NULL;
460: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpg[f] = NULL;
461: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgp[f] = NULL;
462: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgt[f] = NULL;
463: for (f = Nf; f < NfNew; ++f) tmpr[f] = NULL;
464: for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL;
465: for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL;
466: PetscFree7(prob->obj, prob->f, prob->g, prob->gp, prob->gt, prob->r, prob->ctx);
467: PetscFree(prob->update);
468: prob->obj = tmpobj;
469: prob->f = tmpf;
470: prob->g = tmpg;
471: prob->gp = tmpgp;
472: prob->gt = tmpgt;
473: prob->r = tmpr;
474: prob->update = tmpup;
475: prob->ctx = tmpctx;
476: PetscCalloc3(NfNew*2, &tmpfbd, NfNew*NfNew*4, &tmpgbd, NfNew, &tmpexactSol);
477: for (f = 0; f < Nf*2; ++f) tmpfbd[f] = prob->fBd[f];
478: for (f = 0; f < Nf*Nf*4; ++f) tmpgbd[f] = prob->gBd[f];
479: for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f];
480: for (f = Nf*2; f < NfNew*2; ++f) tmpfbd[f] = NULL;
481: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgbd[f] = NULL;
482: for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL;
483: PetscFree3(prob->fBd, prob->gBd, prob->exactSol);
484: prob->fBd = tmpfbd;
485: prob->gBd = tmpgbd;
486: prob->exactSol = tmpexactSol;
487: return(0);
488: }
490: /*@
491: PetscDSDestroy - Destroys a PetscDS object
493: Collective on PetscDS
495: Input Parameter:
496: . prob - the PetscDS object to destroy
498: Level: developer
500: .seealso PetscDSView()
501: @*/
502: PetscErrorCode PetscDSDestroy(PetscDS *prob)
503: {
504: PetscInt f;
505: DSBoundary next;
509: if (!*prob) return(0);
512: if (--((PetscObject)(*prob))->refct > 0) {*prob = 0; return(0);}
513: ((PetscObject) (*prob))->refct = 0;
514: if ((*prob)->subprobs) {
515: PetscInt dim, d;
517: PetscDSGetSpatialDimension(*prob, &dim);
518: for (d = 0; d < dim; ++d) {PetscDSDestroy(&(*prob)->subprobs[d]);}
519: }
520: PetscFree((*prob)->subprobs);
521: PetscDSDestroyStructs_Static(*prob);
522: for (f = 0; f < (*prob)->Nf; ++f) {
523: PetscObjectDereference((*prob)->disc[f]);
524: }
525: PetscFree2((*prob)->disc, (*prob)->implicit);
526: PetscFree7((*prob)->obj,(*prob)->f,(*prob)->g,(*prob)->gp,(*prob)->gt,(*prob)->r,(*prob)->ctx);
527: PetscFree((*prob)->update);
528: PetscFree3((*prob)->fBd,(*prob)->gBd,(*prob)->exactSol);
529: if ((*prob)->ops->destroy) {(*(*prob)->ops->destroy)(*prob);}
530: next = (*prob)->boundary;
531: while (next) {
532: DSBoundary b = next;
534: next = b->next;
535: PetscFree(b->comps);
536: PetscFree(b->ids);
537: PetscFree(b->name);
538: PetscFree(b->labelname);
539: PetscFree(b);
540: }
541: PetscFree((*prob)->constants);
542: PetscHeaderDestroy(prob);
543: return(0);
544: }
546: /*@
547: PetscDSCreate - Creates an empty PetscDS object. The type can then be set with PetscDSSetType().
549: Collective on MPI_Comm
551: Input Parameter:
552: . comm - The communicator for the PetscDS object
554: Output Parameter:
555: . prob - The PetscDS object
557: Level: beginner
559: .seealso: PetscDSSetType(), PETSCDSBASIC
560: @*/
561: PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *prob)
562: {
563: PetscDS p;
568: *prob = NULL;
569: PetscDSInitializePackage();
571: PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView);
573: p->Nf = 0;
574: p->setup = PETSC_FALSE;
575: p->numConstants = 0;
576: p->constants = NULL;
577: p->dimEmbed = -1;
578: p->useJacPre = PETSC_TRUE;
580: *prob = p;
581: return(0);
582: }
584: /*@
585: PetscDSGetNumFields - Returns the number of fields in the DS
587: Not collective
589: Input Parameter:
590: . prob - The PetscDS object
592: Output Parameter:
593: . Nf - The number of fields
595: Level: beginner
597: .seealso: PetscDSGetSpatialDimension(), PetscDSCreate()
598: @*/
599: PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf)
600: {
604: *Nf = prob->Nf;
605: return(0);
606: }
608: /*@
609: PetscDSGetSpatialDimension - Returns the spatial dimension of the DS, meaning the topological dimension of the discretizations
611: Not collective
613: Input Parameter:
614: . prob - The PetscDS object
616: Output Parameter:
617: . dim - The spatial dimension
619: Level: beginner
621: .seealso: PetscDSGetCoordinateDimension(), PetscDSGetNumFields(), PetscDSCreate()
622: @*/
623: PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim)
624: {
630: *dim = 0;
631: if (prob->Nf) {
632: PetscObject obj;
633: PetscClassId id;
635: PetscDSGetDiscretization(prob, 0, &obj);
636: PetscObjectGetClassId(obj, &id);
637: if (id == PETSCFE_CLASSID) {PetscFEGetSpatialDimension((PetscFE) obj, dim);}
638: else if (id == PETSCFV_CLASSID) {PetscFVGetSpatialDimension((PetscFV) obj, dim);}
639: else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
640: }
641: return(0);
642: }
644: /*@
645: PetscDSGetCoordinateDimension - Returns the coordinate dimension of the DS, meaning the dimension of the space into which the discretiaztions are embedded
647: Not collective
649: Input Parameter:
650: . prob - The PetscDS object
652: Output Parameter:
653: . dimEmbed - The coordinate dimension
655: Level: beginner
657: .seealso: PetscDSSetCoordinateDimension(), PetscDSGetSpatialDimension(), PetscDSGetNumFields(), PetscDSCreate()
658: @*/
659: PetscErrorCode PetscDSGetCoordinateDimension(PetscDS prob, PetscInt *dimEmbed)
660: {
664: if (prob->dimEmbed < 0) SETERRQ(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONGSTATE, "No coordinate dimension set for this DS");
665: *dimEmbed = prob->dimEmbed;
666: return(0);
667: }
669: /*@
670: PetscDSSetCoordinateDimension - Set the coordinate dimension of the DS, meaning the dimension of the space into which the discretiaztions are embedded
672: Logically collective on DS
674: Input Parameters:
675: + prob - The PetscDS object
676: - dimEmbed - The coordinate dimension
678: Level: beginner
680: .seealso: PetscDSGetCoordinateDimension(), PetscDSGetSpatialDimension(), PetscDSGetNumFields(), PetscDSCreate()
681: @*/
682: PetscErrorCode PetscDSSetCoordinateDimension(PetscDS prob, PetscInt dimEmbed)
683: {
686: if (dimEmbed < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Coordinate dimension must be non-negative, not %D", dimEmbed);
687: prob->dimEmbed = dimEmbed;
688: return(0);
689: }
691: /*@
692: PetscDSGetHybrid - Returns the flag for a hybrid (cohesive) cell
694: Not collective
696: Input Parameter:
697: . prob - The PetscDS object
699: Output Parameter:
700: . isHybrid - The flag
702: Level: developer
704: .seealso: PetscDSSetHybrid(), PetscDSCreate()
705: @*/
706: PetscErrorCode PetscDSGetHybrid(PetscDS prob, PetscBool *isHybrid)
707: {
711: *isHybrid = prob->isHybrid;
712: return(0);
713: }
715: /*@
716: PetscDSSetHybrid - Set the flag for a hybrid (cohesive) cell
718: Not collective
720: Input Parameters:
721: + prob - The PetscDS object
722: - isHybrid - The flag
724: Level: developer
726: .seealso: PetscDSGetHybrid(), PetscDSCreate()
727: @*/
728: PetscErrorCode PetscDSSetHybrid(PetscDS prob, PetscBool isHybrid)
729: {
732: prob->isHybrid = isHybrid;
733: return(0);
734: }
736: /*@
737: PetscDSGetTotalDimension - Returns the total size of the approximation space for this system
739: Not collective
741: Input Parameter:
742: . prob - The PetscDS object
744: Output Parameter:
745: . dim - The total problem dimension
747: Level: beginner
749: .seealso: PetscDSGetNumFields(), PetscDSCreate()
750: @*/
751: PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim)
752: {
757: PetscDSSetUp(prob);
759: *dim = prob->totDim;
760: return(0);
761: }
763: /*@
764: PetscDSGetTotalComponents - Returns the total number of components in this system
766: Not collective
768: Input Parameter:
769: . prob - The PetscDS object
771: Output Parameter:
772: . dim - The total number of components
774: Level: beginner
776: .seealso: PetscDSGetNumFields(), PetscDSCreate()
777: @*/
778: PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc)
779: {
784: PetscDSSetUp(prob);
786: *Nc = prob->totComp;
787: return(0);
788: }
790: /*@
791: PetscDSGetDiscretization - Returns the discretization object for the given field
793: Not collective
795: Input Parameters:
796: + prob - The PetscDS object
797: - f - The field number
799: Output Parameter:
800: . disc - The discretization object
802: Level: beginner
804: .seealso: PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
805: @*/
806: PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc)
807: {
811: 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);
812: *disc = prob->disc[f];
813: return(0);
814: }
816: /*@
817: PetscDSSetDiscretization - Sets the discretization object for the given field
819: Not collective
821: Input Parameters:
822: + prob - The PetscDS object
823: . f - The field number
824: - disc - The discretization object
826: Level: beginner
828: .seealso: PetscDSGetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
829: @*/
830: PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc)
831: {
837: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
838: PetscDSEnlarge_Static(prob, f+1);
839: PetscObjectDereference(prob->disc[f]);
840: prob->disc[f] = disc;
841: PetscObjectReference(disc);
842: {
843: PetscClassId id;
845: PetscObjectGetClassId(disc, &id);
846: if (id == PETSCFE_CLASSID) {
847: PetscDSSetImplicit(prob, f, PETSC_TRUE);
848: } else if (id == PETSCFV_CLASSID) {
849: PetscDSSetImplicit(prob, f, PETSC_FALSE);
850: }
851: }
852: return(0);
853: }
855: /*@
856: PetscDSAddDiscretization - Adds a discretization object
858: Not collective
860: Input Parameters:
861: + prob - The PetscDS object
862: - disc - The boundary discretization object
864: Level: beginner
866: .seealso: PetscDSGetDiscretization(), PetscDSSetDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
867: @*/
868: PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc)
869: {
873: PetscDSSetDiscretization(prob, prob->Nf, disc);
874: return(0);
875: }
877: /*@
878: PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for IMEX
880: Not collective
882: Input Parameters:
883: + prob - The PetscDS object
884: - f - The field number
886: Output Parameter:
887: . implicit - The flag indicating what kind of solve to use for this field
889: Level: developer
891: .seealso: PetscDSSetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
892: @*/
893: PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit)
894: {
898: 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);
899: *implicit = prob->implicit[f];
900: return(0);
901: }
903: /*@
904: PetscDSSetImplicit - Set 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
911: - implicit - The flag indicating what kind of solve to use for this field
913: Level: developer
915: .seealso: PetscDSGetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
916: @*/
917: PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit)
918: {
921: 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);
922: prob->implicit[f] = implicit;
923: return(0);
924: }
926: PetscErrorCode PetscDSGetObjective(PetscDS prob, PetscInt f,
927: void (**obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
928: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
929: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
930: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
931: {
935: 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);
936: *obj = prob->obj[f];
937: return(0);
938: }
940: PetscErrorCode PetscDSSetObjective(PetscDS prob, PetscInt f,
941: void (*obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
942: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
943: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
944: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
945: {
951: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
952: PetscDSEnlarge_Static(prob, f+1);
953: prob->obj[f] = obj;
954: return(0);
955: }
957: /*@C
958: PetscDSGetResidual - Get the pointwise residual function for a given test field
960: Not collective
962: Input Parameters:
963: + prob - The PetscDS
964: - f - The test field number
966: Output Parameters:
967: + f0 - integrand for the test function term
968: - f1 - integrand for the test function gradient term
970: Note: We are using a first order FEM model for the weak form:
972: \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)
974: The calling sequence for the callbacks f0 and f1 is given by:
976: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
977: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
978: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
979: $ PetscReal t, const PetscReal x[], PetscScalar f0[])
981: + dim - the spatial dimension
982: . Nf - the number of fields
983: . uOff - the offset into u[] and u_t[] for each field
984: . uOff_x - the offset into u_x[] for each field
985: . u - each field evaluated at the current point
986: . u_t - the time derivative of each field evaluated at the current point
987: . u_x - the gradient of each field evaluated at the current point
988: . aOff - the offset into a[] and a_t[] for each auxiliary field
989: . aOff_x - the offset into a_x[] for each auxiliary field
990: . a - each auxiliary field evaluated at the current point
991: . a_t - the time derivative of each auxiliary field evaluated at the current point
992: . a_x - the gradient of auxiliary each field evaluated at the current point
993: . t - current time
994: . x - coordinates of the current point
995: . numConstants - number of constant parameters
996: . constants - constant parameters
997: - f0 - output values at the current point
999: Level: intermediate
1001: .seealso: PetscDSSetResidual()
1002: @*/
1003: PetscErrorCode PetscDSGetResidual(PetscDS prob, PetscInt f,
1004: void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1005: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1006: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1007: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1008: void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1009: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1010: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1011: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1012: {
1015: 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);
1018: return(0);
1019: }
1021: /*@C
1022: PetscDSSetResidual - Set the pointwise residual function for a given test field
1024: Not collective
1026: Input Parameters:
1027: + prob - The PetscDS
1028: . f - The test field number
1029: . f0 - integrand for the test function term
1030: - f1 - integrand for the test function gradient term
1032: Note: We are using a first order FEM model for the weak form:
1034: \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)
1036: The calling sequence for the callbacks f0 and f1 is given by:
1038: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1039: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1040: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1041: $ PetscReal t, const PetscReal x[], PetscScalar f0[])
1043: + dim - the spatial dimension
1044: . Nf - the number of fields
1045: . uOff - the offset into u[] and u_t[] for each field
1046: . uOff_x - the offset into u_x[] for each field
1047: . u - each field evaluated at the current point
1048: . u_t - the time derivative of each field evaluated at the current point
1049: . u_x - the gradient of each field evaluated at the current point
1050: . aOff - the offset into a[] and a_t[] for each auxiliary field
1051: . aOff_x - the offset into a_x[] for each auxiliary field
1052: . a - each auxiliary field evaluated at the current point
1053: . a_t - the time derivative of each auxiliary field evaluated at the current point
1054: . a_x - the gradient of auxiliary each field evaluated at the current point
1055: . t - current time
1056: . x - coordinates of the current point
1057: . numConstants - number of constant parameters
1058: . constants - constant parameters
1059: - f0 - output values at the current point
1061: Level: intermediate
1063: .seealso: PetscDSGetResidual()
1064: @*/
1065: PetscErrorCode PetscDSSetResidual(PetscDS prob, PetscInt f,
1066: void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1067: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1068: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1069: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1070: void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1071: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1072: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1073: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1074: {
1081: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1082: PetscDSEnlarge_Static(prob, f+1);
1083: prob->f[f*2+0] = f0;
1084: prob->f[f*2+1] = f1;
1085: return(0);
1086: }
1088: /*@C
1089: PetscDSHasJacobian - Signals that Jacobian functions have been set
1091: Not collective
1093: Input Parameter:
1094: . prob - The PetscDS
1096: Output Parameter:
1097: . hasJac - flag that pointwise function for the Jacobian has been set
1099: Level: intermediate
1101: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1102: @*/
1103: PetscErrorCode PetscDSHasJacobian(PetscDS prob, PetscBool *hasJac)
1104: {
1105: PetscInt f, g, h;
1109: *hasJac = PETSC_FALSE;
1110: for (f = 0; f < prob->Nf; ++f) {
1111: for (g = 0; g < prob->Nf; ++g) {
1112: for (h = 0; h < 4; ++h) {
1113: if (prob->g[(f*prob->Nf + g)*4+h]) *hasJac = PETSC_TRUE;
1114: }
1115: }
1116: }
1117: return(0);
1118: }
1120: /*@C
1121: PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field
1123: Not collective
1125: Input Parameters:
1126: + prob - The PetscDS
1127: . f - The test field number
1128: - g - The field number
1130: Output Parameters:
1131: + g0 - integrand for the test and basis function term
1132: . g1 - integrand for the test function and basis function gradient term
1133: . g2 - integrand for the test function gradient and basis function term
1134: - g3 - integrand for the test function gradient and basis function gradient term
1136: Note: We are using a first order FEM model for the weak form:
1138: \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
1140: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1142: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1143: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1144: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1145: $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1147: + dim - the spatial dimension
1148: . Nf - the number of fields
1149: . uOff - the offset into u[] and u_t[] for each field
1150: . uOff_x - the offset into u_x[] for each field
1151: . u - each field evaluated at the current point
1152: . u_t - the time derivative of each field evaluated at the current point
1153: . u_x - the gradient of each field evaluated at the current point
1154: . aOff - the offset into a[] and a_t[] for each auxiliary field
1155: . aOff_x - the offset into a_x[] for each auxiliary field
1156: . a - each auxiliary field evaluated at the current point
1157: . a_t - the time derivative of each auxiliary field evaluated at the current point
1158: . a_x - the gradient of auxiliary each field evaluated at the current point
1159: . t - current time
1160: . u_tShift - the multiplier a for dF/dU_t
1161: . x - coordinates of the current point
1162: . numConstants - number of constant parameters
1163: . constants - constant parameters
1164: - g0 - output values at the current point
1166: Level: intermediate
1168: .seealso: PetscDSSetJacobian()
1169: @*/
1170: PetscErrorCode PetscDSGetJacobian(PetscDS prob, PetscInt f, PetscInt g,
1171: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1172: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1173: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1174: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1175: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1176: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1177: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1178: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1179: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1180: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1181: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1182: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1183: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1184: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1185: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1186: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1187: {
1190: 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);
1191: 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);
1196: return(0);
1197: }
1199: /*@C
1200: PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields
1202: Not collective
1204: Input Parameters:
1205: + prob - The PetscDS
1206: . f - The test field number
1207: . g - The field number
1208: . g0 - integrand for the test and basis function term
1209: . g1 - integrand for the test function and basis function gradient term
1210: . g2 - integrand for the test function gradient and basis function term
1211: - g3 - integrand for the test function gradient and basis function gradient term
1213: Note: We are using a first order FEM model for the weak form:
1215: \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
1217: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1219: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1220: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1221: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1222: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1224: + dim - the spatial dimension
1225: . Nf - the number of fields
1226: . uOff - the offset into u[] and u_t[] for each field
1227: . uOff_x - the offset into u_x[] for each field
1228: . u - each field evaluated at the current point
1229: . u_t - the time derivative of each field evaluated at the current point
1230: . u_x - the gradient of each field evaluated at the current point
1231: . aOff - the offset into a[] and a_t[] for each auxiliary field
1232: . aOff_x - the offset into a_x[] for each auxiliary field
1233: . a - each auxiliary field evaluated at the current point
1234: . a_t - the time derivative of each auxiliary field evaluated at the current point
1235: . a_x - the gradient of auxiliary each field evaluated at the current point
1236: . t - current time
1237: . u_tShift - the multiplier a for dF/dU_t
1238: . x - coordinates of the current point
1239: . numConstants - number of constant parameters
1240: . constants - constant parameters
1241: - g0 - output values at the current point
1243: Level: intermediate
1245: .seealso: PetscDSGetJacobian()
1246: @*/
1247: PetscErrorCode PetscDSSetJacobian(PetscDS prob, PetscInt f, PetscInt g,
1248: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1249: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1250: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1251: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1252: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1253: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1254: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1255: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1256: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1257: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1258: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1259: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1260: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1261: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1262: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1263: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1264: {
1273: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1274: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1275: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1276: prob->g[(f*prob->Nf + g)*4+0] = g0;
1277: prob->g[(f*prob->Nf + g)*4+1] = g1;
1278: prob->g[(f*prob->Nf + g)*4+2] = g2;
1279: prob->g[(f*prob->Nf + g)*4+3] = g3;
1280: return(0);
1281: }
1283: /*@C
1284: PetscDSUseJacobianPreconditioner - Whether to construct a Jacobian preconditioner
1286: Not collective
1288: Input Parameters:
1289: + prob - The PetscDS
1290: - useJacPre - flag that enables construction of a Jacobian preconditioner
1292: Level: intermediate
1294: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1295: @*/
1296: PetscErrorCode PetscDSUseJacobianPreconditioner(PetscDS prob, PetscBool useJacPre)
1297: {
1300: prob->useJacPre = useJacPre;
1301: return(0);
1302: }
1304: /*@C
1305: PetscDSHasJacobianPreconditioner - Signals that a Jacobian preconditioner matrix has been set
1307: Not collective
1309: Input Parameter:
1310: . prob - The PetscDS
1312: Output Parameter:
1313: . hasJacPre - flag that pointwise function for Jacobian preconditioner matrix has been set
1315: Level: intermediate
1317: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1318: @*/
1319: PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS prob, PetscBool *hasJacPre)
1320: {
1321: PetscInt f, g, h;
1325: *hasJacPre = PETSC_FALSE;
1326: if (!prob->useJacPre) return(0);
1327: for (f = 0; f < prob->Nf; ++f) {
1328: for (g = 0; g < prob->Nf; ++g) {
1329: for (h = 0; h < 4; ++h) {
1330: if (prob->gp[(f*prob->Nf + g)*4+h]) *hasJacPre = PETSC_TRUE;
1331: }
1332: }
1333: }
1334: return(0);
1335: }
1337: /*@C
1338: 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.
1340: Not collective
1342: Input Parameters:
1343: + prob - The PetscDS
1344: . f - The test field number
1345: - g - The field number
1347: Output Parameters:
1348: + g0 - integrand for the test and basis function term
1349: . g1 - integrand for the test function and basis function gradient term
1350: . g2 - integrand for the test function gradient and basis function term
1351: - g3 - integrand for the test function gradient and basis function gradient term
1353: Note: We are using a first order FEM model for the weak form:
1355: \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
1357: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1359: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1360: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1361: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1362: $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1364: + dim - the spatial dimension
1365: . Nf - the number of fields
1366: . uOff - the offset into u[] and u_t[] for each field
1367: . uOff_x - the offset into u_x[] for each field
1368: . u - each field evaluated at the current point
1369: . u_t - the time derivative of each field evaluated at the current point
1370: . u_x - the gradient of each field evaluated at the current point
1371: . aOff - the offset into a[] and a_t[] for each auxiliary field
1372: . aOff_x - the offset into a_x[] for each auxiliary field
1373: . a - each auxiliary field evaluated at the current point
1374: . a_t - the time derivative of each auxiliary field evaluated at the current point
1375: . a_x - the gradient of auxiliary each field evaluated at the current point
1376: . t - current time
1377: . u_tShift - the multiplier a for dF/dU_t
1378: . x - coordinates of the current point
1379: . numConstants - number of constant parameters
1380: . constants - constant parameters
1381: - g0 - output values at the current point
1383: Level: intermediate
1385: .seealso: PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1386: @*/
1387: PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g,
1388: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1389: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1390: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1391: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1392: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1393: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1394: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1395: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1396: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1397: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1398: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1399: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1400: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1401: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1402: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1403: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1404: {
1407: 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);
1408: 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);
1413: return(0);
1414: }
1416: /*@C
1417: 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.
1419: Not collective
1421: Input Parameters:
1422: + prob - The PetscDS
1423: . f - The test field number
1424: . g - The field number
1425: . g0 - integrand for the test and basis function term
1426: . g1 - integrand for the test function and basis function gradient term
1427: . g2 - integrand for the test function gradient and basis function term
1428: - g3 - integrand for the test function gradient and basis function gradient term
1430: Note: We are using a first order FEM model for the weak form:
1432: \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
1434: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1436: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1437: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1438: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1439: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1441: + dim - the spatial dimension
1442: . Nf - the number of fields
1443: . uOff - the offset into u[] and u_t[] for each field
1444: . uOff_x - the offset into u_x[] for each field
1445: . u - each field evaluated at the current point
1446: . u_t - the time derivative of each field evaluated at the current point
1447: . u_x - the gradient of each field evaluated at the current point
1448: . aOff - the offset into a[] and a_t[] for each auxiliary field
1449: . aOff_x - the offset into a_x[] for each auxiliary field
1450: . a - each auxiliary field evaluated at the current point
1451: . a_t - the time derivative of each auxiliary field evaluated at the current point
1452: . a_x - the gradient of auxiliary each field evaluated at the current point
1453: . t - current time
1454: . u_tShift - the multiplier a for dF/dU_t
1455: . x - coordinates of the current point
1456: . numConstants - number of constant parameters
1457: . constants - constant parameters
1458: - g0 - output values at the current point
1460: Level: intermediate
1462: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobian()
1463: @*/
1464: PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g,
1465: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1466: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1467: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1468: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1469: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1470: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1471: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1472: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1473: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1474: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1475: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1476: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1477: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1478: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1479: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1480: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1481: {
1490: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1491: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1492: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1493: prob->gp[(f*prob->Nf + g)*4+0] = g0;
1494: prob->gp[(f*prob->Nf + g)*4+1] = g1;
1495: prob->gp[(f*prob->Nf + g)*4+2] = g2;
1496: prob->gp[(f*prob->Nf + g)*4+3] = g3;
1497: return(0);
1498: }
1500: /*@C
1501: PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, dF/du_t, has been set
1503: Not collective
1505: Input Parameter:
1506: . prob - The PetscDS
1508: Output Parameter:
1509: . hasDynJac - flag that pointwise function for dynamic Jacobian has been set
1511: Level: intermediate
1513: .seealso: PetscDSGetDynamicJacobian(), PetscDSSetDynamicJacobian(), PetscDSGetJacobian()
1514: @*/
1515: PetscErrorCode PetscDSHasDynamicJacobian(PetscDS prob, PetscBool *hasDynJac)
1516: {
1517: PetscInt f, g, h;
1521: *hasDynJac = PETSC_FALSE;
1522: for (f = 0; f < prob->Nf; ++f) {
1523: for (g = 0; g < prob->Nf; ++g) {
1524: for (h = 0; h < 4; ++h) {
1525: if (prob->gt[(f*prob->Nf + g)*4+h]) *hasDynJac = PETSC_TRUE;
1526: }
1527: }
1528: }
1529: return(0);
1530: }
1532: /*@C
1533: PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, dF/du_t, function for given test and basis field
1535: Not collective
1537: Input Parameters:
1538: + prob - The PetscDS
1539: . f - The test field number
1540: - g - The field number
1542: Output Parameters:
1543: + g0 - integrand for the test and basis function term
1544: . g1 - integrand for the test function and basis function gradient term
1545: . g2 - integrand for the test function gradient and basis function term
1546: - g3 - integrand for the test function gradient and basis function gradient term
1548: Note: We are using a first order FEM model for the weak form:
1550: \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
1552: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1554: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1555: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1556: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1557: $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1559: + dim - the spatial dimension
1560: . Nf - the number of fields
1561: . uOff - the offset into u[] and u_t[] for each field
1562: . uOff_x - the offset into u_x[] for each field
1563: . u - each field evaluated at the current point
1564: . u_t - the time derivative of each field evaluated at the current point
1565: . u_x - the gradient of each field evaluated at the current point
1566: . aOff - the offset into a[] and a_t[] for each auxiliary field
1567: . aOff_x - the offset into a_x[] for each auxiliary field
1568: . a - each auxiliary field evaluated at the current point
1569: . a_t - the time derivative of each auxiliary field evaluated at the current point
1570: . a_x - the gradient of auxiliary each field evaluated at the current point
1571: . t - current time
1572: . u_tShift - the multiplier a for dF/dU_t
1573: . x - coordinates of the current point
1574: . numConstants - number of constant parameters
1575: . constants - constant parameters
1576: - g0 - output values at the current point
1578: Level: intermediate
1580: .seealso: PetscDSSetJacobian()
1581: @*/
1582: PetscErrorCode PetscDSGetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g,
1583: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1584: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1585: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1586: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1587: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1588: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1589: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1590: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1591: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1592: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1593: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1594: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1595: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1596: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1597: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1598: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1599: {
1602: 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);
1603: 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);
1608: return(0);
1609: }
1611: /*@C
1612: PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, dF/du_t, function for given test and basis fields
1614: Not collective
1616: Input Parameters:
1617: + prob - The PetscDS
1618: . f - The test field number
1619: . g - The field number
1620: . g0 - integrand for the test and basis function term
1621: . g1 - integrand for the test function and basis function gradient term
1622: . g2 - integrand for the test function gradient and basis function term
1623: - g3 - integrand for the test function gradient and basis function gradient term
1625: Note: We are using a first order FEM model for the weak form:
1627: \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
1629: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1631: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1632: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1633: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1634: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1636: + dim - the spatial dimension
1637: . Nf - the number of fields
1638: . uOff - the offset into u[] and u_t[] for each field
1639: . uOff_x - the offset into u_x[] for each field
1640: . u - each field evaluated at the current point
1641: . u_t - the time derivative of each field evaluated at the current point
1642: . u_x - the gradient of each field evaluated at the current point
1643: . aOff - the offset into a[] and a_t[] for each auxiliary field
1644: . aOff_x - the offset into a_x[] for each auxiliary field
1645: . a - each auxiliary field evaluated at the current point
1646: . a_t - the time derivative of each auxiliary field evaluated at the current point
1647: . a_x - the gradient of auxiliary each field evaluated at the current point
1648: . t - current time
1649: . u_tShift - the multiplier a for dF/dU_t
1650: . x - coordinates of the current point
1651: . numConstants - number of constant parameters
1652: . constants - constant parameters
1653: - g0 - output values at the current point
1655: Level: intermediate
1657: .seealso: PetscDSGetJacobian()
1658: @*/
1659: PetscErrorCode PetscDSSetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g,
1660: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1661: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1662: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1663: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1664: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1665: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1666: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1667: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1668: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1669: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1670: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1671: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1672: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1673: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1674: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1675: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1676: {
1685: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1686: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1687: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1688: prob->gt[(f*prob->Nf + g)*4+0] = g0;
1689: prob->gt[(f*prob->Nf + g)*4+1] = g1;
1690: prob->gt[(f*prob->Nf + g)*4+2] = g2;
1691: prob->gt[(f*prob->Nf + g)*4+3] = g3;
1692: return(0);
1693: }
1695: /*@C
1696: PetscDSGetRiemannSolver - Returns the Riemann solver for the given field
1698: Not collective
1700: Input Arguments:
1701: + prob - The PetscDS object
1702: - f - The field number
1704: Output Argument:
1705: . r - Riemann solver
1707: Calling sequence for r:
1709: $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)
1711: + dim - The spatial dimension
1712: . Nf - The number of fields
1713: . x - The coordinates at a point on the interface
1714: . n - The normal vector to the interface
1715: . uL - The state vector to the left of the interface
1716: . uR - The state vector to the right of the interface
1717: . flux - output array of flux through the interface
1718: . numConstants - number of constant parameters
1719: . constants - constant parameters
1720: - ctx - optional user context
1722: Level: intermediate
1724: .seealso: PetscDSSetRiemannSolver()
1725: @*/
1726: PetscErrorCode PetscDSGetRiemannSolver(PetscDS prob, PetscInt f,
1727: 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))
1728: {
1731: 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);
1733: *r = prob->r[f];
1734: return(0);
1735: }
1737: /*@C
1738: PetscDSSetRiemannSolver - Sets the Riemann solver for the given field
1740: Not collective
1742: Input Arguments:
1743: + prob - The PetscDS object
1744: . f - The field number
1745: - r - Riemann solver
1747: Calling sequence for r:
1749: $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)
1751: + dim - The spatial dimension
1752: . Nf - The number of fields
1753: . x - The coordinates at a point on the interface
1754: . n - The normal vector to the interface
1755: . uL - The state vector to the left of the interface
1756: . uR - The state vector to the right of the interface
1757: . flux - output array of flux through the interface
1758: . numConstants - number of constant parameters
1759: . constants - constant parameters
1760: - ctx - optional user context
1762: Level: intermediate
1764: .seealso: PetscDSGetRiemannSolver()
1765: @*/
1766: PetscErrorCode PetscDSSetRiemannSolver(PetscDS prob, PetscInt f,
1767: 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))
1768: {
1774: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1775: PetscDSEnlarge_Static(prob, f+1);
1776: prob->r[f] = r;
1777: return(0);
1778: }
1780: /*@C
1781: PetscDSGetUpdate - Get the pointwise update function for a given field
1783: Not collective
1785: Input Parameters:
1786: + prob - The PetscDS
1787: - f - The field number
1789: Output Parameters:
1790: + update - update function
1792: Note: The calling sequence for the callback update is given by:
1794: $ update(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1795: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1796: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1797: $ PetscReal t, const PetscReal x[], PetscScalar uNew[])
1799: + dim - the spatial dimension
1800: . Nf - the number of fields
1801: . uOff - the offset into u[] and u_t[] for each field
1802: . uOff_x - the offset into u_x[] for each field
1803: . u - each field evaluated at the current point
1804: . u_t - the time derivative of each field evaluated at the current point
1805: . u_x - the gradient of each field evaluated at the current point
1806: . aOff - the offset into a[] and a_t[] for each auxiliary field
1807: . aOff_x - the offset into a_x[] for each auxiliary field
1808: . a - each auxiliary field evaluated at the current point
1809: . a_t - the time derivative of each auxiliary field evaluated at the current point
1810: . a_x - the gradient of auxiliary each field evaluated at the current point
1811: . t - current time
1812: . x - coordinates of the current point
1813: - uNew - new value for field at the current point
1815: Level: intermediate
1817: .seealso: PetscDSSetUpdate(), PetscDSSetResidual()
1818: @*/
1819: PetscErrorCode PetscDSGetUpdate(PetscDS prob, PetscInt f,
1820: void (**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[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
1824: {
1827: 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);
1829: return(0);
1830: }
1832: /*@C
1833: PetscDSSetUpdate - Set the pointwise update function for a given field
1835: Not collective
1837: Input Parameters:
1838: + prob - The PetscDS
1839: . f - The field number
1840: - update - update function
1842: Note: The calling sequence for the callback update is given by:
1844: $ update(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1845: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1846: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1847: $ PetscReal t, const PetscReal x[], PetscScalar uNew[])
1849: + dim - the spatial dimension
1850: . Nf - the number of fields
1851: . uOff - the offset into u[] and u_t[] for each field
1852: . uOff_x - the offset into u_x[] for each field
1853: . u - each field evaluated at the current point
1854: . u_t - the time derivative of each field evaluated at the current point
1855: . u_x - the gradient of each field evaluated at the current point
1856: . aOff - the offset into a[] and a_t[] for each auxiliary field
1857: . aOff_x - the offset into a_x[] for each auxiliary field
1858: . a - each auxiliary field evaluated at the current point
1859: . a_t - the time derivative of each auxiliary field evaluated at the current point
1860: . a_x - the gradient of auxiliary each field evaluated at the current point
1861: . t - current time
1862: . x - coordinates of the current point
1863: - uNew - new field values at the current point
1865: Level: intermediate
1867: .seealso: PetscDSGetResidual()
1868: @*/
1869: PetscErrorCode PetscDSSetUpdate(PetscDS prob, PetscInt f,
1870: void (*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[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
1874: {
1880: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1881: PetscDSEnlarge_Static(prob, f+1);
1882: prob->update[f] = update;
1883: return(0);
1884: }
1886: PetscErrorCode PetscDSGetContext(PetscDS prob, PetscInt f, void **ctx)
1887: {
1890: 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);
1892: *ctx = prob->ctx[f];
1893: return(0);
1894: }
1896: PetscErrorCode PetscDSSetContext(PetscDS prob, PetscInt f, void *ctx)
1897: {
1902: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1903: PetscDSEnlarge_Static(prob, f+1);
1904: prob->ctx[f] = ctx;
1905: return(0);
1906: }
1908: /*@C
1909: PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field
1911: Not collective
1913: Input Parameters:
1914: + prob - The PetscDS
1915: - f - The test field number
1917: Output Parameters:
1918: + f0 - boundary integrand for the test function term
1919: - f1 - boundary integrand for the test function gradient term
1921: Note: We are using a first order FEM model for the weak form:
1923: \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
1925: The calling sequence for the callbacks f0 and f1 is given by:
1927: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1928: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1929: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1930: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[])
1932: + dim - the spatial dimension
1933: . Nf - the number of fields
1934: . uOff - the offset into u[] and u_t[] for each field
1935: . uOff_x - the offset into u_x[] for each field
1936: . u - each field evaluated at the current point
1937: . u_t - the time derivative of each field evaluated at the current point
1938: . u_x - the gradient of each field evaluated at the current point
1939: . aOff - the offset into a[] and a_t[] for each auxiliary field
1940: . aOff_x - the offset into a_x[] for each auxiliary field
1941: . a - each auxiliary field evaluated at the current point
1942: . a_t - the time derivative of each auxiliary field evaluated at the current point
1943: . a_x - the gradient of auxiliary each field evaluated at the current point
1944: . t - current time
1945: . x - coordinates of the current point
1946: . n - unit normal at the current point
1947: . numConstants - number of constant parameters
1948: . constants - constant parameters
1949: - f0 - output values at the current point
1951: Level: intermediate
1953: .seealso: PetscDSSetBdResidual()
1954: @*/
1955: PetscErrorCode PetscDSGetBdResidual(PetscDS prob, PetscInt f,
1956: void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1957: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1958: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1959: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1960: void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1961: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1962: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1963: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1964: {
1967: 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);
1970: return(0);
1971: }
1973: /*@C
1974: PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field
1976: Not collective
1978: Input Parameters:
1979: + prob - The PetscDS
1980: . f - The test field number
1981: . f0 - boundary integrand for the test function term
1982: - f1 - boundary integrand for the test function gradient term
1984: Note: We are using a first order FEM model for the weak form:
1986: \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
1988: The calling sequence for the callbacks f0 and f1 is given by:
1990: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1991: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1992: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1993: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[])
1995: + dim - the spatial dimension
1996: . Nf - the number of fields
1997: . uOff - the offset into u[] and u_t[] for each field
1998: . uOff_x - the offset into u_x[] for each field
1999: . u - each field evaluated at the current point
2000: . u_t - the time derivative of each field evaluated at the current point
2001: . u_x - the gradient of each field evaluated at the current point
2002: . aOff - the offset into a[] and a_t[] for each auxiliary field
2003: . aOff_x - the offset into a_x[] for each auxiliary field
2004: . a - each auxiliary field evaluated at the current point
2005: . a_t - the time derivative of each auxiliary field evaluated at the current point
2006: . a_x - the gradient of auxiliary each field evaluated at the current point
2007: . t - current time
2008: . x - coordinates of the current point
2009: . n - unit normal at the current point
2010: . numConstants - number of constant parameters
2011: . constants - constant parameters
2012: - f0 - output values at the current point
2014: Level: intermediate
2016: .seealso: PetscDSGetBdResidual()
2017: @*/
2018: PetscErrorCode PetscDSSetBdResidual(PetscDS prob, PetscInt f,
2019: void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2020: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2021: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2022: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
2023: void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2024: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2025: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2026: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
2027: {
2032: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2033: PetscDSEnlarge_Static(prob, f+1);
2036: return(0);
2037: }
2039: /*@C
2040: PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field
2042: Not collective
2044: Input Parameters:
2045: + prob - The PetscDS
2046: . f - The test field number
2047: - g - The field number
2049: Output Parameters:
2050: + g0 - integrand for the test and basis function term
2051: . g1 - integrand for the test function and basis function gradient term
2052: . g2 - integrand for the test function gradient and basis function term
2053: - g3 - integrand for the test function gradient and basis function gradient term
2055: Note: We are using a first order FEM model for the weak form:
2057: \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
2059: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
2061: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2062: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2063: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2064: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[])
2066: + dim - the spatial dimension
2067: . Nf - the number of fields
2068: . uOff - the offset into u[] and u_t[] for each field
2069: . uOff_x - the offset into u_x[] for each field
2070: . u - each field evaluated at the current point
2071: . u_t - the time derivative of each field evaluated at the current point
2072: . u_x - the gradient of each field evaluated at the current point
2073: . aOff - the offset into a[] and a_t[] for each auxiliary field
2074: . aOff_x - the offset into a_x[] for each auxiliary field
2075: . a - each auxiliary field evaluated at the current point
2076: . a_t - the time derivative of each auxiliary field evaluated at the current point
2077: . a_x - the gradient of auxiliary each field evaluated at the current point
2078: . t - current time
2079: . u_tShift - the multiplier a for dF/dU_t
2080: . x - coordinates of the current point
2081: . n - normal at the current point
2082: . numConstants - number of constant parameters
2083: . constants - constant parameters
2084: - g0 - output values at the current point
2086: Level: intermediate
2088: .seealso: PetscDSSetBdJacobian()
2089: @*/
2090: PetscErrorCode PetscDSGetBdJacobian(PetscDS prob, PetscInt f, PetscInt g,
2091: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2092: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2093: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2094: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
2095: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2096: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2097: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2098: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
2099: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2100: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2101: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2102: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
2103: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2104: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2105: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2106: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2107: {
2110: 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);
2111: 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);
2116: return(0);
2117: }
2119: /*@C
2120: PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field
2122: Not collective
2124: Input Parameters:
2125: + prob - The PetscDS
2126: . f - The test field number
2127: . g - The field number
2128: . g0 - integrand for the test and basis function term
2129: . g1 - integrand for the test function and basis function gradient term
2130: . g2 - integrand for the test function gradient and basis function term
2131: - g3 - integrand for the test function gradient and basis function gradient term
2133: Note: We are using a first order FEM model for the weak form:
2135: \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
2137: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
2139: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2140: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2141: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2142: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[])
2144: + dim - the spatial dimension
2145: . Nf - the number of fields
2146: . uOff - the offset into u[] and u_t[] for each field
2147: . uOff_x - the offset into u_x[] for each field
2148: . u - each field evaluated at the current point
2149: . u_t - the time derivative of each field evaluated at the current point
2150: . u_x - the gradient of each field evaluated at the current point
2151: . aOff - the offset into a[] and a_t[] for each auxiliary field
2152: . aOff_x - the offset into a_x[] for each auxiliary field
2153: . a - each auxiliary field evaluated at the current point
2154: . a_t - the time derivative of each auxiliary field evaluated at the current point
2155: . a_x - the gradient of auxiliary each field evaluated at the current point
2156: . t - current time
2157: . u_tShift - the multiplier a for dF/dU_t
2158: . x - coordinates of the current point
2159: . n - normal at the current point
2160: . numConstants - number of constant parameters
2161: . constants - constant parameters
2162: - g0 - output values at the current point
2164: Level: intermediate
2166: .seealso: PetscDSGetBdJacobian()
2167: @*/
2168: PetscErrorCode PetscDSSetBdJacobian(PetscDS prob, PetscInt f, PetscInt g,
2169: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2170: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2171: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2172: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
2173: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2174: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2175: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2176: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
2177: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2178: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2179: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2180: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
2181: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2182: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2183: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2184: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2185: {
2194: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2195: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
2196: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
2197: prob->gBd[(f*prob->Nf + g)*4+0] = g0;
2198: prob->gBd[(f*prob->Nf + g)*4+1] = g1;
2199: prob->gBd[(f*prob->Nf + g)*4+2] = g2;
2200: prob->gBd[(f*prob->Nf + g)*4+3] = g3;
2201: return(0);
2202: }
2204: /*@C
2205: PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field
2207: Not collective
2209: Input Parameters:
2210: + prob - The PetscDS
2211: - f - The test field number
2213: Output Parameter:
2214: . exactSol - exact solution for the test field
2216: Note: The calling sequence for the solution functions is given by:
2218: $ sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2220: + dim - the spatial dimension
2221: . t - current time
2222: . x - coordinates of the current point
2223: . Nc - the number of field components
2224: . u - the solution field evaluated at the current point
2225: - ctx - a user context
2227: Level: intermediate
2229: .seealso: PetscDSSetExactSolution()
2230: @*/
2231: PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx))
2232: {
2235: 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);
2237: return(0);
2238: }
2240: /*@C
2241: PetscDSSetExactSolution - Get the pointwise exact solution function for a given test field
2243: Not collective
2245: Input Parameters:
2246: + prob - The PetscDS
2247: . f - The test field number
2248: - sol - solution function for the test fields
2250: Note: The calling sequence for solution functions is given by:
2252: $ sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2254: + dim - the spatial dimension
2255: . t - current time
2256: . x - coordinates of the current point
2257: . Nc - the number of field components
2258: . u - the solution field evaluated at the current point
2259: - ctx - a user context
2261: Level: intermediate
2263: .seealso: PetscDSGetExactSolution()
2264: @*/
2265: PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx))
2266: {
2271: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2272: PetscDSEnlarge_Static(prob, f+1);
2274: return(0);
2275: }
2277: /*@C
2278: PetscDSGetConstants - Returns the array of constants passed to point functions
2280: Not collective
2282: Input Parameter:
2283: . prob - The PetscDS object
2285: Output Parameters:
2286: + numConstants - The number of constants
2287: - constants - The array of constants, NULL if there are none
2289: Level: intermediate
2291: .seealso: PetscDSSetConstants(), PetscDSCreate()
2292: @*/
2293: PetscErrorCode PetscDSGetConstants(PetscDS prob, PetscInt *numConstants, const PetscScalar *constants[])
2294: {
2299: return(0);
2300: }
2302: /*@C
2303: PetscDSSetConstants - Set the array of constants passed to point functions
2305: Not collective
2307: Input Parameters:
2308: + prob - The PetscDS object
2309: . numConstants - The number of constants
2310: - constants - The array of constants, NULL if there are none
2312: Level: intermediate
2314: .seealso: PetscDSGetConstants(), PetscDSCreate()
2315: @*/
2316: PetscErrorCode PetscDSSetConstants(PetscDS prob, PetscInt numConstants, PetscScalar constants[])
2317: {
2322: if (numConstants != prob->numConstants) {
2323: PetscFree(prob->constants);
2324: prob->numConstants = numConstants;
2325: if (prob->numConstants) {
2326: PetscMalloc1(prob->numConstants, &prob->constants);
2327: } else {
2328: prob->constants = NULL;
2329: }
2330: }
2331: if (prob->numConstants) {
2333: PetscMemcpy(prob->constants, constants, prob->numConstants * sizeof(PetscScalar));
2334: }
2335: return(0);
2336: }
2338: /*@
2339: PetscDSGetFieldIndex - Returns the index of the given field
2341: Not collective
2343: Input Parameters:
2344: + prob - The PetscDS object
2345: - disc - The discretization object
2347: Output Parameter:
2348: . f - The field number
2350: Level: beginner
2352: .seealso: PetscGetDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
2353: @*/
2354: PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f)
2355: {
2356: PetscInt g;
2361: *f = -1;
2362: for (g = 0; g < prob->Nf; ++g) {if (disc == prob->disc[g]) break;}
2363: if (g == prob->Nf) SETERRQ(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS.");
2364: *f = g;
2365: return(0);
2366: }
2368: /*@
2369: PetscDSGetFieldSize - Returns the size of the given field in the full space basis
2371: Not collective
2373: Input Parameters:
2374: + prob - The PetscDS object
2375: - f - The field number
2377: Output Parameter:
2378: . size - The size
2380: Level: beginner
2382: .seealso: PetscDSGetFieldOffset(), PetscDSGetNumFields(), PetscDSCreate()
2383: @*/
2384: PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size)
2385: {
2391: 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);
2392: PetscDSSetUp(prob);
2393: *size = prob->Nb[f];
2394: return(0);
2395: }
2397: /*@
2398: PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis
2400: Not collective
2402: Input Parameters:
2403: + prob - The PetscDS object
2404: - f - The field number
2406: Output Parameter:
2407: . off - The offset
2409: Level: beginner
2411: .seealso: PetscDSGetFieldSize(), PetscDSGetNumFields(), PetscDSCreate()
2412: @*/
2413: PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off)
2414: {
2415: PetscInt size, g;
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: *off = 0;
2423: for (g = 0; g < f; ++g) {
2424: PetscDSGetFieldSize(prob, g, &size);
2425: *off += size;
2426: }
2427: return(0);
2428: }
2430: /*@
2431: PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point
2433: Not collective
2435: Input Parameter:
2436: . prob - The PetscDS object
2438: Output Parameter:
2439: . dimensions - The number of dimensions
2441: Level: beginner
2443: .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate()
2444: @*/
2445: PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[])
2446: {
2451: PetscDSSetUp(prob);
2453: *dimensions = prob->Nb;
2454: return(0);
2455: }
2457: /*@
2458: PetscDSGetComponents - Returns the number of components for each field on an evaluation point
2460: Not collective
2462: Input Parameter:
2463: . prob - The PetscDS object
2465: Output Parameter:
2466: . components - The number of components
2468: Level: beginner
2470: .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate()
2471: @*/
2472: PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[])
2473: {
2478: PetscDSSetUp(prob);
2480: *components = prob->Nc;
2481: return(0);
2482: }
2484: /*@
2485: PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point
2487: Not collective
2489: Input Parameters:
2490: + prob - The PetscDS object
2491: - f - The field number
2493: Output Parameter:
2494: . off - The offset
2496: Level: beginner
2498: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2499: @*/
2500: PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off)
2501: {
2505: 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);
2506: *off = prob->off[f];
2507: return(0);
2508: }
2510: /*@
2511: PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point
2513: Not collective
2515: Input Parameter:
2516: . prob - The PetscDS object
2518: Output Parameter:
2519: . offsets - The offsets
2521: Level: beginner
2523: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2524: @*/
2525: PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[])
2526: {
2530: *offsets = prob->off;
2531: return(0);
2532: }
2534: /*@
2535: PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point
2537: Not collective
2539: Input Parameter:
2540: . prob - The PetscDS object
2542: Output Parameter:
2543: . offsets - The offsets
2545: Level: beginner
2547: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2548: @*/
2549: PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[])
2550: {
2554: *offsets = prob->offDer;
2555: return(0);
2556: }
2558: /*@C
2559: PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization
2561: Not collective
2563: Input Parameter:
2564: . prob - The PetscDS object
2566: Output Parameters:
2567: + basis - The basis function tabulation at quadrature points
2568: - basisDer - The basis function derivative tabulation at quadrature points
2570: Level: intermediate
2572: .seealso: PetscDSCreate()
2573: @*/
2574: PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer)
2575: {
2580: PetscDSSetUp(prob);
2583: return(0);
2584: }
2586: /*@C
2587: PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces
2589: Not collective
2591: Input Parameter:
2592: . prob - The PetscDS object
2594: Output Parameters:
2595: + basisFace - The basis function tabulation at quadrature points
2596: - basisDerFace - The basis function derivative tabulation at quadrature points
2598: Level: intermediate
2600: .seealso: PetscDSGetTabulation(), PetscDSCreate()
2601: @*/
2602: PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer)
2603: {
2608: PetscDSSetUp(prob);
2611: return(0);
2612: }
2614: PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x)
2615: {
2620: PetscDSSetUp(prob);
2624: return(0);
2625: }
2627: PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3)
2628: {
2633: PetscDSSetUp(prob);
2640: return(0);
2641: }
2643: PetscErrorCode PetscDSGetRefCoordArrays(PetscDS prob, PetscReal **x, PetscScalar **refSpaceDer)
2644: {
2649: PetscDSSetUp(prob);
2652: return(0);
2653: }
2655: /*@C
2656: PetscDSAddBoundary - Add a boundary condition to the model
2658: Input Parameters:
2659: + ds - The PetscDS object
2660: . type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2661: . name - The BC name
2662: . labelname - The label defining constrained points
2663: . field - The field to constrain
2664: . numcomps - The number of constrained field components (0 will constrain all fields)
2665: . comps - An array of constrained component numbers
2666: . bcFunc - A pointwise function giving boundary values
2667: . numids - The number of DMLabel ids for constrained points
2668: . ids - An array of ids for constrained points
2669: - ctx - An optional user context for bcFunc
2671: Options Database Keys:
2672: + -bc_<boundary name> <num> - Overrides the boundary ids
2673: - -bc_<boundary name>_comp <num> - Overrides the boundary components
2675: Level: developer
2677: .seealso: PetscDSGetBoundary()
2678: @*/
2679: 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)
2680: {
2681: DSBoundary b;
2686: PetscNew(&b);
2687: PetscStrallocpy(name, (char **) &b->name);
2688: PetscStrallocpy(labelname, (char **) &b->labelname);
2689: PetscMalloc1(numcomps, &b->comps);
2690: if (numcomps) {PetscMemcpy(b->comps, comps, numcomps*sizeof(PetscInt));}
2691: PetscMalloc1(numids, &b->ids);
2692: if (numids) {PetscMemcpy(b->ids, ids, numids*sizeof(PetscInt));}
2693: b->type = type;
2694: b->field = field;
2695: b->numcomps = numcomps;
2696: b->func = bcFunc;
2697: b->numids = numids;
2698: b->ctx = ctx;
2699: b->next = ds->boundary;
2700: ds->boundary = b;
2701: return(0);
2702: }
2704: /*@C
2705: PetscDSUpdateBoundary - Change a boundary condition for the model
2707: Input Parameters:
2708: + ds - The PetscDS object
2709: . bd - The boundary condition number
2710: . type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2711: . name - The BC name
2712: . labelname - The label defining constrained points
2713: . field - The field to constrain
2714: . numcomps - The number of constrained field components
2715: . comps - An array of constrained component numbers
2716: . bcFunc - A pointwise function giving boundary values
2717: . numids - The number of DMLabel ids for constrained points
2718: . ids - An array of ids for constrained points
2719: - ctx - An optional user context for bcFunc
2721: Note: The boundary condition number is the order in which it was registered. The user can get the number of boundary conditions from PetscDSGetNumBoundary().
2723: Level: developer
2725: .seealso: PetscDSAddBoundary(), PetscDSGetBoundary(), PetscDSGetNumBoundary()
2726: @*/
2727: 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)
2728: {
2729: DSBoundary b = ds->boundary;
2730: PetscInt n = 0;
2735: while (b) {
2736: if (n == bd) break;
2737: b = b->next;
2738: ++n;
2739: }
2740: if (!b) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %d is not in [0, %d)", bd, n);
2741: if (name) {
2742: PetscFree(b->name);
2743: PetscStrallocpy(name, (char **) &b->name);
2744: }
2745: if (labelname) {
2746: PetscFree(b->labelname);
2747: PetscStrallocpy(labelname, (char **) &b->labelname);
2748: }
2749: if (numcomps >= 0 && numcomps != b->numcomps) {
2750: b->numcomps = numcomps;
2751: PetscFree(b->comps);
2752: PetscMalloc1(numcomps, &b->comps);
2753: if (numcomps) {PetscMemcpy(b->comps, comps, numcomps*sizeof(PetscInt));}
2754: }
2755: if (numids >= 0 && numids != b->numids) {
2756: b->numids = numids;
2757: PetscFree(b->ids);
2758: PetscMalloc1(numids, &b->ids);
2759: if (numids) {PetscMemcpy(b->ids, ids, numids*sizeof(PetscInt));}
2760: }
2761: b->type = type;
2762: if (field >= 0) {b->field = field;}
2763: if (bcFunc) {b->func = bcFunc;}
2764: if (ctx) {b->ctx = ctx;}
2765: return(0);
2766: }
2768: /*@
2769: PetscDSGetNumBoundary - Get the number of registered BC
2771: Input Parameters:
2772: . ds - The PetscDS object
2774: Output Parameters:
2775: . numBd - The number of BC
2777: Level: intermediate
2779: .seealso: PetscDSAddBoundary(), PetscDSGetBoundary()
2780: @*/
2781: PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd)
2782: {
2783: DSBoundary b = ds->boundary;
2788: *numBd = 0;
2789: while (b) {++(*numBd); b = b->next;}
2790: return(0);
2791: }
2793: /*@C
2794: PetscDSGetBoundary - Gets a boundary condition to the model
2796: Input Parameters:
2797: + ds - The PetscDS object
2798: - bd - The BC number
2800: Output Parameters:
2801: + type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2802: . name - The BC name
2803: . labelname - The label defining constrained points
2804: . field - The field to constrain
2805: . numcomps - The number of constrained field components
2806: . comps - An array of constrained component numbers
2807: . bcFunc - A pointwise function giving boundary values
2808: . numids - The number of DMLabel ids for constrained points
2809: . ids - An array of ids for constrained points
2810: - ctx - An optional user context for bcFunc
2812: Options Database Keys:
2813: + -bc_<boundary name> <num> - Overrides the boundary ids
2814: - -bc_<boundary name>_comp <num> - Overrides the boundary components
2816: Level: developer
2818: .seealso: PetscDSAddBoundary()
2819: @*/
2820: 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)
2821: {
2822: DSBoundary b = ds->boundary;
2823: PetscInt n = 0;
2827: while (b) {
2828: if (n == bd) break;
2829: b = b->next;
2830: ++n;
2831: }
2832: if (!b) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %d is not in [0, %d)", bd, n);
2833: if (type) {
2835: *type = b->type;
2836: }
2837: if (name) {
2839: *name = b->name;
2840: }
2841: if (labelname) {
2843: *labelname = b->labelname;
2844: }
2845: if (field) {
2847: *field = b->field;
2848: }
2849: if (numcomps) {
2851: *numcomps = b->numcomps;
2852: }
2853: if (comps) {
2855: *comps = b->comps;
2856: }
2857: if (func) {
2859: *func = b->func;
2860: }
2861: if (numids) {
2863: *numids = b->numids;
2864: }
2865: if (ids) {
2867: *ids = b->ids;
2868: }
2869: if (ctx) {
2871: *ctx = b->ctx;
2872: }
2873: return(0);
2874: }
2876: /*@
2877: PetscDSCopyBoundary - Copy all boundary condition objects to the new problem
2879: Not collective
2881: Input Parameter:
2882: . prob - The PetscDS object
2884: Output Parameter:
2885: . newprob - The PetscDS copy
2887: Level: intermediate
2889: .seealso: PetscDSCopyEquations(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
2890: @*/
2891: PetscErrorCode PetscDSCopyBoundary(PetscDS probA, PetscDS probB)
2892: {
2893: DSBoundary b, next, *lastnext;
2899: if (probA == probB) return(0);
2900: next = probB->boundary;
2901: while (next) {
2902: DSBoundary b = next;
2904: next = b->next;
2905: PetscFree(b->comps);
2906: PetscFree(b->ids);
2907: PetscFree(b->name);
2908: PetscFree(b->labelname);
2909: PetscFree(b);
2910: }
2911: lastnext = &(probB->boundary);
2912: for (b = probA->boundary; b; b = b->next) {
2913: DSBoundary bNew;
2915: PetscNew(&bNew);
2916: bNew->numcomps = b->numcomps;
2917: PetscMalloc1(bNew->numcomps, &bNew->comps);
2918: PetscMemcpy(bNew->comps, b->comps, bNew->numcomps*sizeof(PetscInt));
2919: bNew->numids = b->numids;
2920: PetscMalloc1(bNew->numids, &bNew->ids);
2921: PetscMemcpy(bNew->ids, b->ids, bNew->numids*sizeof(PetscInt));
2922: PetscStrallocpy(b->labelname,(char **) &(bNew->labelname));
2923: PetscStrallocpy(b->name,(char **) &(bNew->name));
2924: bNew->ctx = b->ctx;
2925: bNew->type = b->type;
2926: bNew->field = b->field;
2927: bNew->func = b->func;
2929: *lastnext = bNew;
2930: lastnext = &(bNew->next);
2931: }
2932: return(0);
2933: }
2935: /*@C
2936: PetscDSSelectEquations - Copy pointwise function pointers to the new problem with different field layout
2938: Not collective
2940: Input Parameter:
2941: + prob - The PetscDS object
2942: . numFields - Number of new fields
2943: - fields - Old field number for each new field
2945: Output Parameter:
2946: . newprob - The PetscDS copy
2948: Level: intermediate
2950: .seealso: PetscDSCopyBoundary(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
2951: @*/
2952: PetscErrorCode PetscDSSelectEquations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
2953: {
2954: PetscInt Nf, Nfn, fn, gn;
2961: PetscDSGetNumFields(prob, &Nf);
2962: PetscDSGetNumFields(newprob, &Nfn);
2963: 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);
2964: for (fn = 0; fn < numFields; ++fn) {
2965: const PetscInt f = fields ? fields[fn] : fn;
2966: PetscPointFunc obj;
2967: PetscPointFunc f0, f1;
2968: PetscBdPointFunc f0Bd, f1Bd;
2969: PetscRiemannFunc r;
2971: if (f >= Nf) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Field %D must be in [0, %D)", f, Nf);
2972: PetscDSGetObjective(prob, f, &obj);
2973: PetscDSGetResidual(prob, f, &f0, &f1);
2974: PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd);
2975: PetscDSGetRiemannSolver(prob, f, &r);
2976: PetscDSSetObjective(newprob, fn, obj);
2977: PetscDSSetResidual(newprob, fn, f0, f1);
2978: PetscDSSetBdResidual(newprob, fn, f0Bd, f1Bd);
2979: PetscDSSetRiemannSolver(newprob, fn, r);
2980: for (gn = 0; gn < numFields; ++gn) {
2981: const PetscInt g = fields ? fields[gn] : gn;
2982: PetscPointJac g0, g1, g2, g3;
2983: PetscPointJac g0p, g1p, g2p, g3p;
2984: PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd;
2986: if (g >= Nf) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Field %D must be in [0, %D)", g, Nf);
2987: PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3);
2988: PetscDSGetJacobianPreconditioner(prob, f, g, &g0p, &g1p, &g2p, &g3p);
2989: PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd);
2990: PetscDSSetJacobian(newprob, fn, gn, g0, g1, g2, g3);
2991: PetscDSSetJacobianPreconditioner(prob, fn, gn, g0p, g1p, g2p, g3p);
2992: PetscDSSetBdJacobian(newprob, fn, gn, g0Bd, g1Bd, g2Bd, g3Bd);
2993: }
2994: }
2995: return(0);
2996: }
2998: /*@
2999: PetscDSCopyEquations - Copy all pointwise function pointers to the new problem
3001: Not collective
3003: Input Parameter:
3004: . prob - The PetscDS object
3006: Output Parameter:
3007: . newprob - The PetscDS copy
3009: Level: intermediate
3011: .seealso: PetscDSCopyBoundary(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
3012: @*/
3013: PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob)
3014: {
3015: PetscInt Nf, Ng;
3021: PetscDSGetNumFields(prob, &Nf);
3022: PetscDSGetNumFields(newprob, &Ng);
3023: if (Nf != Ng) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %D != %D", Nf, Ng);
3024: PetscDSSelectEquations(prob, Nf, NULL, newprob);
3025: return(0);
3026: }
3027: /*@
3028: PetscDSCopyConstants - Copy all constants to the new problem
3030: Not collective
3032: Input Parameter:
3033: . prob - The PetscDS object
3035: Output Parameter:
3036: . newprob - The PetscDS copy
3038: Level: intermediate
3040: .seealso: PetscDSCopyBoundary(), PetscDSCopyEquations(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
3041: @*/
3042: PetscErrorCode PetscDSCopyConstants(PetscDS prob, PetscDS newprob)
3043: {
3044: PetscInt Nc;
3045: const PetscScalar *constants;
3046: PetscErrorCode ierr;
3051: PetscDSGetConstants(prob, &Nc, &constants);
3052: PetscDSSetConstants(newprob, Nc, (PetscScalar *) constants);
3053: return(0);
3054: }
3056: PetscErrorCode PetscDSGetHeightSubspace(PetscDS prob, PetscInt height, PetscDS *subprob)
3057: {
3058: PetscInt dim, Nf, f;
3064: if (height == 0) {*subprob = prob; return(0);}
3065: PetscDSGetNumFields(prob, &Nf);
3066: PetscDSGetSpatialDimension(prob, &dim);
3067: if (height > dim) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_OUTOFRANGE, "DS can only handle height in [0, %D], not %D", dim, height);
3068: if (!prob->subprobs) {PetscCalloc1(dim, &prob->subprobs);}
3069: if (!prob->subprobs[height-1]) {
3070: PetscInt cdim;
3072: PetscDSCreate(PetscObjectComm((PetscObject) prob), &prob->subprobs[height-1]);
3073: PetscDSGetCoordinateDimension(prob, &cdim);
3074: PetscDSSetCoordinateDimension(prob->subprobs[height-1], cdim);
3075: for (f = 0; f < Nf; ++f) {
3076: PetscFE subfe;
3077: PetscObject obj;
3078: PetscClassId id;
3080: PetscDSGetDiscretization(prob, f, &obj);
3081: PetscObjectGetClassId(obj, &id);
3082: if (id == PETSCFE_CLASSID) {PetscFEGetHeightSubspace((PetscFE) obj, height, &subfe);}
3083: else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unsupported discretization type for field %d", f);
3084: PetscDSSetDiscretization(prob->subprobs[height-1], f, (PetscObject) subfe);
3085: }
3086: }
3087: *subprob = prob->subprobs[height-1];
3088: return(0);
3089: }
3091: static PetscErrorCode PetscDSDestroy_Basic(PetscDS prob)
3092: {
3093: PetscErrorCode ierr;
3096: PetscFree(prob->data);
3097: return(0);
3098: }
3100: static PetscErrorCode PetscDSInitialize_Basic(PetscDS prob)
3101: {
3103: prob->ops->setfromoptions = NULL;
3104: prob->ops->setup = NULL;
3105: prob->ops->view = NULL;
3106: prob->ops->destroy = PetscDSDestroy_Basic;
3107: return(0);
3108: }
3110: /*MC
3111: PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions
3113: Level: intermediate
3115: .seealso: PetscDSType, PetscDSCreate(), PetscDSSetType()
3116: M*/
3118: PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS prob)
3119: {
3120: PetscDS_Basic *b;
3121: PetscErrorCode ierr;
3125: PetscNewLog(prob, &b);
3126: prob->data = b;
3128: PetscDSInitialize_Basic(prob);
3129: return(0);
3130: }