Actual source code: dtds.c
petsc-3.10.5 2019-03-28
1: #include <petsc/private/petscdsimpl.h>
3: PetscClassId PETSCDS_CLASSID = 0;
5: PetscFunctionList PetscDSList = NULL;
6: PetscBool PetscDSRegisterAllCalled = PETSC_FALSE;
8: /*@C
9: PetscDSRegister - Adds a new PetscDS implementation
11: Not Collective
13: Input Parameters:
14: + name - The name of a new user-defined creation routine
15: - create_func - The creation routine itself
17: Notes:
18: PetscDSRegister() may be called multiple times to add several user-defined PetscDSs
20: Sample usage:
21: .vb
22: PetscDSRegister("my_ds", MyPetscDSCreate);
23: .ve
25: Then, your PetscDS type can be chosen with the procedural interface via
26: .vb
27: PetscDSCreate(MPI_Comm, PetscDS *);
28: PetscDSSetType(PetscDS, "my_ds");
29: .ve
30: or at runtime via the option
31: .vb
32: -petscds_type my_ds
33: .ve
35: Level: advanced
37: Not available from Fortran
39: .keywords: PetscDS, register
40: .seealso: PetscDSRegisterAll(), PetscDSRegisterDestroy()
42: @*/
43: PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS))
44: {
48: PetscFunctionListAdd(&PetscDSList, sname, function);
49: return(0);
50: }
52: /*@C
53: PetscDSSetType - Builds a particular PetscDS
55: Collective on PetscDS
57: Input Parameters:
58: + prob - The PetscDS object
59: - name - The kind of system
61: Options Database Key:
62: . -petscds_type <type> - Sets the PetscDS type; use -help for a list of available types
64: Level: intermediate
66: Not available from Fortran
68: .keywords: PetscDS, set, type
69: .seealso: PetscDSGetType(), PetscDSCreate()
70: @*/
71: PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name)
72: {
73: PetscErrorCode (*r)(PetscDS);
74: PetscBool match;
79: PetscObjectTypeCompare((PetscObject) prob, name, &match);
80: if (match) return(0);
82: PetscDSRegisterAll();
83: PetscFunctionListFind(PetscDSList, name, &r);
84: if (!r) SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscDS type: %s", name);
86: if (prob->ops->destroy) {
87: (*prob->ops->destroy)(prob);
88: prob->ops->destroy = NULL;
89: }
90: (*r)(prob);
91: PetscObjectChangeTypeName((PetscObject) prob, name);
92: return(0);
93: }
95: /*@C
96: PetscDSGetType - Gets the PetscDS type name (as a string) from the object.
98: Not Collective
100: Input Parameter:
101: . prob - The PetscDS
103: Output Parameter:
104: . name - The PetscDS type name
106: Level: intermediate
108: Not available from Fortran
110: .keywords: PetscDS, get, type, name
111: .seealso: PetscDSSetType(), PetscDSCreate()
112: @*/
113: PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name)
114: {
120: PetscDSRegisterAll();
121: *name = ((PetscObject) prob)->type_name;
122: return(0);
123: }
125: static PetscErrorCode PetscDSView_Ascii(PetscDS prob, PetscViewer viewer)
126: {
127: PetscViewerFormat format;
128: const PetscScalar *constants;
129: PetscInt numConstants, f;
130: PetscErrorCode ierr;
133: PetscViewerGetFormat(viewer, &format);
134: PetscViewerASCIIPrintf(viewer, "Discrete System with %d fields\n", prob->Nf);
135: PetscViewerASCIIPushTab(viewer);
136: for (f = 0; f < prob->Nf; ++f) {
137: PetscObject obj;
138: PetscClassId id;
139: PetscQuadrature q;
140: const char *name;
141: PetscInt Nc, Nq, Nqc;
143: PetscDSGetDiscretization(prob, f, &obj);
144: PetscObjectGetClassId(obj, &id);
145: PetscObjectGetName(obj, &name);
146: PetscViewerASCIIPrintf(viewer, "Field %s", name ? name : "<unknown>");
147: if (id == PETSCFE_CLASSID) {
148: PetscFEGetNumComponents((PetscFE) obj, &Nc);
149: PetscFEGetQuadrature((PetscFE) obj, &q);
150: PetscViewerASCIIPrintf(viewer, " FEM");
151: } else if (id == PETSCFV_CLASSID) {
152: PetscFVGetNumComponents((PetscFV) obj, &Nc);
153: PetscFVGetQuadrature((PetscFV) obj, &q);
154: PetscViewerASCIIPrintf(viewer, " FVM");
155: }
156: else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %D", f);
157: if (Nc > 1) {PetscViewerASCIIPrintf(viewer, "%D components", Nc);}
158: else {PetscViewerASCIIPrintf(viewer, "%D component ", Nc);}
159: if (prob->implicit[f]) {PetscViewerASCIIPrintf(viewer, " (implicit)");}
160: else {PetscViewerASCIIPrintf(viewer, " (explicit)");}
161: if (prob->adjacency[f*2+0]) {
162: if (prob->adjacency[f*2+1]) {PetscViewerASCIIPrintf(viewer, " (adj FVM++)");}
163: else {PetscViewerASCIIPrintf(viewer, " (adj FVM)");}
164: } else {
165: if (prob->adjacency[f*2+1]) {PetscViewerASCIIPrintf(viewer, " (adj FEM)");}
166: else {PetscViewerASCIIPrintf(viewer, " (adj FUNKY)");}
167: }
168: if (q) {
169: PetscQuadratureGetData(q, NULL, &Nqc, &Nq, NULL, NULL);
170: PetscViewerASCIIPrintf(viewer, " (Nq %D Nqc %D)", Nq, Nqc);
171: }
172: PetscViewerASCIIPrintf(viewer, "\n");
173: PetscViewerASCIIPushTab(viewer);
174: if (id == PETSCFE_CLASSID) {PetscFEView((PetscFE) obj, viewer);}
175: else if (id == PETSCFV_CLASSID) {PetscFVView((PetscFV) obj, viewer);}
176: PetscViewerASCIIPopTab(viewer);
177: }
178: PetscDSGetConstants(prob, &numConstants, &constants);
179: if (numConstants) {
180: PetscViewerASCIIPrintf(viewer, "%D constants\n", numConstants);
181: PetscViewerASCIIPushTab(viewer);
182: for (f = 0; f < numConstants; ++f) {PetscViewerASCIIPrintf(viewer, "%g\n", (double) PetscRealPart(constants[f]));}
183: PetscViewerASCIIPopTab(viewer);
184: }
185: PetscViewerASCIIPopTab(viewer);
186: return(0);
187: }
189: /*@C
190: PetscDSView - Views a PetscDS
192: Collective on PetscDS
194: Input Parameter:
195: + prob - the PetscDS object to view
196: - v - the viewer
198: Level: developer
200: .seealso PetscDSDestroy()
201: @*/
202: PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v)
203: {
204: PetscBool iascii;
209: if (!v) {PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject) prob), &v);}
211: PetscObjectTypeCompare((PetscObject) v, PETSCVIEWERASCII, &iascii);
212: if (iascii) {PetscDSView_Ascii(prob, v);}
213: if (prob->ops->view) {(*prob->ops->view)(prob, v);}
214: return(0);
215: }
217: /*@
218: PetscDSSetFromOptions - sets parameters in a PetscDS from the options database
220: Collective on PetscDS
222: Input Parameter:
223: . prob - the PetscDS object to set options for
225: Options Database:
226: + -petscds_type <type> : Set the DS type
227: . -petscds_view <view opt> : View the DS
228: . -petscds_jac_pre : Turn formation of a separate Jacobian preconditioner on and off
229: . -bc_<name> <ids> : Specify a list of label ids for a boundary condition
230: - -bc_<name>_comp <comps> : Specify a list of field components to constrain for a boundary condition
232: Level: developer
234: .seealso PetscDSView()
235: @*/
236: PetscErrorCode PetscDSSetFromOptions(PetscDS prob)
237: {
238: DSBoundary b;
239: const char *defaultType;
240: char name[256];
241: PetscBool flg;
246: if (!((PetscObject) prob)->type_name) {
247: defaultType = PETSCDSBASIC;
248: } else {
249: defaultType = ((PetscObject) prob)->type_name;
250: }
251: PetscDSRegisterAll();
253: PetscObjectOptionsBegin((PetscObject) prob);
254: for (b = prob->boundary; b; b = b->next) {
255: char optname[1024];
256: PetscInt ids[1024], len = 1024;
257: PetscBool flg;
259: PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name);
260: PetscMemzero(ids, sizeof(ids));
261: PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg);
262: if (flg) {
263: b->numids = len;
264: PetscFree(b->ids);
265: PetscMalloc1(len, &b->ids);
266: PetscMemcpy(b->ids, ids, len*sizeof(PetscInt));
267: }
268: len = 1024;
269: PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name);
270: PetscMemzero(ids, sizeof(ids));
271: PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg);
272: if (flg) {
273: b->numcomps = len;
274: PetscFree(b->comps);
275: PetscMalloc1(len, &b->comps);
276: PetscMemcpy(b->comps, ids, len*sizeof(PetscInt));
277: }
278: }
279: PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg);
280: if (flg) {
281: PetscDSSetType(prob, name);
282: } else if (!((PetscObject) prob)->type_name) {
283: PetscDSSetType(prob, defaultType);
284: }
285: PetscOptionsBool("-petscds_jac_pre", "Discrete System", "PetscDSUseJacobianPreconditioner", prob->useJacPre, &prob->useJacPre, &flg);
286: if (prob->ops->setfromoptions) {(*prob->ops->setfromoptions)(prob);}
287: /* process any options handlers added with PetscObjectAddOptionsHandler() */
288: PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject) prob);
289: PetscOptionsEnd();
290: if (prob->Nf) {PetscDSViewFromOptions(prob, NULL, "-petscds_view");}
291: return(0);
292: }
294: /*@C
295: PetscDSSetUp - Construct data structures for the PetscDS
297: Collective on PetscDS
299: Input Parameter:
300: . prob - the PetscDS object to setup
302: Level: developer
304: .seealso PetscDSView(), PetscDSDestroy()
305: @*/
306: PetscErrorCode PetscDSSetUp(PetscDS prob)
307: {
308: const PetscInt Nf = prob->Nf;
309: PetscInt dim, dimEmbed, work, NcMax = 0, NqMax = 0, f;
314: if (prob->setup) return(0);
315: /* Calculate sizes */
316: PetscDSGetSpatialDimension(prob, &dim);
317: PetscDSGetCoordinateDimension(prob, &dimEmbed);
318: prob->totDim = prob->totComp = 0;
319: PetscMalloc2(Nf,&prob->Nc,Nf,&prob->Nb);
320: PetscCalloc2(Nf+1,&prob->off,Nf+1,&prob->offDer);
321: PetscMalloc4(Nf,&prob->basis,Nf,&prob->basisDer,Nf,&prob->basisFace,Nf,&prob->basisDerFace);
322: for (f = 0; f < Nf; ++f) {
323: PetscObject obj;
324: PetscClassId id;
325: PetscQuadrature q;
326: PetscInt Nq = 0, Nb, Nc;
328: PetscDSGetDiscretization(prob, f, &obj);
329: PetscObjectGetClassId(obj, &id);
330: if (id == PETSCFE_CLASSID) {
331: PetscFE fe = (PetscFE) obj;
333: PetscFEGetQuadrature(fe, &q);
334: PetscFEGetDimension(fe, &Nb);
335: PetscFEGetNumComponents(fe, &Nc);
336: PetscFEGetDefaultTabulation(fe, &prob->basis[f], &prob->basisDer[f], NULL);
337: PetscFEGetFaceTabulation(fe, &prob->basisFace[f], &prob->basisDerFace[f], NULL);
338: } else if (id == PETSCFV_CLASSID) {
339: PetscFV fv = (PetscFV) obj;
341: PetscFVGetQuadrature(fv, &q);
342: PetscFVGetNumComponents(fv, &Nc);
343: Nb = Nc;
344: PetscFVGetDefaultTabulation(fv, &prob->basis[f], &prob->basisDer[f], NULL);
345: /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */
346: } else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", f);
347: prob->Nc[f] = Nc;
348: prob->Nb[f] = Nb;
349: prob->off[f+1] = Nc + prob->off[f];
350: prob->offDer[f+1] = Nc*dim + prob->offDer[f];
351: if (q) {PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL);}
352: NqMax = PetscMax(NqMax, Nq);
353: NcMax = PetscMax(NcMax, Nc);
354: prob->totDim += Nb;
355: prob->totComp += Nc;
356: }
357: work = PetscMax(prob->totComp*dim, PetscSqr(NcMax*dim));
358: /* Allocate works space */
359: PetscMalloc5(prob->totComp,&prob->u,prob->totComp,&prob->u_t,prob->totComp*dimEmbed,&prob->u_x,dimEmbed,&prob->x,work,&prob->refSpaceDer);
360: PetscMalloc6(NqMax*NcMax,&prob->f0,NqMax*NcMax*dim,&prob->f1,NqMax*NcMax*NcMax,&prob->g0,NqMax*NcMax*NcMax*dim,&prob->g1,NqMax*NcMax*NcMax*dim,&prob->g2,NqMax*NcMax*NcMax*dim*dim,&prob->g3);
361: if (prob->ops->setup) {(*prob->ops->setup)(prob);}
362: prob->setup = PETSC_TRUE;
363: return(0);
364: }
366: static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob)
367: {
371: PetscFree2(prob->Nc,prob->Nb);
372: PetscFree2(prob->off,prob->offDer);
373: PetscFree4(prob->basis,prob->basisDer,prob->basisFace,prob->basisDerFace);
374: PetscFree5(prob->u,prob->u_t,prob->u_x,prob->x,prob->refSpaceDer);
375: PetscFree6(prob->f0,prob->f1,prob->g0,prob->g1,prob->g2,prob->g3);
376: return(0);
377: }
379: static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew)
380: {
381: PetscObject *tmpd;
382: PetscBool *tmpi, *tmpa;
383: PetscPointFunc *tmpobj, *tmpf, *tmpup;
384: PetscPointJac *tmpg, *tmpgp, *tmpgt;
385: PetscBdPointFunc *tmpfbd;
386: PetscBdPointJac *tmpgbd;
387: PetscRiemannFunc *tmpr;
388: PetscSimplePointFunc *tmpexactSol;
389: void **tmpctx;
390: PetscInt Nf = prob->Nf, f, i;
391: PetscErrorCode ierr;
394: if (Nf >= NfNew) return(0);
395: prob->setup = PETSC_FALSE;
396: PetscDSDestroyStructs_Static(prob);
397: PetscMalloc3(NfNew, &tmpd, NfNew, &tmpi, NfNew*2, &tmpa);
398: for (f = 0; f < Nf; ++f) {tmpd[f] = prob->disc[f]; tmpi[f] = prob->implicit[f]; for (i = 0; i < 2; ++i) tmpa[f*2+i] = prob->adjacency[f*2+i];}
399: for (f = Nf; f < NfNew; ++f) {tmpd[f] = NULL; tmpi[f] = PETSC_TRUE; tmpa[f*2+0] = PETSC_FALSE; tmpa[f*2+1] = PETSC_TRUE;}
400: PetscFree3(prob->disc, prob->implicit, prob->adjacency);
401: prob->Nf = NfNew;
402: prob->disc = tmpd;
403: prob->implicit = tmpi;
404: prob->adjacency = tmpa;
405: PetscCalloc7(NfNew, &tmpobj, NfNew*2, &tmpf, NfNew*NfNew*4, &tmpg, NfNew*NfNew*4, &tmpgp, NfNew*NfNew*4, &tmpgt, NfNew, &tmpr, NfNew, &tmpctx);
406: PetscCalloc1(NfNew, &tmpup);
407: for (f = 0; f < Nf; ++f) tmpobj[f] = prob->obj[f];
408: for (f = 0; f < Nf*2; ++f) tmpf[f] = prob->f[f];
409: for (f = 0; f < Nf*Nf*4; ++f) tmpg[f] = prob->g[f];
410: for (f = 0; f < Nf*Nf*4; ++f) tmpgp[f] = prob->gp[f];
411: for (f = 0; f < Nf; ++f) tmpr[f] = prob->r[f];
412: for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f];
413: for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f];
414: for (f = Nf; f < NfNew; ++f) tmpobj[f] = NULL;
415: for (f = Nf*2; f < NfNew*2; ++f) tmpf[f] = NULL;
416: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpg[f] = NULL;
417: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgp[f] = NULL;
418: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgt[f] = NULL;
419: for (f = Nf; f < NfNew; ++f) tmpr[f] = NULL;
420: for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL;
421: for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL;
422: PetscFree7(prob->obj, prob->f, prob->g, prob->gp, prob->gt, prob->r, prob->ctx);
423: PetscFree(prob->update);
424: prob->obj = tmpobj;
425: prob->f = tmpf;
426: prob->g = tmpg;
427: prob->gp = tmpgp;
428: prob->gt = tmpgt;
429: prob->r = tmpr;
430: prob->update = tmpup;
431: prob->ctx = tmpctx;
432: PetscCalloc3(NfNew*2, &tmpfbd, NfNew*NfNew*4, &tmpgbd, NfNew, &tmpexactSol);
433: for (f = 0; f < Nf*2; ++f) tmpfbd[f] = prob->fBd[f];
434: for (f = 0; f < Nf*Nf*4; ++f) tmpgbd[f] = prob->gBd[f];
435: for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f];
436: for (f = Nf*2; f < NfNew*2; ++f) tmpfbd[f] = NULL;
437: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgbd[f] = NULL;
438: for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL;
439: PetscFree3(prob->fBd, prob->gBd, prob->exactSol);
440: prob->fBd = tmpfbd;
441: prob->gBd = tmpgbd;
442: prob->exactSol = tmpexactSol;
443: return(0);
444: }
446: /*@
447: PetscDSDestroy - Destroys a PetscDS object
449: Collective on PetscDS
451: Input Parameter:
452: . prob - the PetscDS object to destroy
454: Level: developer
456: .seealso PetscDSView()
457: @*/
458: PetscErrorCode PetscDSDestroy(PetscDS *prob)
459: {
460: PetscInt f;
461: DSBoundary next;
465: if (!*prob) return(0);
468: if (--((PetscObject)(*prob))->refct > 0) {*prob = 0; return(0);}
469: ((PetscObject) (*prob))->refct = 0;
470: if ((*prob)->subprobs) {
471: PetscInt dim, d;
473: PetscDSGetSpatialDimension(*prob, &dim);
474: for (d = 0; d < dim; ++d) {PetscDSDestroy(&(*prob)->subprobs[d]);}
475: }
476: PetscFree((*prob)->subprobs);
477: PetscDSDestroyStructs_Static(*prob);
478: for (f = 0; f < (*prob)->Nf; ++f) {
479: PetscObjectDereference((*prob)->disc[f]);
480: }
481: PetscFree3((*prob)->disc, (*prob)->implicit, (*prob)->adjacency);
482: PetscFree7((*prob)->obj,(*prob)->f,(*prob)->g,(*prob)->gp,(*prob)->gt,(*prob)->r,(*prob)->ctx);
483: PetscFree((*prob)->update);
484: PetscFree3((*prob)->fBd,(*prob)->gBd,(*prob)->exactSol);
485: if ((*prob)->ops->destroy) {(*(*prob)->ops->destroy)(*prob);}
486: next = (*prob)->boundary;
487: while (next) {
488: DSBoundary b = next;
490: next = b->next;
491: PetscFree(b->comps);
492: PetscFree(b->ids);
493: PetscFree(b->name);
494: PetscFree(b->labelname);
495: PetscFree(b);
496: }
497: PetscFree((*prob)->constants);
498: PetscHeaderDestroy(prob);
499: return(0);
500: }
502: /*@
503: PetscDSCreate - Creates an empty PetscDS object. The type can then be set with PetscDSSetType().
505: Collective on MPI_Comm
507: Input Parameter:
508: . comm - The communicator for the PetscDS object
510: Output Parameter:
511: . prob - The PetscDS object
513: Level: beginner
515: .seealso: PetscDSSetType(), PETSCDSBASIC
516: @*/
517: PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *prob)
518: {
519: PetscDS p;
524: *prob = NULL;
525: PetscDSInitializePackage();
527: PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView);
529: p->Nf = 0;
530: p->setup = PETSC_FALSE;
531: p->numConstants = 0;
532: p->constants = NULL;
533: p->dimEmbed = -1;
534: p->defaultAdj[0] = PETSC_FALSE;
535: p->defaultAdj[1] = PETSC_TRUE;
536: p->useJacPre = PETSC_TRUE;
538: *prob = p;
539: return(0);
540: }
542: /*@
543: PetscDSGetNumFields - Returns the number of fields in the DS
545: Not collective
547: Input Parameter:
548: . prob - The PetscDS object
550: Output Parameter:
551: . Nf - The number of fields
553: Level: beginner
555: .seealso: PetscDSGetSpatialDimension(), PetscDSCreate()
556: @*/
557: PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf)
558: {
562: *Nf = prob->Nf;
563: return(0);
564: }
566: /*@
567: PetscDSGetSpatialDimension - Returns the spatial dimension of the DS, meaning the topological dimension of the discretizations
569: Not collective
571: Input Parameter:
572: . prob - The PetscDS object
574: Output Parameter:
575: . dim - The spatial dimension
577: Level: beginner
579: .seealso: PetscDSGetCoordinateDimension(), PetscDSGetNumFields(), PetscDSCreate()
580: @*/
581: PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim)
582: {
588: *dim = 0;
589: if (prob->Nf) {
590: PetscObject obj;
591: PetscClassId id;
593: PetscDSGetDiscretization(prob, 0, &obj);
594: PetscObjectGetClassId(obj, &id);
595: if (id == PETSCFE_CLASSID) {PetscFEGetSpatialDimension((PetscFE) obj, dim);}
596: else if (id == PETSCFV_CLASSID) {PetscFVGetSpatialDimension((PetscFV) obj, dim);}
597: else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
598: }
599: return(0);
600: }
602: /*@
603: PetscDSGetCoordinateDimension - Returns the coordinate dimension of the DS, meaning the dimension of the space into which the discretiaztions are embedded
605: Not collective
607: Input Parameter:
608: . prob - The PetscDS object
610: Output Parameter:
611: . dimEmbed - The coordinate dimension
613: Level: beginner
615: .seealso: PetscDSSetCoordinateDimension(), PetscDSGetSpatialDimension(), PetscDSGetNumFields(), PetscDSCreate()
616: @*/
617: PetscErrorCode PetscDSGetCoordinateDimension(PetscDS prob, PetscInt *dimEmbed)
618: {
622: if (prob->dimEmbed < 0) SETERRQ(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONGSTATE, "No coordinate dimension set for this DS");
623: *dimEmbed = prob->dimEmbed;
624: return(0);
625: }
627: /*@
628: PetscDSSetCoordinateDimension - Set the coordinate dimension of the DS, meaning the dimension of the space into which the discretiaztions are embedded
630: Not collective
632: Input Parameters:
633: + prob - The PetscDS object
634: - dimEmbed - The coordinate dimension
636: Level: beginner
638: .seealso: PetscDSGetCoordinateDimension(), PetscDSGetSpatialDimension(), PetscDSGetNumFields(), PetscDSCreate()
639: @*/
640: PetscErrorCode PetscDSSetCoordinateDimension(PetscDS prob, PetscInt dimEmbed)
641: {
644: prob->dimEmbed = dimEmbed;
645: return(0);
646: }
648: /*@
649: PetscDSGetTotalDimension - Returns the total size of the approximation space for this system
651: Not collective
653: Input Parameter:
654: . prob - The PetscDS object
656: Output Parameter:
657: . dim - The total problem dimension
659: Level: beginner
661: .seealso: PetscDSGetNumFields(), PetscDSCreate()
662: @*/
663: PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim)
664: {
669: PetscDSSetUp(prob);
671: *dim = prob->totDim;
672: return(0);
673: }
675: /*@
676: PetscDSGetTotalComponents - Returns the total number of components in this system
678: Not collective
680: Input Parameter:
681: . prob - The PetscDS object
683: Output Parameter:
684: . dim - The total number of components
686: Level: beginner
688: .seealso: PetscDSGetNumFields(), PetscDSCreate()
689: @*/
690: PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc)
691: {
696: PetscDSSetUp(prob);
698: *Nc = prob->totComp;
699: return(0);
700: }
702: /*@
703: PetscDSGetDiscretization - Returns the discretization object for the given field
705: Not collective
707: Input Parameters:
708: + prob - The PetscDS object
709: - f - The field number
711: Output Parameter:
712: . disc - The discretization object
714: Level: beginner
716: .seealso: PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
717: @*/
718: PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc)
719: {
723: 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);
724: *disc = prob->disc[f];
725: return(0);
726: }
728: /*@
729: PetscDSSetDiscretization - Sets the discretization object for the given field
731: Not collective
733: Input Parameters:
734: + prob - The PetscDS object
735: . f - The field number
736: - disc - The discretization object
738: Level: beginner
740: .seealso: PetscDSGetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
741: @*/
742: PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc)
743: {
749: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
750: PetscDSEnlarge_Static(prob, f+1);
751: if (prob->disc[f]) {PetscObjectDereference(prob->disc[f]);}
752: prob->disc[f] = disc;
753: PetscObjectReference(disc);
754: {
755: PetscClassId id;
757: PetscObjectGetClassId(disc, &id);
758: if (id == PETSCFE_CLASSID) {
759: PetscDSSetImplicit(prob, f, PETSC_TRUE);
760: PetscDSSetAdjacency(prob, f, PETSC_FALSE, PETSC_TRUE);
761: } else if (id == PETSCFV_CLASSID) {
762: PetscDSSetImplicit(prob, f, PETSC_FALSE);
763: PetscDSSetAdjacency(prob, f, PETSC_TRUE, PETSC_FALSE);
764: }
765: }
766: return(0);
767: }
769: /*@
770: PetscDSAddDiscretization - Adds a discretization object
772: Not collective
774: Input Parameters:
775: + prob - The PetscDS object
776: - disc - The boundary discretization object
778: Level: beginner
780: .seealso: PetscDSGetDiscretization(), PetscDSSetDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
781: @*/
782: PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc)
783: {
787: PetscDSSetDiscretization(prob, prob->Nf, disc);
788: return(0);
789: }
791: /*@
792: PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for IMEX
794: Not collective
796: Input Parameters:
797: + prob - The PetscDS object
798: - f - The field number
800: Output Parameter:
801: . implicit - The flag indicating what kind of solve to use for this field
803: Level: developer
805: .seealso: PetscDSSetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
806: @*/
807: PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit)
808: {
812: 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);
813: *implicit = prob->implicit[f];
814: return(0);
815: }
817: /*@
818: PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for IMEX
820: Not collective
822: Input Parameters:
823: + prob - The PetscDS object
824: . f - The field number
825: - implicit - The flag indicating what kind of solve to use for this field
827: Level: developer
829: .seealso: PetscDSGetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
830: @*/
831: PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit)
832: {
835: 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);
836: prob->implicit[f] = implicit;
837: return(0);
838: }
840: /*@
841: PetscDSGetAdjacency - Returns the flags for determining variable influence
843: Not collective
845: Input Parameters:
846: + prob - The PetscDS object
847: - f - The field number
849: Output Parameter:
850: + useCone - Flag for variable influence starting with the cone operation
851: - useClosure - Flag for variable influence using transitive closure
853: Note: See the discussion in DMPlexGetAdjacencyUseCone() and DMPlexGetAdjacencyUseClosure()
855: Level: developer
857: .seealso: PetscDSSetAdjacency(), DMPlexGetAdjacencyUseCone(), DMPlexGetAdjacencyUseClosure(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
858: @*/
859: PetscErrorCode PetscDSGetAdjacency(PetscDS prob, PetscInt f, PetscBool *useCone, PetscBool *useClosure)
860: {
865: if (f < 0) {
866: if (useCone) *useCone = prob->defaultAdj[0];
867: if (useClosure) *useClosure = prob->defaultAdj[1];
868: } else {
869: if (f >= prob->Nf) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
870: if (useCone) *useCone = prob->adjacency[f*2+0];
871: if (useClosure) *useClosure = prob->adjacency[f*2+1];
872: }
873: return(0);
874: }
876: /*@
877: PetscDSSetAdjacency - Set the flags for determining variable influence
879: Not collective
881: Input Parameters:
882: + prob - The PetscDS object
883: . f - The field number
884: . useCone - Flag for variable influence starting with the cone operation
885: - useClosure - Flag for variable influence using transitive closure
887: Note: See the discussion in DMPlexGetAdjacencyUseCone() and DMPlexGetAdjacencyUseClosure()
889: Level: developer
891: .seealso: PetscDSGetAdjacency(), DMPlexGetAdjacencyUseCone(), DMPlexGetAdjacencyUseClosure(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
892: @*/
893: PetscErrorCode PetscDSSetAdjacency(PetscDS prob, PetscInt f, PetscBool useCone, PetscBool useClosure)
894: {
897: if (f < 0) {
898: prob->defaultAdj[0] = useCone;
899: prob->defaultAdj[1] = useClosure;
900: } else {
901: if (f >= prob->Nf) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
902: prob->adjacency[f*2+0] = useCone;
903: prob->adjacency[f*2+1] = useClosure;
904: }
905: return(0);
906: }
908: PetscErrorCode PetscDSGetObjective(PetscDS prob, PetscInt f,
909: void (**obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
910: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
911: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
912: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
913: {
917: 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);
918: *obj = prob->obj[f];
919: return(0);
920: }
922: PetscErrorCode PetscDSSetObjective(PetscDS prob, PetscInt f,
923: void (*obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
924: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
925: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
926: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
927: {
933: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
934: PetscDSEnlarge_Static(prob, f+1);
935: prob->obj[f] = obj;
936: return(0);
937: }
939: /*@C
940: PetscDSGetResidual - Get the pointwise residual function for a given test field
942: Not collective
944: Input Parameters:
945: + prob - The PetscDS
946: - f - The test field number
948: Output Parameters:
949: + f0 - integrand for the test function term
950: - f1 - integrand for the test function gradient term
952: Note: We are using a first order FEM model for the weak form:
954: \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)
956: The calling sequence for the callbacks f0 and f1 is given by:
958: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
959: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
960: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
961: $ PetscReal t, const PetscReal x[], PetscScalar f0[])
963: + dim - the spatial dimension
964: . Nf - the number of fields
965: . uOff - the offset into u[] and u_t[] for each field
966: . uOff_x - the offset into u_x[] for each field
967: . u - each field evaluated at the current point
968: . u_t - the time derivative of each field evaluated at the current point
969: . u_x - the gradient of each field evaluated at the current point
970: . aOff - the offset into a[] and a_t[] for each auxiliary field
971: . aOff_x - the offset into a_x[] for each auxiliary field
972: . a - each auxiliary field evaluated at the current point
973: . a_t - the time derivative of each auxiliary field evaluated at the current point
974: . a_x - the gradient of auxiliary each field evaluated at the current point
975: . t - current time
976: . x - coordinates of the current point
977: . numConstants - number of constant parameters
978: . constants - constant parameters
979: - f0 - output values at the current point
981: Level: intermediate
983: .seealso: PetscDSSetResidual()
984: @*/
985: PetscErrorCode PetscDSGetResidual(PetscDS prob, PetscInt f,
986: void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
987: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
988: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
989: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
990: void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
991: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
992: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
993: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
994: {
997: 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);
1000: return(0);
1001: }
1003: /*@C
1004: PetscDSSetResidual - Set the pointwise residual function for a given test field
1006: Not collective
1008: Input Parameters:
1009: + prob - The PetscDS
1010: . f - The test field number
1011: . f0 - integrand for the test function term
1012: - f1 - integrand for the test function gradient term
1014: Note: We are using a first order FEM model for the weak form:
1016: \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)
1018: The calling sequence for the callbacks f0 and f1 is given by:
1020: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1021: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1022: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1023: $ PetscReal t, const PetscReal x[], PetscScalar f0[])
1025: + dim - the spatial dimension
1026: . Nf - the number of fields
1027: . uOff - the offset into u[] and u_t[] for each field
1028: . uOff_x - the offset into u_x[] for each field
1029: . u - each field evaluated at the current point
1030: . u_t - the time derivative of each field evaluated at the current point
1031: . u_x - the gradient of each field evaluated at the current point
1032: . aOff - the offset into a[] and a_t[] for each auxiliary field
1033: . aOff_x - the offset into a_x[] for each auxiliary field
1034: . a - each auxiliary field evaluated at the current point
1035: . a_t - the time derivative of each auxiliary field evaluated at the current point
1036: . a_x - the gradient of auxiliary each field evaluated at the current point
1037: . t - current time
1038: . x - coordinates of the current point
1039: . numConstants - number of constant parameters
1040: . constants - constant parameters
1041: - f0 - output values at the current point
1043: Level: intermediate
1045: .seealso: PetscDSGetResidual()
1046: @*/
1047: PetscErrorCode PetscDSSetResidual(PetscDS prob, PetscInt f,
1048: void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1049: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1050: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1051: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1052: void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1053: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1054: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1055: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1056: {
1063: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1064: PetscDSEnlarge_Static(prob, f+1);
1065: prob->f[f*2+0] = f0;
1066: prob->f[f*2+1] = f1;
1067: return(0);
1068: }
1070: /*@C
1071: PetscDSHasJacobian - Signals that Jacobian functions have been set
1073: Not collective
1075: Input Parameter:
1076: . prob - The PetscDS
1078: Output Parameter:
1079: . hasJac - flag that pointwise function for the Jacobian has been set
1081: Level: intermediate
1083: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1084: @*/
1085: PetscErrorCode PetscDSHasJacobian(PetscDS prob, PetscBool *hasJac)
1086: {
1087: PetscInt f, g, h;
1091: *hasJac = PETSC_FALSE;
1092: for (f = 0; f < prob->Nf; ++f) {
1093: for (g = 0; g < prob->Nf; ++g) {
1094: for (h = 0; h < 4; ++h) {
1095: if (prob->g[(f*prob->Nf + g)*4+h]) *hasJac = PETSC_TRUE;
1096: }
1097: }
1098: }
1099: return(0);
1100: }
1102: /*@C
1103: PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field
1105: Not collective
1107: Input Parameters:
1108: + prob - The PetscDS
1109: . f - The test field number
1110: - g - The field number
1112: Output Parameters:
1113: + g0 - integrand for the test and basis function term
1114: . g1 - integrand for the test function and basis function gradient term
1115: . g2 - integrand for the test function gradient and basis function term
1116: - g3 - integrand for the test function gradient and basis function gradient term
1118: Note: We are using a first order FEM model for the weak form:
1120: \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
1122: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1124: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1125: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1126: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1127: $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1129: + dim - the spatial dimension
1130: . Nf - the number of fields
1131: . uOff - the offset into u[] and u_t[] for each field
1132: . uOff_x - the offset into u_x[] for each field
1133: . u - each field evaluated at the current point
1134: . u_t - the time derivative of each field evaluated at the current point
1135: . u_x - the gradient of each field evaluated at the current point
1136: . aOff - the offset into a[] and a_t[] for each auxiliary field
1137: . aOff_x - the offset into a_x[] for each auxiliary field
1138: . a - each auxiliary field evaluated at the current point
1139: . a_t - the time derivative of each auxiliary field evaluated at the current point
1140: . a_x - the gradient of auxiliary each field evaluated at the current point
1141: . t - current time
1142: . u_tShift - the multiplier a for dF/dU_t
1143: . x - coordinates of the current point
1144: . numConstants - number of constant parameters
1145: . constants - constant parameters
1146: - g0 - output values at the current point
1148: Level: intermediate
1150: .seealso: PetscDSSetJacobian()
1151: @*/
1152: PetscErrorCode PetscDSGetJacobian(PetscDS prob, PetscInt f, PetscInt g,
1153: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1154: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1155: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1156: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1157: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1158: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1159: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1160: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1161: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1162: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1163: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1164: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1165: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1166: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1167: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1168: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1169: {
1172: 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);
1173: 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);
1178: return(0);
1179: }
1181: /*@C
1182: PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields
1184: Not collective
1186: Input Parameters:
1187: + prob - The PetscDS
1188: . f - The test field number
1189: . g - The field number
1190: . g0 - integrand for the test and basis function term
1191: . g1 - integrand for the test function and basis function gradient term
1192: . g2 - integrand for the test function gradient and basis function term
1193: - g3 - integrand for the test function gradient and basis function gradient term
1195: Note: We are using a first order FEM model for the weak form:
1197: \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
1199: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1201: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1202: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1203: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1204: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1206: + dim - the spatial dimension
1207: . Nf - the number of fields
1208: . uOff - the offset into u[] and u_t[] for each field
1209: . uOff_x - the offset into u_x[] for each field
1210: . u - each field evaluated at the current point
1211: . u_t - the time derivative of each field evaluated at the current point
1212: . u_x - the gradient of each field evaluated at the current point
1213: . aOff - the offset into a[] and a_t[] for each auxiliary field
1214: . aOff_x - the offset into a_x[] for each auxiliary field
1215: . a - each auxiliary field evaluated at the current point
1216: . a_t - the time derivative of each auxiliary field evaluated at the current point
1217: . a_x - the gradient of auxiliary each field evaluated at the current point
1218: . t - current time
1219: . u_tShift - the multiplier a for dF/dU_t
1220: . x - coordinates of the current point
1221: . numConstants - number of constant parameters
1222: . constants - constant parameters
1223: - g0 - output values at the current point
1225: Level: intermediate
1227: .seealso: PetscDSGetJacobian()
1228: @*/
1229: PetscErrorCode PetscDSSetJacobian(PetscDS prob, PetscInt f, PetscInt g,
1230: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1231: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1232: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1233: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1234: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1235: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1236: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1237: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1238: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1239: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1240: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1241: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1242: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1243: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1244: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1245: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1246: {
1255: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1256: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1257: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1258: prob->g[(f*prob->Nf + g)*4+0] = g0;
1259: prob->g[(f*prob->Nf + g)*4+1] = g1;
1260: prob->g[(f*prob->Nf + g)*4+2] = g2;
1261: prob->g[(f*prob->Nf + g)*4+3] = g3;
1262: return(0);
1263: }
1265: /*@C
1266: PetscDSUseJacobianPreconditioner - Whether to construct a Jacobian preconditioner
1268: Not collective
1270: Input Parameters:
1271: + prob - The PetscDS
1272: - useJacPre - flag that enables construction of a Jacobian preconditioner
1274: Level: intermediate
1276: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1277: @*/
1278: PetscErrorCode PetscDSUseJacobianPreconditioner(PetscDS prob, PetscBool useJacPre)
1279: {
1282: prob->useJacPre = useJacPre;
1283: return(0);
1284: }
1286: /*@C
1287: PetscDSHasJacobianPreconditioner - Signals that a Jacobian preconditioner matrix has been set
1289: Not collective
1291: Input Parameter:
1292: . prob - The PetscDS
1294: Output Parameter:
1295: . hasJacPre - flag that pointwise function for Jacobian preconditioner matrix has been set
1297: Level: intermediate
1299: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1300: @*/
1301: PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS prob, PetscBool *hasJacPre)
1302: {
1303: PetscInt f, g, h;
1307: *hasJacPre = PETSC_FALSE;
1308: if (!prob->useJacPre) return(0);
1309: for (f = 0; f < prob->Nf; ++f) {
1310: for (g = 0; g < prob->Nf; ++g) {
1311: for (h = 0; h < 4; ++h) {
1312: if (prob->gp[(f*prob->Nf + g)*4+h]) *hasJacPre = PETSC_TRUE;
1313: }
1314: }
1315: }
1316: return(0);
1317: }
1319: /*@C
1320: 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.
1322: Not collective
1324: Input Parameters:
1325: + prob - The PetscDS
1326: . f - The test field number
1327: - g - The field number
1329: Output Parameters:
1330: + g0 - integrand for the test and basis function term
1331: . g1 - integrand for the test function and basis function gradient term
1332: . g2 - integrand for the test function gradient and basis function term
1333: - g3 - integrand for the test function gradient and basis function gradient term
1335: Note: We are using a first order FEM model for the weak form:
1337: \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
1339: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1341: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1342: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1343: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1344: $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1346: + dim - the spatial dimension
1347: . Nf - the number of fields
1348: . uOff - the offset into u[] and u_t[] for each field
1349: . uOff_x - the offset into u_x[] for each field
1350: . u - each field evaluated at the current point
1351: . u_t - the time derivative of each field evaluated at the current point
1352: . u_x - the gradient of each field evaluated at the current point
1353: . aOff - the offset into a[] and a_t[] for each auxiliary field
1354: . aOff_x - the offset into a_x[] for each auxiliary field
1355: . a - each auxiliary field evaluated at the current point
1356: . a_t - the time derivative of each auxiliary field evaluated at the current point
1357: . a_x - the gradient of auxiliary each field evaluated at the current point
1358: . t - current time
1359: . u_tShift - the multiplier a for dF/dU_t
1360: . x - coordinates of the current point
1361: . numConstants - number of constant parameters
1362: . constants - constant parameters
1363: - g0 - output values at the current point
1365: Level: intermediate
1367: .seealso: PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1368: @*/
1369: PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g,
1370: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1371: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1372: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1373: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1374: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1375: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1376: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1377: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1378: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1379: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1380: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1381: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1382: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1383: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1384: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1385: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1386: {
1389: 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);
1390: 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);
1395: return(0);
1396: }
1398: /*@C
1399: 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.
1401: Not collective
1403: Input Parameters:
1404: + prob - The PetscDS
1405: . f - The test field number
1406: . g - The field number
1407: . g0 - integrand for the test and basis function term
1408: . g1 - integrand for the test function and basis function gradient term
1409: . g2 - integrand for the test function gradient and basis function term
1410: - g3 - integrand for the test function gradient and basis function gradient term
1412: Note: We are using a first order FEM model for the weak form:
1414: \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
1416: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1418: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1419: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1420: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1421: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1423: + dim - the spatial dimension
1424: . Nf - the number of fields
1425: . uOff - the offset into u[] and u_t[] for each field
1426: . uOff_x - the offset into u_x[] for each field
1427: . u - each field evaluated at the current point
1428: . u_t - the time derivative of each field evaluated at the current point
1429: . u_x - the gradient of each field evaluated at the current point
1430: . aOff - the offset into a[] and a_t[] for each auxiliary field
1431: . aOff_x - the offset into a_x[] for each auxiliary field
1432: . a - each auxiliary field evaluated at the current point
1433: . a_t - the time derivative of each auxiliary field evaluated at the current point
1434: . a_x - the gradient of auxiliary each field evaluated at the current point
1435: . t - current time
1436: . u_tShift - the multiplier a for dF/dU_t
1437: . x - coordinates of the current point
1438: . numConstants - number of constant parameters
1439: . constants - constant parameters
1440: - g0 - output values at the current point
1442: Level: intermediate
1444: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobian()
1445: @*/
1446: PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g,
1447: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1448: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1449: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1450: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1451: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1452: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1453: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1454: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1455: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1456: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1457: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1458: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1459: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1460: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1461: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1462: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1463: {
1472: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1473: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1474: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1475: prob->gp[(f*prob->Nf + g)*4+0] = g0;
1476: prob->gp[(f*prob->Nf + g)*4+1] = g1;
1477: prob->gp[(f*prob->Nf + g)*4+2] = g2;
1478: prob->gp[(f*prob->Nf + g)*4+3] = g3;
1479: return(0);
1480: }
1482: /*@C
1483: PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, dF/du_t, has been set
1485: Not collective
1487: Input Parameter:
1488: . prob - The PetscDS
1490: Output Parameter:
1491: . hasDynJac - flag that pointwise function for dynamic Jacobian has been set
1493: Level: intermediate
1495: .seealso: PetscDSGetDynamicJacobian(), PetscDSSetDynamicJacobian(), PetscDSGetJacobian()
1496: @*/
1497: PetscErrorCode PetscDSHasDynamicJacobian(PetscDS prob, PetscBool *hasDynJac)
1498: {
1499: PetscInt f, g, h;
1503: *hasDynJac = PETSC_FALSE;
1504: for (f = 0; f < prob->Nf; ++f) {
1505: for (g = 0; g < prob->Nf; ++g) {
1506: for (h = 0; h < 4; ++h) {
1507: if (prob->gt[(f*prob->Nf + g)*4+h]) *hasDynJac = PETSC_TRUE;
1508: }
1509: }
1510: }
1511: return(0);
1512: }
1514: /*@C
1515: PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, dF/du_t, function for given test and basis field
1517: Not collective
1519: Input Parameters:
1520: + prob - The PetscDS
1521: . f - The test field number
1522: - g - The field number
1524: Output Parameters:
1525: + g0 - integrand for the test and basis function term
1526: . g1 - integrand for the test function and basis function gradient term
1527: . g2 - integrand for the test function gradient and basis function term
1528: - g3 - integrand for the test function gradient and basis function gradient term
1530: Note: We are using a first order FEM model for the weak form:
1532: \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
1534: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1536: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1537: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1538: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1539: $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1541: + dim - the spatial dimension
1542: . Nf - the number of fields
1543: . uOff - the offset into u[] and u_t[] for each field
1544: . uOff_x - the offset into u_x[] for each field
1545: . u - each field evaluated at the current point
1546: . u_t - the time derivative of each field evaluated at the current point
1547: . u_x - the gradient of each field evaluated at the current point
1548: . aOff - the offset into a[] and a_t[] for each auxiliary field
1549: . aOff_x - the offset into a_x[] for each auxiliary field
1550: . a - each auxiliary field evaluated at the current point
1551: . a_t - the time derivative of each auxiliary field evaluated at the current point
1552: . a_x - the gradient of auxiliary each field evaluated at the current point
1553: . t - current time
1554: . u_tShift - the multiplier a for dF/dU_t
1555: . x - coordinates of the current point
1556: . numConstants - number of constant parameters
1557: . constants - constant parameters
1558: - g0 - output values at the current point
1560: Level: intermediate
1562: .seealso: PetscDSSetJacobian()
1563: @*/
1564: PetscErrorCode PetscDSGetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g,
1565: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1566: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1567: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1568: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1569: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1570: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1571: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1572: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1573: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1574: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1575: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1576: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1577: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1578: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1579: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1580: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1581: {
1584: 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);
1585: 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);
1590: return(0);
1591: }
1593: /*@C
1594: PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, dF/du_t, function for given test and basis fields
1596: Not collective
1598: Input Parameters:
1599: + prob - The PetscDS
1600: . f - The test field number
1601: . g - The field number
1602: . g0 - integrand for the test and basis function term
1603: . g1 - integrand for the test function and basis function gradient term
1604: . g2 - integrand for the test function gradient and basis function term
1605: - g3 - integrand for the test function gradient and basis function gradient term
1607: Note: We are using a first order FEM model for the weak form:
1609: \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
1611: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1613: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1614: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1615: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1616: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1618: + dim - the spatial dimension
1619: . Nf - the number of fields
1620: . uOff - the offset into u[] and u_t[] for each field
1621: . uOff_x - the offset into u_x[] for each field
1622: . u - each field evaluated at the current point
1623: . u_t - the time derivative of each field evaluated at the current point
1624: . u_x - the gradient of each field evaluated at the current point
1625: . aOff - the offset into a[] and a_t[] for each auxiliary field
1626: . aOff_x - the offset into a_x[] for each auxiliary field
1627: . a - each auxiliary field evaluated at the current point
1628: . a_t - the time derivative of each auxiliary field evaluated at the current point
1629: . a_x - the gradient of auxiliary each field evaluated at the current point
1630: . t - current time
1631: . u_tShift - the multiplier a for dF/dU_t
1632: . x - coordinates of the current point
1633: . numConstants - number of constant parameters
1634: . constants - constant parameters
1635: - g0 - output values at the current point
1637: Level: intermediate
1639: .seealso: PetscDSGetJacobian()
1640: @*/
1641: PetscErrorCode PetscDSSetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g,
1642: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1643: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1644: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1645: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1646: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1647: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1648: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1649: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1650: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1651: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1652: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1653: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1654: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1655: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1656: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1657: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1658: {
1667: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1668: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1669: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1670: prob->gt[(f*prob->Nf + g)*4+0] = g0;
1671: prob->gt[(f*prob->Nf + g)*4+1] = g1;
1672: prob->gt[(f*prob->Nf + g)*4+2] = g2;
1673: prob->gt[(f*prob->Nf + g)*4+3] = g3;
1674: return(0);
1675: }
1677: /*@C
1678: PetscDSGetRiemannSolver - Returns the Riemann solver for the given field
1680: Not collective
1682: Input Arguments:
1683: + prob - The PetscDS object
1684: - f - The field number
1686: Output Argument:
1687: . r - Riemann solver
1689: Calling sequence for r:
1691: $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)
1693: + dim - The spatial dimension
1694: . Nf - The number of fields
1695: . x - The coordinates at a point on the interface
1696: . n - The normal vector to the interface
1697: . uL - The state vector to the left of the interface
1698: . uR - The state vector to the right of the interface
1699: . flux - output array of flux through the interface
1700: . numConstants - number of constant parameters
1701: . constants - constant parameters
1702: - ctx - optional user context
1704: Level: intermediate
1706: .seealso: PetscDSSetRiemannSolver()
1707: @*/
1708: PetscErrorCode PetscDSGetRiemannSolver(PetscDS prob, PetscInt f,
1709: 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))
1710: {
1713: 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);
1715: *r = prob->r[f];
1716: return(0);
1717: }
1719: /*@C
1720: PetscDSSetRiemannSolver - Sets the Riemann solver for the given field
1722: Not collective
1724: Input Arguments:
1725: + prob - The PetscDS object
1726: . f - The field number
1727: - r - Riemann solver
1729: Calling sequence for r:
1731: $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)
1733: + dim - The spatial dimension
1734: . Nf - The number of fields
1735: . x - The coordinates at a point on the interface
1736: . n - The normal vector to the interface
1737: . uL - The state vector to the left of the interface
1738: . uR - The state vector to the right of the interface
1739: . flux - output array of flux through the interface
1740: . numConstants - number of constant parameters
1741: . constants - constant parameters
1742: - ctx - optional user context
1744: Level: intermediate
1746: .seealso: PetscDSGetRiemannSolver()
1747: @*/
1748: PetscErrorCode PetscDSSetRiemannSolver(PetscDS prob, PetscInt f,
1749: 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))
1750: {
1756: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1757: PetscDSEnlarge_Static(prob, f+1);
1758: prob->r[f] = r;
1759: return(0);
1760: }
1762: /*@C
1763: PetscDSGetUpdate - Get the pointwise update function for a given field
1765: Not collective
1767: Input Parameters:
1768: + prob - The PetscDS
1769: - f - The field number
1771: Output Parameters:
1772: + update - update function
1774: Note: The calling sequence for the callback update is given by:
1776: $ update(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1777: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1778: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1779: $ PetscReal t, const PetscReal x[], PetscScalar uNew[])
1781: + dim - the spatial dimension
1782: . Nf - the number of fields
1783: . uOff - the offset into u[] and u_t[] for each field
1784: . uOff_x - the offset into u_x[] for each field
1785: . u - each field evaluated at the current point
1786: . u_t - the time derivative of each field evaluated at the current point
1787: . u_x - the gradient of each field evaluated at the current point
1788: . aOff - the offset into a[] and a_t[] for each auxiliary field
1789: . aOff_x - the offset into a_x[] for each auxiliary field
1790: . a - each auxiliary field evaluated at the current point
1791: . a_t - the time derivative of each auxiliary field evaluated at the current point
1792: . a_x - the gradient of auxiliary each field evaluated at the current point
1793: . t - current time
1794: . x - coordinates of the current point
1795: - uNew - new value for field at the current point
1797: Level: intermediate
1799: .seealso: PetscDSSetUpdate(), PetscDSSetResidual()
1800: @*/
1801: PetscErrorCode PetscDSGetUpdate(PetscDS prob, PetscInt f,
1802: void (**update)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1803: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1804: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1805: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
1806: {
1809: 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);
1811: return(0);
1812: }
1814: /*@C
1815: PetscDSSetUpdate - Set the pointwise update function for a given field
1817: Not collective
1819: Input Parameters:
1820: + prob - The PetscDS
1821: . f - The field number
1822: - update - update function
1824: Note: The calling sequence for the callback update is given by:
1826: $ update(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1827: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1828: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1829: $ PetscReal t, const PetscReal x[], PetscScalar uNew[])
1831: + dim - the spatial dimension
1832: . Nf - the number of fields
1833: . uOff - the offset into u[] and u_t[] for each field
1834: . uOff_x - the offset into u_x[] for each field
1835: . u - each field evaluated at the current point
1836: . u_t - the time derivative of each field evaluated at the current point
1837: . u_x - the gradient of each field evaluated at the current point
1838: . aOff - the offset into a[] and a_t[] for each auxiliary field
1839: . aOff_x - the offset into a_x[] for each auxiliary field
1840: . a - each auxiliary field evaluated at the current point
1841: . a_t - the time derivative of each auxiliary field evaluated at the current point
1842: . a_x - the gradient of auxiliary each field evaluated at the current point
1843: . t - current time
1844: . x - coordinates of the current point
1845: - uNew - new field values at the current point
1847: Level: intermediate
1849: .seealso: PetscDSGetResidual()
1850: @*/
1851: PetscErrorCode PetscDSSetUpdate(PetscDS prob, PetscInt f,
1852: void (*update)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1853: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1854: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1855: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
1856: {
1862: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1863: PetscDSEnlarge_Static(prob, f+1);
1864: prob->update[f] = update;
1865: return(0);
1866: }
1868: PetscErrorCode PetscDSGetContext(PetscDS prob, PetscInt f, void **ctx)
1869: {
1872: 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);
1874: *ctx = prob->ctx[f];
1875: return(0);
1876: }
1878: PetscErrorCode PetscDSSetContext(PetscDS prob, PetscInt f, void *ctx)
1879: {
1884: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1885: PetscDSEnlarge_Static(prob, f+1);
1886: prob->ctx[f] = ctx;
1887: return(0);
1888: }
1890: /*@C
1891: PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field
1893: Not collective
1895: Input Parameters:
1896: + prob - The PetscDS
1897: - f - The test field number
1899: Output Parameters:
1900: + f0 - boundary integrand for the test function term
1901: - f1 - boundary integrand for the test function gradient term
1903: Note: We are using a first order FEM model for the weak form:
1905: \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
1907: The calling sequence for the callbacks f0 and f1 is given by:
1909: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1910: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1911: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1912: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[])
1914: + dim - the spatial dimension
1915: . Nf - the number of fields
1916: . uOff - the offset into u[] and u_t[] for each field
1917: . uOff_x - the offset into u_x[] for each field
1918: . u - each field evaluated at the current point
1919: . u_t - the time derivative of each field evaluated at the current point
1920: . u_x - the gradient of each field evaluated at the current point
1921: . aOff - the offset into a[] and a_t[] for each auxiliary field
1922: . aOff_x - the offset into a_x[] for each auxiliary field
1923: . a - each auxiliary field evaluated at the current point
1924: . a_t - the time derivative of each auxiliary field evaluated at the current point
1925: . a_x - the gradient of auxiliary each field evaluated at the current point
1926: . t - current time
1927: . x - coordinates of the current point
1928: . n - unit normal at the current point
1929: . numConstants - number of constant parameters
1930: . constants - constant parameters
1931: - f0 - output values at the current point
1933: Level: intermediate
1935: .seealso: PetscDSSetBdResidual()
1936: @*/
1937: PetscErrorCode PetscDSGetBdResidual(PetscDS prob, PetscInt f,
1938: void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1939: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1940: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1941: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1942: void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1943: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1944: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1945: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1946: {
1949: 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);
1952: return(0);
1953: }
1955: /*@C
1956: PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field
1958: Not collective
1960: Input Parameters:
1961: + prob - The PetscDS
1962: . f - The test field number
1963: . f0 - boundary integrand for the test function term
1964: - f1 - boundary integrand for the test function gradient term
1966: Note: We are using a first order FEM model for the weak form:
1968: \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
1970: The calling sequence for the callbacks f0 and f1 is given by:
1972: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1973: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1974: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1975: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[])
1977: + dim - the spatial dimension
1978: . Nf - the number of fields
1979: . uOff - the offset into u[] and u_t[] for each field
1980: . uOff_x - the offset into u_x[] for each field
1981: . u - each field evaluated at the current point
1982: . u_t - the time derivative of each field evaluated at the current point
1983: . u_x - the gradient of each field evaluated at the current point
1984: . aOff - the offset into a[] and a_t[] for each auxiliary field
1985: . aOff_x - the offset into a_x[] for each auxiliary field
1986: . a - each auxiliary field evaluated at the current point
1987: . a_t - the time derivative of each auxiliary field evaluated at the current point
1988: . a_x - the gradient of auxiliary each field evaluated at the current point
1989: . t - current time
1990: . x - coordinates of the current point
1991: . n - unit normal at the current point
1992: . numConstants - number of constant parameters
1993: . constants - constant parameters
1994: - f0 - output values at the current point
1996: Level: intermediate
1998: .seealso: PetscDSGetBdResidual()
1999: @*/
2000: PetscErrorCode PetscDSSetBdResidual(PetscDS prob, PetscInt f,
2001: void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2002: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2003: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2004: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
2005: void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2006: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2007: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2008: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
2009: {
2014: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2015: PetscDSEnlarge_Static(prob, f+1);
2018: return(0);
2019: }
2021: /*@C
2022: PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field
2024: Not collective
2026: Input Parameters:
2027: + prob - The PetscDS
2028: . f - The test field number
2029: - g - The field number
2031: Output Parameters:
2032: + g0 - integrand for the test and basis function term
2033: . g1 - integrand for the test function and basis function gradient term
2034: . g2 - integrand for the test function gradient and basis function term
2035: - g3 - integrand for the test function gradient and basis function gradient term
2037: Note: We are using a first order FEM model for the weak form:
2039: \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
2041: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
2043: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2044: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2045: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2046: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[])
2048: + dim - the spatial dimension
2049: . Nf - the number of fields
2050: . uOff - the offset into u[] and u_t[] for each field
2051: . uOff_x - the offset into u_x[] for each field
2052: . u - each field evaluated at the current point
2053: . u_t - the time derivative of each field evaluated at the current point
2054: . u_x - the gradient of each field evaluated at the current point
2055: . aOff - the offset into a[] and a_t[] for each auxiliary field
2056: . aOff_x - the offset into a_x[] for each auxiliary field
2057: . a - each auxiliary field evaluated at the current point
2058: . a_t - the time derivative of each auxiliary field evaluated at the current point
2059: . a_x - the gradient of auxiliary each field evaluated at the current point
2060: . t - current time
2061: . u_tShift - the multiplier a for dF/dU_t
2062: . x - coordinates of the current point
2063: . n - normal at the current point
2064: . numConstants - number of constant parameters
2065: . constants - constant parameters
2066: - g0 - output values at the current point
2068: Level: intermediate
2070: .seealso: PetscDSSetBdJacobian()
2071: @*/
2072: PetscErrorCode PetscDSGetBdJacobian(PetscDS prob, PetscInt f, PetscInt g,
2073: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2074: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2075: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2076: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
2077: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2078: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2079: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2080: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
2081: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2082: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2083: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2084: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
2085: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2086: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2087: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2088: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2089: {
2092: 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);
2093: 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);
2098: return(0);
2099: }
2101: /*@C
2102: PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field
2104: Not collective
2106: Input Parameters:
2107: + prob - The PetscDS
2108: . f - The test field number
2109: . g - The field number
2110: . g0 - integrand for the test and basis function term
2111: . g1 - integrand for the test function and basis function gradient term
2112: . g2 - integrand for the test function gradient and basis function term
2113: - g3 - integrand for the test function gradient and basis function gradient term
2115: Note: We are using a first order FEM model for the weak form:
2117: \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
2119: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
2121: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2122: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2123: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2124: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[])
2126: + dim - the spatial dimension
2127: . Nf - the number of fields
2128: . uOff - the offset into u[] and u_t[] for each field
2129: . uOff_x - the offset into u_x[] for each field
2130: . u - each field evaluated at the current point
2131: . u_t - the time derivative of each field evaluated at the current point
2132: . u_x - the gradient of each field evaluated at the current point
2133: . aOff - the offset into a[] and a_t[] for each auxiliary field
2134: . aOff_x - the offset into a_x[] for each auxiliary field
2135: . a - each auxiliary field evaluated at the current point
2136: . a_t - the time derivative of each auxiliary field evaluated at the current point
2137: . a_x - the gradient of auxiliary each field evaluated at the current point
2138: . t - current time
2139: . u_tShift - the multiplier a for dF/dU_t
2140: . x - coordinates of the current point
2141: . n - normal at the current point
2142: . numConstants - number of constant parameters
2143: . constants - constant parameters
2144: - g0 - output values at the current point
2146: Level: intermediate
2148: .seealso: PetscDSGetBdJacobian()
2149: @*/
2150: PetscErrorCode PetscDSSetBdJacobian(PetscDS prob, PetscInt f, PetscInt g,
2151: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2152: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2153: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2154: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
2155: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2156: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2157: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2158: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
2159: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2160: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2161: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2162: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
2163: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2164: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2165: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2166: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2167: {
2176: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2177: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
2178: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
2179: prob->gBd[(f*prob->Nf + g)*4+0] = g0;
2180: prob->gBd[(f*prob->Nf + g)*4+1] = g1;
2181: prob->gBd[(f*prob->Nf + g)*4+2] = g2;
2182: prob->gBd[(f*prob->Nf + g)*4+3] = g3;
2183: return(0);
2184: }
2186: /*@C
2187: PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field
2189: Not collective
2191: Input Parameters:
2192: + prob - The PetscDS
2193: - f - The test field number
2195: Output Parameter:
2196: . exactSol - exact solution for the test field
2198: Note: The calling sequence for the solution functions is given by:
2200: $ sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2202: + dim - the spatial dimension
2203: . t - current time
2204: . x - coordinates of the current point
2205: . Nc - the number of field components
2206: . u - the solution field evaluated at the current point
2207: - ctx - a user context
2209: Level: intermediate
2211: .seealso: PetscDSSetExactSolution()
2212: @*/
2213: PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx))
2214: {
2217: 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);
2219: return(0);
2220: }
2222: /*@C
2223: PetscDSSetExactSolution - Get the pointwise exact solution function for a given test field
2225: Not collective
2227: Input Parameters:
2228: + prob - The PetscDS
2229: . f - The test field number
2230: - sol - solution function for the test fields
2232: Note: The calling sequence for solution functions is given by:
2234: $ sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2236: + dim - the spatial dimension
2237: . t - current time
2238: . x - coordinates of the current point
2239: . Nc - the number of field components
2240: . u - the solution field evaluated at the current point
2241: - ctx - a user context
2243: Level: intermediate
2245: .seealso: PetscDSGetExactSolution()
2246: @*/
2247: PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx))
2248: {
2253: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2254: PetscDSEnlarge_Static(prob, f+1);
2256: return(0);
2257: }
2259: /*@C
2260: PetscDSGetConstants - Returns the array of constants passed to point functions
2262: Not collective
2264: Input Parameter:
2265: . prob - The PetscDS object
2267: Output Parameters:
2268: + numConstants - The number of constants
2269: - constants - The array of constants, NULL if there are none
2271: Level: intermediate
2273: .seealso: PetscDSSetConstants(), PetscDSCreate()
2274: @*/
2275: PetscErrorCode PetscDSGetConstants(PetscDS prob, PetscInt *numConstants, const PetscScalar *constants[])
2276: {
2281: return(0);
2282: }
2284: /*@C
2285: PetscDSSetConstants - Set the array of constants passed to point functions
2287: Not collective
2289: Input Parameters:
2290: + prob - The PetscDS object
2291: . numConstants - The number of constants
2292: - constants - The array of constants, NULL if there are none
2294: Level: intermediate
2296: .seealso: PetscDSGetConstants(), PetscDSCreate()
2297: @*/
2298: PetscErrorCode PetscDSSetConstants(PetscDS prob, PetscInt numConstants, PetscScalar constants[])
2299: {
2304: if (numConstants != prob->numConstants) {
2305: PetscFree(prob->constants);
2306: prob->numConstants = numConstants;
2307: if (prob->numConstants) {
2308: PetscMalloc1(prob->numConstants, &prob->constants);
2309: } else {
2310: prob->constants = NULL;
2311: }
2312: }
2313: if (prob->numConstants) {
2315: PetscMemcpy(prob->constants, constants, prob->numConstants * sizeof(PetscScalar));
2316: }
2317: return(0);
2318: }
2320: /*@
2321: PetscDSGetFieldIndex - Returns the index of the given field
2323: Not collective
2325: Input Parameters:
2326: + prob - The PetscDS object
2327: - disc - The discretization object
2329: Output Parameter:
2330: . f - The field number
2332: Level: beginner
2334: .seealso: PetscGetDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
2335: @*/
2336: PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f)
2337: {
2338: PetscInt g;
2343: *f = -1;
2344: for (g = 0; g < prob->Nf; ++g) {if (disc == prob->disc[g]) break;}
2345: if (g == prob->Nf) SETERRQ(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS.");
2346: *f = g;
2347: return(0);
2348: }
2350: /*@
2351: PetscDSGetFieldSize - Returns the size of the given field in the full space basis
2353: Not collective
2355: Input Parameters:
2356: + prob - The PetscDS object
2357: - f - The field number
2359: Output Parameter:
2360: . size - The size
2362: Level: beginner
2364: .seealso: PetscDSGetFieldOffset(), PetscDSGetNumFields(), PetscDSCreate()
2365: @*/
2366: PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size)
2367: {
2373: 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);
2374: PetscDSSetUp(prob);
2375: *size = prob->Nb[f];
2376: return(0);
2377: }
2379: /*@
2380: PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis
2382: Not collective
2384: Input Parameters:
2385: + prob - The PetscDS object
2386: - f - The field number
2388: Output Parameter:
2389: . off - The offset
2391: Level: beginner
2393: .seealso: PetscDSGetFieldSize(), PetscDSGetNumFields(), PetscDSCreate()
2394: @*/
2395: PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off)
2396: {
2397: PetscInt size, g;
2403: 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);
2404: *off = 0;
2405: for (g = 0; g < f; ++g) {
2406: PetscDSGetFieldSize(prob, g, &size);
2407: *off += size;
2408: }
2409: return(0);
2410: }
2412: /*@
2413: PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point
2415: Not collective
2417: Input Parameter:
2418: . prob - The PetscDS object
2420: Output Parameter:
2421: . dimensions - The number of dimensions
2423: Level: beginner
2425: .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate()
2426: @*/
2427: PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[])
2428: {
2433: PetscDSSetUp(prob);
2435: *dimensions = prob->Nb;
2436: return(0);
2437: }
2439: /*@
2440: PetscDSGetComponents - Returns the number of components for each field on an evaluation point
2442: Not collective
2444: Input Parameter:
2445: . prob - The PetscDS object
2447: Output Parameter:
2448: . components - The number of components
2450: Level: beginner
2452: .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate()
2453: @*/
2454: PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[])
2455: {
2460: PetscDSSetUp(prob);
2462: *components = prob->Nc;
2463: return(0);
2464: }
2466: /*@
2467: PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point
2469: Not collective
2471: Input Parameters:
2472: + prob - The PetscDS object
2473: - f - The field number
2475: Output Parameter:
2476: . off - The offset
2478: Level: beginner
2480: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2481: @*/
2482: PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off)
2483: {
2487: 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);
2488: *off = prob->off[f];
2489: return(0);
2490: }
2492: /*@
2493: PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point
2495: Not collective
2497: Input Parameter:
2498: . prob - The PetscDS object
2500: Output Parameter:
2501: . offsets - The offsets
2503: Level: beginner
2505: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2506: @*/
2507: PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[])
2508: {
2512: *offsets = prob->off;
2513: return(0);
2514: }
2516: /*@
2517: PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point
2519: Not collective
2521: Input Parameter:
2522: . prob - The PetscDS object
2524: Output Parameter:
2525: . offsets - The offsets
2527: Level: beginner
2529: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2530: @*/
2531: PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[])
2532: {
2536: *offsets = prob->offDer;
2537: return(0);
2538: }
2540: /*@C
2541: PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization
2543: Not collective
2545: Input Parameter:
2546: . prob - The PetscDS object
2548: Output Parameters:
2549: + basis - The basis function tabulation at quadrature points
2550: - basisDer - The basis function derivative tabulation at quadrature points
2552: Level: intermediate
2554: .seealso: PetscDSCreate()
2555: @*/
2556: PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer)
2557: {
2562: PetscDSSetUp(prob);
2565: return(0);
2566: }
2568: /*@C
2569: PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces
2571: Not collective
2573: Input Parameter:
2574: . prob - The PetscDS object
2576: Output Parameters:
2577: + basisFace - The basis function tabulation at quadrature points
2578: - basisDerFace - The basis function derivative tabulation at quadrature points
2580: Level: intermediate
2582: .seealso: PetscDSGetTabulation(), PetscDSCreate()
2583: @*/
2584: PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer)
2585: {
2590: PetscDSSetUp(prob);
2593: return(0);
2594: }
2596: PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x)
2597: {
2602: PetscDSSetUp(prob);
2606: return(0);
2607: }
2609: PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3)
2610: {
2615: PetscDSSetUp(prob);
2622: return(0);
2623: }
2625: PetscErrorCode PetscDSGetRefCoordArrays(PetscDS prob, PetscReal **x, PetscScalar **refSpaceDer)
2626: {
2631: PetscDSSetUp(prob);
2634: return(0);
2635: }
2637: /*@C
2638: PetscDSAddBoundary - Add a boundary condition to the model
2640: Input Parameters:
2641: + ds - The PetscDS object
2642: . type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2643: . name - The BC name
2644: . labelname - The label defining constrained points
2645: . field - The field to constrain
2646: . numcomps - The number of constrained field components
2647: . comps - An array of constrained component numbers
2648: . bcFunc - A pointwise function giving boundary values
2649: . numids - The number of DMLabel ids for constrained points
2650: . ids - An array of ids for constrained points
2651: - ctx - An optional user context for bcFunc
2653: Options Database Keys:
2654: + -bc_<boundary name> <num> - Overrides the boundary ids
2655: - -bc_<boundary name>_comp <num> - Overrides the boundary components
2657: Level: developer
2659: .seealso: PetscDSGetBoundary()
2660: @*/
2661: 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)
2662: {
2663: DSBoundary b;
2668: PetscNew(&b);
2669: PetscStrallocpy(name, (char **) &b->name);
2670: PetscStrallocpy(labelname, (char **) &b->labelname);
2671: PetscMalloc1(numcomps, &b->comps);
2672: if (numcomps) {PetscMemcpy(b->comps, comps, numcomps*sizeof(PetscInt));}
2673: PetscMalloc1(numids, &b->ids);
2674: if (numids) {PetscMemcpy(b->ids, ids, numids*sizeof(PetscInt));}
2675: b->type = type;
2676: b->field = field;
2677: b->numcomps = numcomps;
2678: b->func = bcFunc;
2679: b->numids = numids;
2680: b->ctx = ctx;
2681: b->next = ds->boundary;
2682: ds->boundary = b;
2683: return(0);
2684: }
2686: /*@
2687: PetscDSGetNumBoundary - Get the number of registered BC
2689: Input Parameters:
2690: . ds - The PetscDS object
2692: Output Parameters:
2693: . numBd - The number of BC
2695: Level: intermediate
2697: .seealso: PetscDSAddBoundary(), PetscDSGetBoundary()
2698: @*/
2699: PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd)
2700: {
2701: DSBoundary b = ds->boundary;
2706: *numBd = 0;
2707: while (b) {++(*numBd); b = b->next;}
2708: return(0);
2709: }
2711: /*@C
2712: PetscDSGetBoundary - Add a boundary condition to the model
2714: Input Parameters:
2715: + ds - The PetscDS object
2716: - bd - The BC number
2718: Output Parameters:
2719: + type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2720: . name - The BC name
2721: . labelname - The label defining constrained points
2722: . field - The field to constrain
2723: . numcomps - The number of constrained field components
2724: . comps - An array of constrained component numbers
2725: . bcFunc - A pointwise function giving boundary values
2726: . numids - The number of DMLabel ids for constrained points
2727: . ids - An array of ids for constrained points
2728: - ctx - An optional user context for bcFunc
2730: Options Database Keys:
2731: + -bc_<boundary name> <num> - Overrides the boundary ids
2732: - -bc_<boundary name>_comp <num> - Overrides the boundary components
2734: Level: developer
2736: .seealso: PetscDSAddBoundary()
2737: @*/
2738: 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)
2739: {
2740: DSBoundary b = ds->boundary;
2741: PetscInt n = 0;
2745: while (b) {
2746: if (n == bd) break;
2747: b = b->next;
2748: ++n;
2749: }
2750: if (!b) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %d is not in [0, %d)", bd, n);
2751: if (type) {
2753: *type = b->type;
2754: }
2755: if (name) {
2757: *name = b->name;
2758: }
2759: if (labelname) {
2761: *labelname = b->labelname;
2762: }
2763: if (field) {
2765: *field = b->field;
2766: }
2767: if (numcomps) {
2769: *numcomps = b->numcomps;
2770: }
2771: if (comps) {
2773: *comps = b->comps;
2774: }
2775: if (func) {
2777: *func = b->func;
2778: }
2779: if (numids) {
2781: *numids = b->numids;
2782: }
2783: if (ids) {
2785: *ids = b->ids;
2786: }
2787: if (ctx) {
2789: *ctx = b->ctx;
2790: }
2791: return(0);
2792: }
2794: /*@
2795: PetscDSCopyBoundary - Copy all boundary condition objects to the new problem
2797: Not collective
2799: Input Parameter:
2800: . prob - The PetscDS object
2802: Output Parameter:
2803: . newprob - The PetscDS copy
2805: Level: intermediate
2807: .seealso: PetscDSCopyEquations(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
2808: @*/
2809: PetscErrorCode PetscDSCopyBoundary(PetscDS probA, PetscDS probB)
2810: {
2811: DSBoundary b, next, *lastnext;
2817: if (probA == probB) return(0);
2818: next = probB->boundary;
2819: while (next) {
2820: DSBoundary b = next;
2822: next = b->next;
2823: PetscFree(b->comps);
2824: PetscFree(b->ids);
2825: PetscFree(b->name);
2826: PetscFree(b->labelname);
2827: PetscFree(b);
2828: }
2829: lastnext = &(probB->boundary);
2830: for (b = probA->boundary; b; b = b->next) {
2831: DSBoundary bNew;
2833: PetscNew(&bNew);
2834: bNew->numcomps = b->numcomps;
2835: PetscMalloc1(bNew->numcomps, &bNew->comps);
2836: PetscMemcpy(bNew->comps, b->comps, bNew->numcomps*sizeof(PetscInt));
2837: bNew->numids = b->numids;
2838: PetscMalloc1(bNew->numids, &bNew->ids);
2839: PetscMemcpy(bNew->ids, b->ids, bNew->numids*sizeof(PetscInt));
2840: PetscStrallocpy(b->labelname,(char **) &(bNew->labelname));
2841: PetscStrallocpy(b->name,(char **) &(bNew->name));
2842: bNew->ctx = b->ctx;
2843: bNew->type = b->type;
2844: bNew->field = b->field;
2845: bNew->func = b->func;
2847: *lastnext = bNew;
2848: lastnext = &(bNew->next);
2849: }
2850: return(0);
2851: }
2853: /*@C
2854: PetscDSSelectEquations - Copy pointwise function pointers to the new problem with different field layout
2856: Not collective
2858: Input Parameter:
2859: + prob - The PetscDS object
2860: . numFields - Number of new fields
2861: - fields - Old field number for each new field
2863: Output Parameter:
2864: . newprob - The PetscDS copy
2866: Level: intermediate
2868: .seealso: PetscDSCopyBoundary(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
2869: @*/
2870: PetscErrorCode PetscDSSelectEquations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
2871: {
2872: PetscInt Nf, Nfn, fn, gn;
2879: PetscDSGetNumFields(prob, &Nf);
2880: PetscDSGetNumFields(newprob, &Nfn);
2881: 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);
2882: for (fn = 0; fn < numFields; ++fn) {
2883: const PetscInt f = fields ? fields[fn] : fn;
2884: PetscPointFunc obj;
2885: PetscPointFunc f0, f1;
2886: PetscBdPointFunc f0Bd, f1Bd;
2887: PetscRiemannFunc r;
2889: if (f >= Nf) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Field %D must be in [0, %D)", f, Nf);
2890: PetscDSGetObjective(prob, f, &obj);
2891: PetscDSGetResidual(prob, f, &f0, &f1);
2892: PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd);
2893: PetscDSGetRiemannSolver(prob, f, &r);
2894: PetscDSSetObjective(newprob, fn, obj);
2895: PetscDSSetResidual(newprob, fn, f0, f1);
2896: PetscDSSetBdResidual(newprob, fn, f0Bd, f1Bd);
2897: PetscDSSetRiemannSolver(newprob, fn, r);
2898: for (gn = 0; gn < numFields; ++gn) {
2899: const PetscInt g = fields ? fields[gn] : gn;
2900: PetscPointJac g0, g1, g2, g3;
2901: PetscPointJac g0p, g1p, g2p, g3p;
2902: PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd;
2904: if (g >= Nf) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Field %D must be in [0, %D)", g, Nf);
2905: PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3);
2906: PetscDSGetJacobianPreconditioner(prob, f, g, &g0p, &g1p, &g2p, &g3p);
2907: PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd);
2908: PetscDSSetJacobian(newprob, fn, gn, g0, g1, g2, g3);
2909: PetscDSSetJacobianPreconditioner(prob, fn, gn, g0p, g1p, g2p, g3p);
2910: PetscDSSetBdJacobian(newprob, fn, gn, g0Bd, g1Bd, g2Bd, g3Bd);
2911: }
2912: }
2913: return(0);
2914: }
2916: /*@
2917: PetscDSCopyEquations - Copy all pointwise function pointers to the new problem
2919: Not collective
2921: Input Parameter:
2922: . prob - The PetscDS object
2924: Output Parameter:
2925: . newprob - The PetscDS copy
2927: Level: intermediate
2929: .seealso: PetscDSCopyBoundary(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
2930: @*/
2931: PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob)
2932: {
2933: PetscInt Nf, Ng;
2939: PetscDSGetNumFields(prob, &Nf);
2940: PetscDSGetNumFields(newprob, &Ng);
2941: if (Nf != Ng) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %D != %D", Nf, Ng);
2942: PetscDSSelectEquations(prob, Nf, NULL, newprob);
2943: return(0);
2944: }
2945: /*@
2946: PetscDSCopyConstants - Copy all constants to the new problem
2948: Not collective
2950: Input Parameter:
2951: . prob - The PetscDS object
2953: Output Parameter:
2954: . newprob - The PetscDS copy
2956: Level: intermediate
2958: .seealso: PetscDSCopyBoundary(), PetscDSCopyEquations(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
2959: @*/
2960: PetscErrorCode PetscDSCopyConstants(PetscDS prob, PetscDS newprob)
2961: {
2962: PetscInt Nc;
2963: const PetscScalar *constants;
2964: PetscErrorCode ierr;
2969: PetscDSGetConstants(prob, &Nc, &constants);
2970: PetscDSSetConstants(newprob, Nc, (PetscScalar *) constants);
2971: return(0);
2972: }
2974: PetscErrorCode PetscDSGetHeightSubspace(PetscDS prob, PetscInt height, PetscDS *subprob)
2975: {
2976: PetscInt dim, Nf, f;
2982: if (height == 0) {*subprob = prob; return(0);}
2983: PetscDSGetNumFields(prob, &Nf);
2984: PetscDSGetSpatialDimension(prob, &dim);
2985: if (height > dim) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_OUTOFRANGE, "DS can only handle height in [0, %D], not %D", dim, height);
2986: if (!prob->subprobs) {PetscCalloc1(dim, &prob->subprobs);}
2987: if (!prob->subprobs[height-1]) {
2988: PetscInt cdim;
2990: PetscDSCreate(PetscObjectComm((PetscObject) prob), &prob->subprobs[height-1]);
2991: PetscDSGetCoordinateDimension(prob, &cdim);
2992: PetscDSSetCoordinateDimension(prob->subprobs[height-1], cdim);
2993: for (f = 0; f < Nf; ++f) {
2994: PetscFE subfe;
2995: PetscObject obj;
2996: PetscClassId id;
2998: PetscDSGetDiscretization(prob, f, &obj);
2999: PetscObjectGetClassId(obj, &id);
3000: if (id == PETSCFE_CLASSID) {PetscFEGetHeightSubspace((PetscFE) obj, height, &subfe);}
3001: else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unsupported discretization type for field %d", f);
3002: PetscDSSetDiscretization(prob->subprobs[height-1], f, (PetscObject) subfe);
3003: }
3004: }
3005: *subprob = prob->subprobs[height-1];
3006: return(0);
3007: }
3009: static PetscErrorCode PetscDSDestroy_Basic(PetscDS prob)
3010: {
3011: PetscErrorCode ierr;
3014: PetscFree(prob->data);
3015: return(0);
3016: }
3018: static PetscErrorCode PetscDSInitialize_Basic(PetscDS prob)
3019: {
3021: prob->ops->setfromoptions = NULL;
3022: prob->ops->setup = NULL;
3023: prob->ops->view = NULL;
3024: prob->ops->destroy = PetscDSDestroy_Basic;
3025: return(0);
3026: }
3028: /*MC
3029: PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions
3031: Level: intermediate
3033: .seealso: PetscDSType, PetscDSCreate(), PetscDSSetType()
3034: M*/
3036: PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS prob)
3037: {
3038: PetscDS_Basic *b;
3039: PetscErrorCode ierr;
3043: PetscNewLog(prob, &b);
3044: prob->data = b;
3046: PetscDSInitialize_Basic(prob);
3047: return(0);
3048: }