Actual source code: dtds.c
petsc-3.8.4 2018-03-24
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: .keywords: PetscDS, register
38: .seealso: PetscDSRegisterAll(), PetscDSRegisterDestroy()
40: @*/
41: PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS))
42: {
46: PetscFunctionListAdd(&PetscDSList, sname, function);
47: return(0);
48: }
50: /*@C
51: PetscDSSetType - Builds a particular PetscDS
53: Collective on PetscDS
55: Input Parameters:
56: + prob - The PetscDS object
57: - name - The kind of system
59: Options Database Key:
60: . -petscds_type <type> - Sets the PetscDS type; use -help for a list of available types
62: Level: intermediate
64: .keywords: PetscDS, set, type
65: .seealso: PetscDSGetType(), PetscDSCreate()
66: @*/
67: PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name)
68: {
69: PetscErrorCode (*r)(PetscDS);
70: PetscBool match;
75: PetscObjectTypeCompare((PetscObject) prob, name, &match);
76: if (match) return(0);
78: PetscDSRegisterAll();
79: PetscFunctionListFind(PetscDSList, name, &r);
80: if (!r) SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscDS type: %s", name);
82: if (prob->ops->destroy) {
83: (*prob->ops->destroy)(prob);
84: prob->ops->destroy = NULL;
85: }
86: (*r)(prob);
87: PetscObjectChangeTypeName((PetscObject) prob, name);
88: return(0);
89: }
91: /*@C
92: PetscDSGetType - Gets the PetscDS type name (as a string) from the object.
94: Not Collective
96: Input Parameter:
97: . prob - The PetscDS
99: Output Parameter:
100: . name - The PetscDS type name
102: Level: intermediate
104: .keywords: PetscDS, get, type, name
105: .seealso: PetscDSSetType(), PetscDSCreate()
106: @*/
107: PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name)
108: {
114: PetscDSRegisterAll();
115: *name = ((PetscObject) prob)->type_name;
116: return(0);
117: }
119: static PetscErrorCode PetscDSView_Ascii(PetscDS prob, PetscViewer viewer)
120: {
121: PetscViewerFormat format;
122: const PetscScalar *constants;
123: PetscInt numConstants, f;
124: PetscErrorCode ierr;
127: PetscViewerGetFormat(viewer, &format);
128: PetscViewerASCIIPrintf(viewer, "Discrete System with %d fields\n", prob->Nf);
129: PetscViewerASCIIPushTab(viewer);
130: for (f = 0; f < prob->Nf; ++f) {
131: PetscObject obj;
132: PetscClassId id;
133: const char *name;
134: PetscInt Nc;
136: PetscDSGetDiscretization(prob, f, &obj);
137: PetscObjectGetClassId(obj, &id);
138: PetscObjectGetName(obj, &name);
139: PetscViewerASCIIPrintf(viewer, "Field %s", name ? name : "<unknown>");
140: if (id == PETSCFE_CLASSID) {
141: PetscFEGetNumComponents((PetscFE) obj, &Nc);
142: PetscViewerASCIIPrintf(viewer, " FEM");
143: } else if (id == PETSCFV_CLASSID) {
144: PetscFVGetNumComponents((PetscFV) obj, &Nc);
145: PetscViewerASCIIPrintf(viewer, " FVM");
146: }
147: else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %D", f);
148: if (Nc > 1) {PetscViewerASCIIPrintf(viewer, "%D components", Nc);}
149: else {PetscViewerASCIIPrintf(viewer, "%D component ", Nc);}
150: if (prob->implicit[f]) {PetscViewerASCIIPrintf(viewer, " (implicit)");}
151: else {PetscViewerASCIIPrintf(viewer, " (explicit)");}
152: if (prob->adjacency[f*2+0]) {
153: if (prob->adjacency[f*2+1]) {PetscViewerASCIIPrintf(viewer, " (adj FVM++)");}
154: else {PetscViewerASCIIPrintf(viewer, " (adj FVM)");}
155: } else {
156: if (prob->adjacency[f*2+1]) {PetscViewerASCIIPrintf(viewer, " (adj FEM)");}
157: else {PetscViewerASCIIPrintf(viewer, " (adj FUNKY)");}
158: }
159: PetscViewerASCIIPrintf(viewer, "\n");
160: if (format == PETSC_VIEWER_ASCII_INFO_DETAIL) {
161: if (id == PETSCFE_CLASSID) {PetscFEView((PetscFE) obj, viewer);}
162: else if (id == PETSCFV_CLASSID) {PetscFVView((PetscFV) obj, viewer);}
163: }
164: }
165: PetscDSGetConstants(prob, &numConstants, &constants);
166: if (numConstants) {
167: PetscViewerASCIIPrintf(viewer, "%D constants\n", numConstants);
168: PetscViewerASCIIPushTab(viewer);
169: for (f = 0; f < numConstants; ++f) {PetscViewerASCIIPrintf(viewer, "%g\n", (double) PetscRealPart(constants[f]));}
170: PetscViewerASCIIPopTab(viewer);
171: }
172: PetscViewerASCIIPopTab(viewer);
173: return(0);
174: }
176: /*@C
177: PetscDSView - Views a PetscDS
179: Collective on PetscDS
181: Input Parameter:
182: + prob - the PetscDS object to view
183: - v - the viewer
185: Level: developer
187: .seealso PetscDSDestroy()
188: @*/
189: PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v)
190: {
191: PetscBool iascii;
196: if (!v) {PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject) prob), &v);}
198: PetscObjectTypeCompare((PetscObject) v, PETSCVIEWERASCII, &iascii);
199: if (iascii) {PetscDSView_Ascii(prob, v);}
200: if (prob->ops->view) {(*prob->ops->view)(prob, v);}
201: return(0);
202: }
204: /*@
205: PetscDSSetFromOptions - sets parameters in a PetscDS from the options database
207: Collective on PetscDS
209: Input Parameter:
210: . prob - the PetscDS object to set options for
212: Options Database:
214: Level: developer
216: .seealso PetscDSView()
217: @*/
218: PetscErrorCode PetscDSSetFromOptions(PetscDS prob)
219: {
220: DSBoundary b;
221: const char *defaultType;
222: char name[256];
223: PetscBool flg;
228: if (!((PetscObject) prob)->type_name) {
229: defaultType = PETSCDSBASIC;
230: } else {
231: defaultType = ((PetscObject) prob)->type_name;
232: }
233: PetscDSRegisterAll();
235: PetscObjectOptionsBegin((PetscObject) prob);
236: for (b = prob->boundary; b; b = b->next) {
237: char optname[1024];
238: PetscInt ids[1024], len = 1024;
239: PetscBool flg;
241: PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name);
242: PetscMemzero(ids, sizeof(ids));
243: PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg);
244: if (flg) {
245: b->numids = len;
246: PetscFree(b->ids);
247: PetscMalloc1(len, &b->ids);
248: PetscMemcpy(b->ids, ids, len*sizeof(PetscInt));
249: }
250: len = 1024;
251: PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name);
252: PetscMemzero(ids, sizeof(ids));
253: PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg);
254: if (flg) {
255: b->numcomps = len;
256: PetscFree(b->comps);
257: PetscMalloc1(len, &b->comps);
258: PetscMemcpy(b->comps, ids, len*sizeof(PetscInt));
259: }
260: }
261: PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg);
262: if (flg) {
263: PetscDSSetType(prob, name);
264: } else if (!((PetscObject) prob)->type_name) {
265: PetscDSSetType(prob, defaultType);
266: }
267: if (prob->ops->setfromoptions) {(*prob->ops->setfromoptions)(prob);}
268: /* process any options handlers added with PetscObjectAddOptionsHandler() */
269: PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject) prob);
270: PetscOptionsEnd();
271: PetscDSViewFromOptions(prob, NULL, "-petscds_view");
272: return(0);
273: }
275: /*@C
276: PetscDSSetUp - Construct data structures for the PetscDS
278: Collective on PetscDS
280: Input Parameter:
281: . prob - the PetscDS object to setup
283: Level: developer
285: .seealso PetscDSView(), PetscDSDestroy()
286: @*/
287: PetscErrorCode PetscDSSetUp(PetscDS prob)
288: {
289: const PetscInt Nf = prob->Nf;
290: PetscInt dim, work, NcMax = 0, NqMax = 0, f;
295: if (prob->setup) return(0);
296: /* Calculate sizes */
297: PetscDSGetSpatialDimension(prob, &dim);
298: prob->totDim = prob->totComp = 0;
299: PetscMalloc2(Nf,&prob->Nc,Nf,&prob->Nb);
300: PetscCalloc2(Nf+1,&prob->off,Nf+1,&prob->offDer);
301: PetscMalloc4(Nf,&prob->basis,Nf,&prob->basisDer,Nf,&prob->basisFace,Nf,&prob->basisDerFace);
302: for (f = 0; f < Nf; ++f) {
303: PetscObject obj;
304: PetscClassId id;
305: PetscQuadrature q;
306: PetscInt Nq = 0, Nb, Nc;
308: PetscDSGetDiscretization(prob, f, &obj);
309: PetscObjectGetClassId(obj, &id);
310: if (id == PETSCFE_CLASSID) {
311: PetscFE fe = (PetscFE) obj;
313: PetscFEGetQuadrature(fe, &q);
314: PetscFEGetDimension(fe, &Nb);
315: PetscFEGetNumComponents(fe, &Nc);
316: PetscFEGetDefaultTabulation(fe, &prob->basis[f], &prob->basisDer[f], NULL);
317: PetscFEGetFaceTabulation(fe, &prob->basisFace[f], &prob->basisDerFace[f], NULL);
318: } else if (id == PETSCFV_CLASSID) {
319: PetscFV fv = (PetscFV) obj;
321: PetscFVGetQuadrature(fv, &q);
322: PetscFVGetNumComponents(fv, &Nc);
323: Nb = Nc;
324: PetscFVGetDefaultTabulation(fv, &prob->basis[f], &prob->basisDer[f], NULL);
325: /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */
326: } else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", f);
327: prob->Nc[f] = Nc;
328: prob->Nb[f] = Nb;
329: prob->off[f+1] = Nc + prob->off[f];
330: prob->offDer[f+1] = Nc*dim + prob->offDer[f];
331: if (q) {PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL);}
332: NqMax = PetscMax(NqMax, Nq);
333: NcMax = PetscMax(NcMax, Nc);
334: prob->totDim += Nb;
335: prob->totComp += Nc;
336: }
337: work = PetscMax(prob->totComp*dim, PetscSqr(NcMax*dim));
338: /* Allocate works space */
339: PetscMalloc5(prob->totComp,&prob->u,prob->totComp,&prob->u_t,prob->totComp*dim,&prob->u_x,dim,&prob->x,work,&prob->refSpaceDer);
340: 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);
341: if (prob->ops->setup) {(*prob->ops->setup)(prob);}
342: prob->setup = PETSC_TRUE;
343: return(0);
344: }
346: static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob)
347: {
351: PetscFree2(prob->Nc,prob->Nb);
352: PetscFree2(prob->off,prob->offDer);
353: PetscFree4(prob->basis,prob->basisDer,prob->basisFace,prob->basisDerFace);
354: PetscFree5(prob->u,prob->u_t,prob->u_x,prob->x,prob->refSpaceDer);
355: PetscFree6(prob->f0,prob->f1,prob->g0,prob->g1,prob->g2,prob->g3);
356: return(0);
357: }
359: static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew)
360: {
361: PetscObject *tmpd;
362: PetscBool *tmpi, *tmpa;
363: PetscPointFunc *tmpobj, *tmpf, *tmpup;
364: PetscPointJac *tmpg, *tmpgp, *tmpgt;
365: PetscBdPointFunc *tmpfbd;
366: PetscBdPointJac *tmpgbd;
367: PetscRiemannFunc *tmpr;
368: PetscSimplePointFunc *tmpexactSol;
369: void **tmpctx;
370: PetscInt Nf = prob->Nf, f, i;
371: PetscErrorCode ierr;
374: if (Nf >= NfNew) return(0);
375: prob->setup = PETSC_FALSE;
376: PetscDSDestroyStructs_Static(prob);
377: PetscMalloc3(NfNew, &tmpd, NfNew, &tmpi, NfNew*2, &tmpa);
378: 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];}
379: 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;}
380: PetscFree3(prob->disc, prob->implicit, prob->adjacency);
381: prob->Nf = NfNew;
382: prob->disc = tmpd;
383: prob->implicit = tmpi;
384: prob->adjacency = tmpa;
385: PetscCalloc7(NfNew, &tmpobj, NfNew*2, &tmpf, NfNew*NfNew*4, &tmpg, NfNew*NfNew*4, &tmpgp, NfNew*NfNew*4, &tmpgt, NfNew, &tmpr, NfNew, &tmpctx);
386: PetscCalloc1(NfNew, &tmpup);
387: for (f = 0; f < Nf; ++f) tmpobj[f] = prob->obj[f];
388: for (f = 0; f < Nf*2; ++f) tmpf[f] = prob->f[f];
389: for (f = 0; f < Nf*Nf*4; ++f) tmpg[f] = prob->g[f];
390: for (f = 0; f < Nf*Nf*4; ++f) tmpgp[f] = prob->gp[f];
391: for (f = 0; f < Nf; ++f) tmpr[f] = prob->r[f];
392: for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f];
393: for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f];
394: for (f = Nf; f < NfNew; ++f) tmpobj[f] = NULL;
395: for (f = Nf*2; f < NfNew*2; ++f) tmpf[f] = NULL;
396: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpg[f] = NULL;
397: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgp[f] = NULL;
398: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgt[f] = NULL;
399: for (f = Nf; f < NfNew; ++f) tmpr[f] = NULL;
400: for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL;
401: for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL;
402: PetscFree7(prob->obj, prob->f, prob->g, prob->gp, prob->gt, prob->r, prob->ctx);
403: PetscFree(prob->update);
404: prob->obj = tmpobj;
405: prob->f = tmpf;
406: prob->g = tmpg;
407: prob->gp = tmpgp;
408: prob->gt = tmpgt;
409: prob->r = tmpr;
410: prob->update = tmpup;
411: prob->ctx = tmpctx;
412: PetscCalloc3(NfNew*2, &tmpfbd, NfNew*NfNew*4, &tmpgbd, NfNew, &tmpexactSol);
413: for (f = 0; f < Nf*2; ++f) tmpfbd[f] = prob->fBd[f];
414: for (f = 0; f < Nf*Nf*4; ++f) tmpgbd[f] = prob->gBd[f];
415: for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f];
416: for (f = Nf*2; f < NfNew*2; ++f) tmpfbd[f] = NULL;
417: for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgbd[f] = NULL;
418: for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL;
419: PetscFree3(prob->fBd, prob->gBd, prob->exactSol);
420: prob->fBd = tmpfbd;
421: prob->gBd = tmpgbd;
422: prob->exactSol = tmpexactSol;
423: return(0);
424: }
426: /*@
427: PetscDSDestroy - Destroys a PetscDS object
429: Collective on PetscDS
431: Input Parameter:
432: . prob - the PetscDS object to destroy
434: Level: developer
436: .seealso PetscDSView()
437: @*/
438: PetscErrorCode PetscDSDestroy(PetscDS *prob)
439: {
440: PetscInt f;
441: DSBoundary next;
445: if (!*prob) return(0);
448: if (--((PetscObject)(*prob))->refct > 0) {*prob = 0; return(0);}
449: ((PetscObject) (*prob))->refct = 0;
450: PetscDSDestroyStructs_Static(*prob);
451: for (f = 0; f < (*prob)->Nf; ++f) {
452: PetscObjectDereference((*prob)->disc[f]);
453: }
454: PetscFree3((*prob)->disc, (*prob)->implicit, (*prob)->adjacency);
455: PetscFree7((*prob)->obj,(*prob)->f,(*prob)->g,(*prob)->gp,(*prob)->gt,(*prob)->r,(*prob)->ctx);
456: PetscFree((*prob)->update);
457: PetscFree3((*prob)->fBd,(*prob)->gBd,(*prob)->exactSol);
458: if ((*prob)->ops->destroy) {(*(*prob)->ops->destroy)(*prob);}
459: next = (*prob)->boundary;
460: while (next) {
461: DSBoundary b = next;
463: next = b->next;
464: PetscFree(b->comps);
465: PetscFree(b->ids);
466: PetscFree(b->name);
467: PetscFree(b->labelname);
468: PetscFree(b);
469: }
470: PetscFree((*prob)->constants);
471: PetscHeaderDestroy(prob);
472: return(0);
473: }
475: /*@
476: PetscDSCreate - Creates an empty PetscDS object. The type can then be set with PetscDSSetType().
478: Collective on MPI_Comm
480: Input Parameter:
481: . comm - The communicator for the PetscDS object
483: Output Parameter:
484: . prob - The PetscDS object
486: Level: beginner
488: .seealso: PetscDSSetType(), PETSCDSBASIC
489: @*/
490: PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *prob)
491: {
492: PetscDS p;
497: *prob = NULL;
498: PetscDSInitializePackage();
500: PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView);
502: p->Nf = 0;
503: p->setup = PETSC_FALSE;
504: p->numConstants = 0;
505: p->constants = NULL;
507: *prob = p;
508: return(0);
509: }
511: /*@
512: PetscDSGetNumFields - Returns the number of fields in the DS
514: Not collective
516: Input Parameter:
517: . prob - The PetscDS object
519: Output Parameter:
520: . Nf - The number of fields
522: Level: beginner
524: .seealso: PetscDSGetSpatialDimension(), PetscDSCreate()
525: @*/
526: PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf)
527: {
531: *Nf = prob->Nf;
532: return(0);
533: }
535: /*@
536: PetscDSGetSpatialDimension - Returns the spatial dimension of the DS
538: Not collective
540: Input Parameter:
541: . prob - The PetscDS object
543: Output Parameter:
544: . dim - The spatial dimension
546: Level: beginner
548: .seealso: PetscDSGetNumFields(), PetscDSCreate()
549: @*/
550: PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim)
551: {
557: *dim = 0;
558: if (prob->Nf) {
559: PetscObject obj;
560: PetscClassId id;
562: PetscDSGetDiscretization(prob, 0, &obj);
563: PetscObjectGetClassId(obj, &id);
564: if (id == PETSCFE_CLASSID) {PetscFEGetSpatialDimension((PetscFE) obj, dim);}
565: else if (id == PETSCFV_CLASSID) {PetscFVGetSpatialDimension((PetscFV) obj, dim);}
566: else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
567: }
568: return(0);
569: }
571: /*@
572: PetscDSGetTotalDimension - Returns the total size of the approximation space for this system
574: Not collective
576: Input Parameter:
577: . prob - The PetscDS object
579: Output Parameter:
580: . dim - The total problem dimension
582: Level: beginner
584: .seealso: PetscDSGetNumFields(), PetscDSCreate()
585: @*/
586: PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim)
587: {
592: PetscDSSetUp(prob);
594: *dim = prob->totDim;
595: return(0);
596: }
598: /*@
599: PetscDSGetTotalComponents - Returns the total number of components in this system
601: Not collective
603: Input Parameter:
604: . prob - The PetscDS object
606: Output Parameter:
607: . dim - The total number of components
609: Level: beginner
611: .seealso: PetscDSGetNumFields(), PetscDSCreate()
612: @*/
613: PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc)
614: {
619: PetscDSSetUp(prob);
621: *Nc = prob->totComp;
622: return(0);
623: }
625: /*@
626: PetscDSGetDiscretization - Returns the discretization object for the given field
628: Not collective
630: Input Parameters:
631: + prob - The PetscDS object
632: - f - The field number
634: Output Parameter:
635: . disc - The discretization object
637: Level: beginner
639: .seealso: PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
640: @*/
641: PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc)
642: {
646: 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);
647: *disc = prob->disc[f];
648: return(0);
649: }
651: /*@
652: PetscDSSetDiscretization - Sets the discretization object for the given field
654: Not collective
656: Input Parameters:
657: + prob - The PetscDS object
658: . f - The field number
659: - disc - The discretization object
661: Level: beginner
663: .seealso: PetscDSGetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
664: @*/
665: PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc)
666: {
672: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
673: PetscDSEnlarge_Static(prob, f+1);
674: if (prob->disc[f]) {PetscObjectDereference(prob->disc[f]);}
675: prob->disc[f] = disc;
676: PetscObjectReference(disc);
677: {
678: PetscClassId id;
680: PetscObjectGetClassId(disc, &id);
681: if (id == PETSCFV_CLASSID) {
682: prob->implicit[f] = PETSC_FALSE;
683: prob->adjacency[f*2+0] = PETSC_TRUE;
684: prob->adjacency[f*2+1] = PETSC_FALSE;
685: }
686: }
687: return(0);
688: }
690: /*@
691: PetscDSAddDiscretization - Adds a discretization object
693: Not collective
695: Input Parameters:
696: + prob - The PetscDS object
697: - disc - The boundary discretization object
699: Level: beginner
701: .seealso: PetscDSGetDiscretization(), PetscDSSetDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
702: @*/
703: PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc)
704: {
708: PetscDSSetDiscretization(prob, prob->Nf, disc);
709: return(0);
710: }
712: /*@
713: PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for IMEX
715: Not collective
717: Input Parameters:
718: + prob - The PetscDS object
719: - f - The field number
721: Output Parameter:
722: . implicit - The flag indicating what kind of solve to use for this field
724: Level: developer
726: .seealso: PetscDSSetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
727: @*/
728: PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit)
729: {
733: 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);
734: *implicit = prob->implicit[f];
735: return(0);
736: }
738: /*@
739: PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for IMEX
741: Not collective
743: Input Parameters:
744: + prob - The PetscDS object
745: . f - The field number
746: - implicit - The flag indicating what kind of solve to use for this field
748: Level: developer
750: .seealso: PetscDSGetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
751: @*/
752: PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit)
753: {
756: 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);
757: prob->implicit[f] = implicit;
758: return(0);
759: }
761: /*@
762: PetscDSGetAdjacency - Returns the flags for determining variable influence
764: Not collective
766: Input Parameters:
767: + prob - The PetscDS object
768: - f - The field number
770: Output Parameter:
771: + useCone - Flag for variable influence starting with the cone operation
772: - useClosure - Flag for variable influence using transitive closure
774: Note: See the discussion in DMPlexGetAdjacencyUseCone() and DMPlexGetAdjacencyUseClosure()
776: Level: developer
778: .seealso: PetscDSSetAdjacency(), DMPlexGetAdjacencyUseCone(), DMPlexGetAdjacencyUseClosure(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
779: @*/
780: PetscErrorCode PetscDSGetAdjacency(PetscDS prob, PetscInt f, PetscBool *useCone, PetscBool *useClosure)
781: {
786: 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);
787: *useCone = prob->adjacency[f*2+0];
788: *useClosure = prob->adjacency[f*2+1];
789: return(0);
790: }
792: /*@
793: PetscDSSetAdjacency - Set the flags for determining variable influence
795: Not collective
797: Input Parameters:
798: + prob - The PetscDS object
799: . f - The field number
800: . useCone - Flag for variable influence starting with the cone operation
801: - useClosure - Flag for variable influence using transitive closure
803: Note: See the discussion in DMPlexGetAdjacencyUseCone() and DMPlexGetAdjacencyUseClosure()
805: Level: developer
807: .seealso: PetscDSGetAdjacency(), DMPlexGetAdjacencyUseCone(), DMPlexGetAdjacencyUseClosure(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
808: @*/
809: PetscErrorCode PetscDSSetAdjacency(PetscDS prob, PetscInt f, PetscBool useCone, PetscBool useClosure)
810: {
813: 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);
814: prob->adjacency[f*2+0] = useCone;
815: prob->adjacency[f*2+1] = useClosure;
816: return(0);
817: }
819: PetscErrorCode PetscDSGetObjective(PetscDS prob, PetscInt f,
820: void (**obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
821: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
822: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
823: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
824: {
828: 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);
829: *obj = prob->obj[f];
830: return(0);
831: }
833: PetscErrorCode PetscDSSetObjective(PetscDS prob, PetscInt f,
834: void (*obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
835: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
836: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
837: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
838: {
844: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
845: PetscDSEnlarge_Static(prob, f+1);
846: prob->obj[f] = obj;
847: return(0);
848: }
850: /*@C
851: PetscDSGetResidual - Get the pointwise residual function for a given test field
853: Not collective
855: Input Parameters:
856: + prob - The PetscDS
857: - f - The test field number
859: Output Parameters:
860: + f0 - integrand for the test function term
861: - f1 - integrand for the test function gradient term
863: Note: We are using a first order FEM model for the weak form:
865: \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)
867: The calling sequence for the callbacks f0 and f1 is given by:
869: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
870: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
871: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
872: $ PetscReal t, const PetscReal x[], PetscScalar f0[])
874: + dim - the spatial dimension
875: . Nf - the number of fields
876: . uOff - the offset into u[] and u_t[] for each field
877: . uOff_x - the offset into u_x[] for each field
878: . u - each field evaluated at the current point
879: . u_t - the time derivative of each field evaluated at the current point
880: . u_x - the gradient of each field evaluated at the current point
881: . aOff - the offset into a[] and a_t[] for each auxiliary field
882: . aOff_x - the offset into a_x[] for each auxiliary field
883: . a - each auxiliary field evaluated at the current point
884: . a_t - the time derivative of each auxiliary field evaluated at the current point
885: . a_x - the gradient of auxiliary each field evaluated at the current point
886: . t - current time
887: . x - coordinates of the current point
888: . numConstants - number of constant parameters
889: . constants - constant parameters
890: - f0 - output values at the current point
892: Level: intermediate
894: .seealso: PetscDSSetResidual()
895: @*/
896: PetscErrorCode PetscDSGetResidual(PetscDS prob, PetscInt f,
897: void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
898: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
899: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
900: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
901: void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
902: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
903: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
904: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
905: {
908: 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);
911: return(0);
912: }
914: /*@C
915: PetscDSSetResidual - Set the pointwise residual function for a given test field
917: Not collective
919: Input Parameters:
920: + prob - The PetscDS
921: . f - The test field number
922: . f0 - integrand for the test function term
923: - f1 - integrand for the test function gradient term
925: Note: We are using a first order FEM model for the weak form:
927: \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)
929: The calling sequence for the callbacks f0 and f1 is given by:
931: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
932: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
933: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
934: $ PetscReal t, const PetscReal x[], PetscScalar f0[])
936: + dim - the spatial dimension
937: . Nf - the number of fields
938: . uOff - the offset into u[] and u_t[] for each field
939: . uOff_x - the offset into u_x[] for each field
940: . u - each field evaluated at the current point
941: . u_t - the time derivative of each field evaluated at the current point
942: . u_x - the gradient of each field evaluated at the current point
943: . aOff - the offset into a[] and a_t[] for each auxiliary field
944: . aOff_x - the offset into a_x[] for each auxiliary field
945: . a - each auxiliary field evaluated at the current point
946: . a_t - the time derivative of each auxiliary field evaluated at the current point
947: . a_x - the gradient of auxiliary each field evaluated at the current point
948: . t - current time
949: . x - coordinates of the current point
950: . numConstants - number of constant parameters
951: . constants - constant parameters
952: - f0 - output values at the current point
954: Level: intermediate
956: .seealso: PetscDSGetResidual()
957: @*/
958: PetscErrorCode PetscDSSetResidual(PetscDS prob, PetscInt f,
959: void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
960: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
961: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
962: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
963: void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
964: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
965: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
966: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
967: {
974: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
975: PetscDSEnlarge_Static(prob, f+1);
976: prob->f[f*2+0] = f0;
977: prob->f[f*2+1] = f1;
978: return(0);
979: }
981: /*@C
982: PetscDSHasJacobian - Signals that Jacobian functions have been set
984: Not collective
986: Input Parameter:
987: . prob - The PetscDS
989: Output Parameter:
990: . hasJac - flag that pointwise function for the Jacobian has been set
992: Level: intermediate
994: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
995: @*/
996: PetscErrorCode PetscDSHasJacobian(PetscDS prob, PetscBool *hasJac)
997: {
998: PetscInt f, g, h;
1002: *hasJac = PETSC_FALSE;
1003: for (f = 0; f < prob->Nf; ++f) {
1004: for (g = 0; g < prob->Nf; ++g) {
1005: for (h = 0; h < 4; ++h) {
1006: if (prob->g[(f*prob->Nf + g)*4+h]) *hasJac = PETSC_TRUE;
1007: }
1008: }
1009: }
1010: return(0);
1011: }
1013: /*@C
1014: PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field
1016: Not collective
1018: Input Parameters:
1019: + prob - The PetscDS
1020: . f - The test field number
1021: - g - The field number
1023: Output Parameters:
1024: + g0 - integrand for the test and basis function term
1025: . g1 - integrand for the test function and basis function gradient term
1026: . g2 - integrand for the test function gradient and basis function term
1027: - g3 - integrand for the test function gradient and basis function gradient term
1029: Note: We are using a first order FEM model for the weak form:
1031: \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
1033: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1035: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1036: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1037: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1038: $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1040: + dim - the spatial dimension
1041: . Nf - the number of fields
1042: . uOff - the offset into u[] and u_t[] for each field
1043: . uOff_x - the offset into u_x[] for each field
1044: . u - each field evaluated at the current point
1045: . u_t - the time derivative of each field evaluated at the current point
1046: . u_x - the gradient of each field evaluated at the current point
1047: . aOff - the offset into a[] and a_t[] for each auxiliary field
1048: . aOff_x - the offset into a_x[] for each auxiliary field
1049: . a - each auxiliary field evaluated at the current point
1050: . a_t - the time derivative of each auxiliary field evaluated at the current point
1051: . a_x - the gradient of auxiliary each field evaluated at the current point
1052: . t - current time
1053: . u_tShift - the multiplier a for dF/dU_t
1054: . x - coordinates of the current point
1055: . numConstants - number of constant parameters
1056: . constants - constant parameters
1057: - g0 - output values at the current point
1059: Level: intermediate
1061: .seealso: PetscDSSetJacobian()
1062: @*/
1063: PetscErrorCode PetscDSGetJacobian(PetscDS prob, PetscInt f, PetscInt g,
1064: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1065: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1066: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1067: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1068: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1069: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1070: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1071: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1072: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1073: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1074: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1075: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1076: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1077: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1078: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1079: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1080: {
1083: 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);
1084: 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);
1089: return(0);
1090: }
1092: /*@C
1093: PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields
1095: Not collective
1097: Input Parameters:
1098: + prob - The PetscDS
1099: . f - The test field number
1100: . g - The field number
1101: . g0 - integrand for the test and basis function term
1102: . g1 - integrand for the test function and basis function gradient term
1103: . g2 - integrand for the test function gradient and basis function term
1104: - g3 - integrand for the test function gradient and basis function gradient term
1106: Note: We are using a first order FEM model for the weak form:
1108: \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
1110: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1112: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1113: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1114: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1115: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1117: + dim - the spatial dimension
1118: . Nf - the number of fields
1119: . uOff - the offset into u[] and u_t[] for each field
1120: . uOff_x - the offset into u_x[] for each field
1121: . u - each field evaluated at the current point
1122: . u_t - the time derivative of each field evaluated at the current point
1123: . u_x - the gradient of each field evaluated at the current point
1124: . aOff - the offset into a[] and a_t[] for each auxiliary field
1125: . aOff_x - the offset into a_x[] for each auxiliary field
1126: . a - each auxiliary field evaluated at the current point
1127: . a_t - the time derivative of each auxiliary field evaluated at the current point
1128: . a_x - the gradient of auxiliary each field evaluated at the current point
1129: . t - current time
1130: . u_tShift - the multiplier a for dF/dU_t
1131: . x - coordinates of the current point
1132: . numConstants - number of constant parameters
1133: . constants - constant parameters
1134: - g0 - output values at the current point
1136: Level: intermediate
1138: .seealso: PetscDSGetJacobian()
1139: @*/
1140: PetscErrorCode PetscDSSetJacobian(PetscDS prob, PetscInt f, PetscInt g,
1141: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1142: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1143: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1144: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1145: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1146: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1147: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1148: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1149: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1150: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1151: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1152: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1153: void (*g3)(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 g3[]))
1157: {
1166: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1167: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1168: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1169: prob->g[(f*prob->Nf + g)*4+0] = g0;
1170: prob->g[(f*prob->Nf + g)*4+1] = g1;
1171: prob->g[(f*prob->Nf + g)*4+2] = g2;
1172: prob->g[(f*prob->Nf + g)*4+3] = g3;
1173: return(0);
1174: }
1176: /*@C
1177: PetscDSHasJacobianPreconditioner - Signals that a Jacobian preconditioner matrix has been set
1179: Not collective
1181: Input Parameter:
1182: . prob - The PetscDS
1184: Output Parameter:
1185: . hasJacPre - flag that pointwise function for Jacobian preconditioner matrix has been set
1187: Level: intermediate
1189: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1190: @*/
1191: PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS prob, PetscBool *hasJacPre)
1192: {
1193: PetscInt f, g, h;
1197: *hasJacPre = PETSC_FALSE;
1198: for (f = 0; f < prob->Nf; ++f) {
1199: for (g = 0; g < prob->Nf; ++g) {
1200: for (h = 0; h < 4; ++h) {
1201: if (prob->gp[(f*prob->Nf + g)*4+h]) *hasJacPre = PETSC_TRUE;
1202: }
1203: }
1204: }
1205: return(0);
1206: }
1208: /*@C
1209: 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.
1211: Not collective
1213: Input Parameters:
1214: + prob - The PetscDS
1215: . f - The test field number
1216: - g - The field number
1218: Output Parameters:
1219: + g0 - integrand for the test and basis function term
1220: . g1 - integrand for the test function and basis function gradient term
1221: . g2 - integrand for the test function gradient and basis function term
1222: - g3 - integrand for the test function gradient and basis function gradient term
1224: Note: We are using a first order FEM model for the weak form:
1226: \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
1228: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1230: $ 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, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1235: + dim - the spatial dimension
1236: . Nf - the number of fields
1237: . uOff - the offset into u[] and u_t[] for each field
1238: . uOff_x - the offset into u_x[] for each field
1239: . u - each field evaluated at the current point
1240: . u_t - the time derivative of each field evaluated at the current point
1241: . u_x - the gradient of each field evaluated at the current point
1242: . aOff - the offset into a[] and a_t[] for each auxiliary field
1243: . aOff_x - the offset into a_x[] for each auxiliary field
1244: . a - each auxiliary field evaluated at the current point
1245: . a_t - the time derivative of each auxiliary field evaluated at the current point
1246: . a_x - the gradient of auxiliary each field evaluated at the current point
1247: . t - current time
1248: . u_tShift - the multiplier a for dF/dU_t
1249: . x - coordinates of the current point
1250: . numConstants - number of constant parameters
1251: . constants - constant parameters
1252: - g0 - output values at the current point
1254: Level: intermediate
1256: .seealso: PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1257: @*/
1258: PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g,
1259: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1260: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1261: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1262: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1263: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1264: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1265: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1266: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1267: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1268: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1269: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1270: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1271: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1272: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1273: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1274: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1275: {
1278: 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);
1279: 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);
1284: return(0);
1285: }
1287: /*@C
1288: 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.
1290: Not collective
1292: Input Parameters:
1293: + prob - The PetscDS
1294: . f - The test field number
1295: . g - The field number
1296: . g0 - integrand for the test and basis function term
1297: . g1 - integrand for the test function and basis function gradient term
1298: . g2 - integrand for the test function gradient and basis function term
1299: - g3 - integrand for the test function gradient and basis function gradient term
1301: Note: We are using a first order FEM model for the weak form:
1303: \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
1305: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1307: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1308: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1309: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1310: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1312: + dim - the spatial dimension
1313: . Nf - the number of fields
1314: . uOff - the offset into u[] and u_t[] for each field
1315: . uOff_x - the offset into u_x[] for each field
1316: . u - each field evaluated at the current point
1317: . u_t - the time derivative of each field evaluated at the current point
1318: . u_x - the gradient of each field evaluated at the current point
1319: . aOff - the offset into a[] and a_t[] for each auxiliary field
1320: . aOff_x - the offset into a_x[] for each auxiliary field
1321: . a - each auxiliary field evaluated at the current point
1322: . a_t - the time derivative of each auxiliary field evaluated at the current point
1323: . a_x - the gradient of auxiliary each field evaluated at the current point
1324: . t - current time
1325: . u_tShift - the multiplier a for dF/dU_t
1326: . x - coordinates of the current point
1327: . numConstants - number of constant parameters
1328: . constants - constant parameters
1329: - g0 - output values at the current point
1331: Level: intermediate
1333: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobian()
1334: @*/
1335: PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g,
1336: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1337: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1338: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1339: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1340: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1341: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1342: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1343: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1344: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1345: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1346: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1347: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1348: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1349: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1350: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1351: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1352: {
1361: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1362: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1363: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1364: prob->gp[(f*prob->Nf + g)*4+0] = g0;
1365: prob->gp[(f*prob->Nf + g)*4+1] = g1;
1366: prob->gp[(f*prob->Nf + g)*4+2] = g2;
1367: prob->gp[(f*prob->Nf + g)*4+3] = g3;
1368: return(0);
1369: }
1371: /*@C
1372: PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, dF/du_t, has been set
1374: Not collective
1376: Input Parameter:
1377: . prob - The PetscDS
1379: Output Parameter:
1380: . hasDynJac - flag that pointwise function for dynamic Jacobian has been set
1382: Level: intermediate
1384: .seealso: PetscDSGetDynamicJacobian(), PetscDSSetDynamicJacobian(), PetscDSGetJacobian()
1385: @*/
1386: PetscErrorCode PetscDSHasDynamicJacobian(PetscDS prob, PetscBool *hasDynJac)
1387: {
1388: PetscInt f, g, h;
1392: *hasDynJac = PETSC_FALSE;
1393: for (f = 0; f < prob->Nf; ++f) {
1394: for (g = 0; g < prob->Nf; ++g) {
1395: for (h = 0; h < 4; ++h) {
1396: if (prob->gt[(f*prob->Nf + g)*4+h]) *hasDynJac = PETSC_TRUE;
1397: }
1398: }
1399: }
1400: return(0);
1401: }
1403: /*@C
1404: PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, dF/du_t, function for given test and basis field
1406: Not collective
1408: Input Parameters:
1409: + prob - The PetscDS
1410: . f - The test field number
1411: - g - The field number
1413: Output Parameters:
1414: + g0 - integrand for the test and basis function term
1415: . g1 - integrand for the test function and basis function gradient term
1416: . g2 - integrand for the test function gradient and basis function term
1417: - g3 - integrand for the test function gradient and basis function gradient term
1419: Note: We are using a first order FEM model for the weak form:
1421: \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
1423: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1425: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1426: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1427: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1428: $ PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])
1430: + dim - the spatial dimension
1431: . Nf - the number of fields
1432: . uOff - the offset into u[] and u_t[] for each field
1433: . uOff_x - the offset into u_x[] for each field
1434: . u - each field evaluated at the current point
1435: . u_t - the time derivative of each field evaluated at the current point
1436: . u_x - the gradient of each field evaluated at the current point
1437: . aOff - the offset into a[] and a_t[] for each auxiliary field
1438: . aOff_x - the offset into a_x[] for each auxiliary field
1439: . a - each auxiliary field evaluated at the current point
1440: . a_t - the time derivative of each auxiliary field evaluated at the current point
1441: . a_x - the gradient of auxiliary each field evaluated at the current point
1442: . t - current time
1443: . u_tShift - the multiplier a for dF/dU_t
1444: . x - coordinates of the current point
1445: . numConstants - number of constant parameters
1446: . constants - constant parameters
1447: - g0 - output values at the current point
1449: Level: intermediate
1451: .seealso: PetscDSSetJacobian()
1452: @*/
1453: PetscErrorCode PetscDSGetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g,
1454: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1455: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1456: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1457: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1458: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1459: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1460: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1461: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1462: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1463: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1464: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1465: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1466: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1467: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1468: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1469: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1470: {
1473: 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);
1474: 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);
1479: return(0);
1480: }
1482: /*@C
1483: PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, dF/du_t, function for given test and basis fields
1485: Not collective
1487: Input Parameters:
1488: + prob - The PetscDS
1489: . f - The test field number
1490: . g - The field number
1491: . g0 - integrand for the test and basis function term
1492: . g1 - integrand for the test function and basis function gradient term
1493: . g2 - integrand for the test function gradient and basis function term
1494: - g3 - integrand for the test function gradient and basis function gradient term
1496: Note: We are using a first order FEM model for the weak form:
1498: \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
1500: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1502: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1503: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1504: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1505: $ PetscReal t, const PetscReal x[], PetscScalar g0[])
1507: + dim - the spatial dimension
1508: . Nf - the number of fields
1509: . uOff - the offset into u[] and u_t[] for each field
1510: . uOff_x - the offset into u_x[] for each field
1511: . u - each field evaluated at the current point
1512: . u_t - the time derivative of each field evaluated at the current point
1513: . u_x - the gradient of each field evaluated at the current point
1514: . aOff - the offset into a[] and a_t[] for each auxiliary field
1515: . aOff_x - the offset into a_x[] for each auxiliary field
1516: . a - each auxiliary field evaluated at the current point
1517: . a_t - the time derivative of each auxiliary field evaluated at the current point
1518: . a_x - the gradient of auxiliary each field evaluated at the current point
1519: . t - current time
1520: . u_tShift - the multiplier a for dF/dU_t
1521: . x - coordinates of the current point
1522: . numConstants - number of constant parameters
1523: . constants - constant parameters
1524: - g0 - output values at the current point
1526: Level: intermediate
1528: .seealso: PetscDSGetJacobian()
1529: @*/
1530: PetscErrorCode PetscDSSetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g,
1531: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1532: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1533: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1534: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1535: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1536: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1537: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1538: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1539: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1540: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1541: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1542: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1543: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1544: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1545: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1546: PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1547: {
1556: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1557: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1558: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1559: prob->gt[(f*prob->Nf + g)*4+0] = g0;
1560: prob->gt[(f*prob->Nf + g)*4+1] = g1;
1561: prob->gt[(f*prob->Nf + g)*4+2] = g2;
1562: prob->gt[(f*prob->Nf + g)*4+3] = g3;
1563: return(0);
1564: }
1566: /*@C
1567: PetscDSGetRiemannSolver - Returns the Riemann solver for the given field
1569: Not collective
1571: Input Arguments:
1572: + prob - The PetscDS object
1573: - f - The field number
1575: Output Argument:
1576: . r - Riemann solver
1578: Calling sequence for r:
1580: $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)
1582: + dim - The spatial dimension
1583: . Nf - The number of fields
1584: . x - The coordinates at a point on the interface
1585: . n - The normal vector to the interface
1586: . uL - The state vector to the left of the interface
1587: . uR - The state vector to the right of the interface
1588: . flux - output array of flux through the interface
1589: . numConstants - number of constant parameters
1590: . constants - constant parameters
1591: - ctx - optional user context
1593: Level: intermediate
1595: .seealso: PetscDSSetRiemannSolver()
1596: @*/
1597: PetscErrorCode PetscDSGetRiemannSolver(PetscDS prob, PetscInt f,
1598: 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))
1599: {
1602: if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1604: *r = prob->r[f];
1605: return(0);
1606: }
1608: /*@C
1609: PetscDSSetRiemannSolver - Sets the Riemann solver for the given field
1611: Not collective
1613: Input Arguments:
1614: + prob - The PetscDS object
1615: . f - The field number
1616: - r - Riemann solver
1618: Calling sequence for r:
1620: $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)
1622: + dim - The spatial dimension
1623: . Nf - The number of fields
1624: . x - The coordinates at a point on the interface
1625: . n - The normal vector to the interface
1626: . uL - The state vector to the left of the interface
1627: . uR - The state vector to the right of the interface
1628: . flux - output array of flux through the interface
1629: . numConstants - number of constant parameters
1630: . constants - constant parameters
1631: - ctx - optional user context
1633: Level: intermediate
1635: .seealso: PetscDSGetRiemannSolver()
1636: @*/
1637: PetscErrorCode PetscDSSetRiemannSolver(PetscDS prob, PetscInt f,
1638: 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))
1639: {
1645: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1646: PetscDSEnlarge_Static(prob, f+1);
1647: prob->r[f] = r;
1648: return(0);
1649: }
1651: /*@C
1652: PetscDSGetUpdate - Get the pointwise update function for a given field
1654: Not collective
1656: Input Parameters:
1657: + prob - The PetscDS
1658: - f - The field number
1660: Output Parameters:
1661: + update - update function
1663: Note: The calling sequence for the callback update is given by:
1665: $ update(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1666: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1667: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1668: $ PetscReal t, const PetscReal x[], PetscScalar uNew[])
1670: + dim - the spatial dimension
1671: . Nf - the number of fields
1672: . uOff - the offset into u[] and u_t[] for each field
1673: . uOff_x - the offset into u_x[] for each field
1674: . u - each field evaluated at the current point
1675: . u_t - the time derivative of each field evaluated at the current point
1676: . u_x - the gradient of each field evaluated at the current point
1677: . aOff - the offset into a[] and a_t[] for each auxiliary field
1678: . aOff_x - the offset into a_x[] for each auxiliary field
1679: . a - each auxiliary field evaluated at the current point
1680: . a_t - the time derivative of each auxiliary field evaluated at the current point
1681: . a_x - the gradient of auxiliary each field evaluated at the current point
1682: . t - current time
1683: . x - coordinates of the current point
1684: - uNew - new value for field at the current point
1686: Level: intermediate
1688: .seealso: PetscDSSetUpdate(), PetscDSSetResidual()
1689: @*/
1690: PetscErrorCode PetscDSGetUpdate(PetscDS prob, PetscInt f,
1691: void (**update)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1692: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1693: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1694: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
1695: {
1698: 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);
1700: return(0);
1701: }
1703: /*@C
1704: PetscDSSetUpdate - Set the pointwise update function for a given field
1706: Not collective
1708: Input Parameters:
1709: + prob - The PetscDS
1710: . f - The field number
1711: - update - update function
1713: Note: The calling sequence for the callback update is given by:
1715: $ update(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1716: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1717: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1718: $ PetscReal t, const PetscReal x[], PetscScalar uNew[])
1720: + dim - the spatial dimension
1721: . Nf - the number of fields
1722: . uOff - the offset into u[] and u_t[] for each field
1723: . uOff_x - the offset into u_x[] for each field
1724: . u - each field evaluated at the current point
1725: . u_t - the time derivative of each field evaluated at the current point
1726: . u_x - the gradient of each field evaluated at the current point
1727: . aOff - the offset into a[] and a_t[] for each auxiliary field
1728: . aOff_x - the offset into a_x[] for each auxiliary field
1729: . a - each auxiliary field evaluated at the current point
1730: . a_t - the time derivative of each auxiliary field evaluated at the current point
1731: . a_x - the gradient of auxiliary each field evaluated at the current point
1732: . t - current time
1733: . x - coordinates of the current point
1734: - uNew - new field values at the current point
1736: Level: intermediate
1738: .seealso: PetscDSGetResidual()
1739: @*/
1740: PetscErrorCode PetscDSSetUpdate(PetscDS prob, PetscInt f,
1741: void (*update)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1742: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1743: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1744: PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
1745: {
1751: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1752: PetscDSEnlarge_Static(prob, f+1);
1753: prob->update[f] = update;
1754: return(0);
1755: }
1757: PetscErrorCode PetscDSGetContext(PetscDS prob, PetscInt f, void **ctx)
1758: {
1761: 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);
1763: *ctx = prob->ctx[f];
1764: return(0);
1765: }
1767: PetscErrorCode PetscDSSetContext(PetscDS prob, PetscInt f, void *ctx)
1768: {
1773: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1774: PetscDSEnlarge_Static(prob, f+1);
1775: prob->ctx[f] = ctx;
1776: return(0);
1777: }
1779: /*@C
1780: PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field
1782: Not collective
1784: Input Parameters:
1785: + prob - The PetscDS
1786: - f - The test field number
1788: Output Parameters:
1789: + f0 - boundary integrand for the test function term
1790: - f1 - boundary integrand for the test function gradient term
1792: Note: We are using a first order FEM model for the weak form:
1794: \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
1796: The calling sequence for the callbacks f0 and f1 is given by:
1798: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1799: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1800: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1801: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[])
1803: + dim - the spatial dimension
1804: . Nf - the number of fields
1805: . uOff - the offset into u[] and u_t[] for each field
1806: . uOff_x - the offset into u_x[] for each field
1807: . u - each field evaluated at the current point
1808: . u_t - the time derivative of each field evaluated at the current point
1809: . u_x - the gradient of each field evaluated at the current point
1810: . aOff - the offset into a[] and a_t[] for each auxiliary field
1811: . aOff_x - the offset into a_x[] for each auxiliary field
1812: . a - each auxiliary field evaluated at the current point
1813: . a_t - the time derivative of each auxiliary field evaluated at the current point
1814: . a_x - the gradient of auxiliary each field evaluated at the current point
1815: . t - current time
1816: . x - coordinates of the current point
1817: . n - unit normal at the current point
1818: . numConstants - number of constant parameters
1819: . constants - constant parameters
1820: - f0 - output values at the current point
1822: Level: intermediate
1824: .seealso: PetscDSSetBdResidual()
1825: @*/
1826: PetscErrorCode PetscDSGetBdResidual(PetscDS prob, PetscInt f,
1827: void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1828: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1829: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1830: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1831: void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1832: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1833: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1834: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1835: {
1838: 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);
1841: return(0);
1842: }
1844: /*@C
1845: PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field
1847: Not collective
1849: Input Parameters:
1850: + prob - The PetscDS
1851: . f - The test field number
1852: . f0 - boundary integrand for the test function term
1853: - f1 - boundary integrand for the test function gradient term
1855: Note: We are using a first order FEM model for the weak form:
1857: \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
1859: The calling sequence for the callbacks f0 and f1 is given by:
1861: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1862: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1863: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1864: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[])
1866: + dim - the spatial dimension
1867: . Nf - the number of fields
1868: . uOff - the offset into u[] and u_t[] for each field
1869: . uOff_x - the offset into u_x[] for each field
1870: . u - each field evaluated at the current point
1871: . u_t - the time derivative of each field evaluated at the current point
1872: . u_x - the gradient of each field evaluated at the current point
1873: . aOff - the offset into a[] and a_t[] for each auxiliary field
1874: . aOff_x - the offset into a_x[] for each auxiliary field
1875: . a - each auxiliary field evaluated at the current point
1876: . a_t - the time derivative of each auxiliary field evaluated at the current point
1877: . a_x - the gradient of auxiliary each field evaluated at the current point
1878: . t - current time
1879: . x - coordinates of the current point
1880: . n - unit normal at the current point
1881: . numConstants - number of constant parameters
1882: . constants - constant parameters
1883: - f0 - output values at the current point
1885: Level: intermediate
1887: .seealso: PetscDSGetBdResidual()
1888: @*/
1889: PetscErrorCode PetscDSSetBdResidual(PetscDS prob, PetscInt f,
1890: void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1891: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1892: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1893: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1894: void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1895: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1896: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1897: PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1898: {
1903: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1904: PetscDSEnlarge_Static(prob, f+1);
1907: return(0);
1908: }
1910: /*@C
1911: PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field
1913: Not collective
1915: Input Parameters:
1916: + prob - The PetscDS
1917: . f - The test field number
1918: - g - The field number
1920: Output Parameters:
1921: + g0 - integrand for the test and basis function term
1922: . g1 - integrand for the test function and basis function gradient term
1923: . g2 - integrand for the test function gradient and basis function term
1924: - g3 - integrand for the test function gradient and basis function gradient term
1926: Note: We are using a first order FEM model for the weak form:
1928: \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
1930: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
1932: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1933: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1934: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1935: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[])
1937: + dim - the spatial dimension
1938: . Nf - the number of fields
1939: . uOff - the offset into u[] and u_t[] for each field
1940: . uOff_x - the offset into u_x[] for each field
1941: . u - each field evaluated at the current point
1942: . u_t - the time derivative of each field evaluated at the current point
1943: . u_x - the gradient of each field evaluated at the current point
1944: . aOff - the offset into a[] and a_t[] for each auxiliary field
1945: . aOff_x - the offset into a_x[] for each auxiliary field
1946: . a - each auxiliary field evaluated at the current point
1947: . a_t - the time derivative of each auxiliary field evaluated at the current point
1948: . a_x - the gradient of auxiliary each field evaluated at the current point
1949: . t - current time
1950: . u_tShift - the multiplier a for dF/dU_t
1951: . x - coordinates of the current point
1952: . n - normal at the current point
1953: . numConstants - number of constant parameters
1954: . constants - constant parameters
1955: - g0 - output values at the current point
1957: Level: intermediate
1959: .seealso: PetscDSSetBdJacobian()
1960: @*/
1961: PetscErrorCode PetscDSGetBdJacobian(PetscDS prob, PetscInt f, PetscInt g,
1962: void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1963: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1964: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1965: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1966: void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1967: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1968: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1969: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1970: void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1971: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1972: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1973: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1974: void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1975: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1976: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1977: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1978: {
1981: 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);
1982: 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);
1987: return(0);
1988: }
1990: /*@C
1991: PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field
1993: Not collective
1995: Input Parameters:
1996: + prob - The PetscDS
1997: . f - The test field number
1998: . g - The field number
1999: . g0 - integrand for the test and basis function term
2000: . g1 - integrand for the test function and basis function gradient term
2001: . g2 - integrand for the test function gradient and basis function term
2002: - g3 - integrand for the test function gradient and basis function gradient term
2004: Note: We are using a first order FEM model for the weak form:
2006: \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
2008: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:
2010: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2011: $ const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2012: $ const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2013: $ PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[])
2015: + dim - the spatial dimension
2016: . Nf - the number of fields
2017: . uOff - the offset into u[] and u_t[] for each field
2018: . uOff_x - the offset into u_x[] for each field
2019: . u - each field evaluated at the current point
2020: . u_t - the time derivative of each field evaluated at the current point
2021: . u_x - the gradient of each field evaluated at the current point
2022: . aOff - the offset into a[] and a_t[] for each auxiliary field
2023: . aOff_x - the offset into a_x[] for each auxiliary field
2024: . a - each auxiliary field evaluated at the current point
2025: . a_t - the time derivative of each auxiliary field evaluated at the current point
2026: . a_x - the gradient of auxiliary each field evaluated at the current point
2027: . t - current time
2028: . u_tShift - the multiplier a for dF/dU_t
2029: . x - coordinates of the current point
2030: . n - normal at the current point
2031: . numConstants - number of constant parameters
2032: . constants - constant parameters
2033: - g0 - output values at the current point
2035: Level: intermediate
2037: .seealso: PetscDSGetBdJacobian()
2038: @*/
2039: PetscErrorCode PetscDSSetBdJacobian(PetscDS prob, PetscInt f, PetscInt g,
2040: void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2041: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2042: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2043: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
2044: void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2045: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2046: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2047: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
2048: void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2049: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2050: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2051: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
2052: void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2053: const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2054: const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2055: PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2056: {
2065: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2066: if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
2067: PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
2068: prob->gBd[(f*prob->Nf + g)*4+0] = g0;
2069: prob->gBd[(f*prob->Nf + g)*4+1] = g1;
2070: prob->gBd[(f*prob->Nf + g)*4+2] = g2;
2071: prob->gBd[(f*prob->Nf + g)*4+3] = g3;
2072: return(0);
2073: }
2075: /*@C
2076: PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field
2078: Not collective
2080: Input Parameters:
2081: + prob - The PetscDS
2082: - f - The test field number
2084: Output Parameter:
2085: . exactSol - exact solution for the test field
2087: Note: The calling sequence for the solution functions is given by:
2089: $ sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2091: + dim - the spatial dimension
2092: . t - current time
2093: . x - coordinates of the current point
2094: . Nc - the number of field components
2095: . u - the solution field evaluated at the current point
2096: - ctx - a user context
2098: Level: intermediate
2100: .seealso: PetscDSSetExactSolution()
2101: @*/
2102: PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx))
2103: {
2106: 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);
2108: return(0);
2109: }
2111: /*@C
2112: PetscDSSetExactSolution - Get the pointwise exact solution function for a given test field
2114: Not collective
2116: Input Parameters:
2117: + prob - The PetscDS
2118: . f - The test field number
2119: - sol - solution function for the test fields
2121: Note: The calling sequence for solution functions is given by:
2123: $ sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)
2125: + dim - the spatial dimension
2126: . t - current time
2127: . x - coordinates of the current point
2128: . Nc - the number of field components
2129: . u - the solution field evaluated at the current point
2130: - ctx - a user context
2132: Level: intermediate
2134: .seealso: PetscDSGetExactSolution()
2135: @*/
2136: PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx))
2137: {
2142: if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2143: PetscDSEnlarge_Static(prob, f+1);
2145: return(0);
2146: }
2148: /*@C
2149: PetscDSGetConstants - Returns the array of constants passed to point functions
2151: Not collective
2153: Input Parameter:
2154: . prob - The PetscDS object
2156: Output Parameters:
2157: + numConstants - The number of constants
2158: - constants - The array of constants, NULL if there are none
2160: Level: intermediate
2162: .seealso: PetscDSSetConstants(), PetscDSCreate()
2163: @*/
2164: PetscErrorCode PetscDSGetConstants(PetscDS prob, PetscInt *numConstants, const PetscScalar *constants[])
2165: {
2170: return(0);
2171: }
2173: /*@C
2174: PetscDSSetConstants - Set the array of constants passed to point functions
2176: Not collective
2178: Input Parameters:
2179: + prob - The PetscDS object
2180: . numConstants - The number of constants
2181: - constants - The array of constants, NULL if there are none
2183: Level: intermediate
2185: .seealso: PetscDSGetConstants(), PetscDSCreate()
2186: @*/
2187: PetscErrorCode PetscDSSetConstants(PetscDS prob, PetscInt numConstants, PetscScalar constants[])
2188: {
2193: if (numConstants != prob->numConstants) {
2194: PetscFree(prob->constants);
2195: prob->constants = NULL;
2196: }
2197: prob->numConstants = numConstants;
2198: if (prob->numConstants) {
2200: PetscMalloc1(prob->numConstants, &prob->constants);
2201: PetscMemcpy(prob->constants, constants, prob->numConstants * sizeof(PetscScalar));
2202: }
2203: return(0);
2204: }
2206: /*@
2207: PetscDSGetFieldIndex - Returns the index of the given field
2209: Not collective
2211: Input Parameters:
2212: + prob - The PetscDS object
2213: - disc - The discretization object
2215: Output Parameter:
2216: . f - The field number
2218: Level: beginner
2220: .seealso: PetscGetDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
2221: @*/
2222: PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f)
2223: {
2224: PetscInt g;
2229: *f = -1;
2230: for (g = 0; g < prob->Nf; ++g) {if (disc == prob->disc[g]) break;}
2231: if (g == prob->Nf) SETERRQ(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS.");
2232: *f = g;
2233: return(0);
2234: }
2236: /*@
2237: PetscDSGetFieldSize - Returns the size of the given field in the full space basis
2239: Not collective
2241: Input Parameters:
2242: + prob - The PetscDS object
2243: - f - The field number
2245: Output Parameter:
2246: . size - The size
2248: Level: beginner
2250: .seealso: PetscDSGetFieldOffset(), PetscDSGetNumFields(), PetscDSCreate()
2251: @*/
2252: PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size)
2253: {
2259: 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);
2260: PetscDSSetUp(prob);
2261: *size = prob->Nb[f];
2262: return(0);
2263: }
2265: /*@
2266: PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis
2268: Not collective
2270: Input Parameters:
2271: + prob - The PetscDS object
2272: - f - The field number
2274: Output Parameter:
2275: . off - The offset
2277: Level: beginner
2279: .seealso: PetscDSGetFieldSize(), PetscDSGetNumFields(), PetscDSCreate()
2280: @*/
2281: PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off)
2282: {
2283: PetscInt size, g;
2289: 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);
2290: *off = 0;
2291: for (g = 0; g < f; ++g) {
2292: PetscDSGetFieldSize(prob, g, &size);
2293: *off += size;
2294: }
2295: return(0);
2296: }
2298: /*@
2299: PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point
2301: Not collective
2303: Input Parameter:
2304: . prob - The PetscDS object
2306: Output Parameter:
2307: . dimensions - The number of dimensions
2309: Level: beginner
2311: .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate()
2312: @*/
2313: PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[])
2314: {
2319: PetscDSSetUp(prob);
2321: *dimensions = prob->Nb;
2322: return(0);
2323: }
2325: /*@
2326: PetscDSGetComponents - Returns the number of components for each field on an evaluation point
2328: Not collective
2330: Input Parameter:
2331: . prob - The PetscDS object
2333: Output Parameter:
2334: . components - The number of components
2336: Level: beginner
2338: .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate()
2339: @*/
2340: PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[])
2341: {
2346: PetscDSSetUp(prob);
2348: *components = prob->Nc;
2349: return(0);
2350: }
2352: /*@
2353: PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point
2355: Not collective
2357: Input Parameters:
2358: + prob - The PetscDS object
2359: - f - The field number
2361: Output Parameter:
2362: . off - The offset
2364: Level: beginner
2366: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2367: @*/
2368: PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off)
2369: {
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: *off = prob->off[f];
2375: return(0);
2376: }
2378: /*@
2379: PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point
2381: Not collective
2383: Input Parameter:
2384: . prob - The PetscDS object
2386: Output Parameter:
2387: . offsets - The offsets
2389: Level: beginner
2391: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2392: @*/
2393: PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[])
2394: {
2398: *offsets = prob->off;
2399: return(0);
2400: }
2402: /*@
2403: PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point
2405: Not collective
2407: Input Parameter:
2408: . prob - The PetscDS object
2410: Output Parameter:
2411: . offsets - The offsets
2413: Level: beginner
2415: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2416: @*/
2417: PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[])
2418: {
2422: *offsets = prob->offDer;
2423: return(0);
2424: }
2426: /*@C
2427: PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization
2429: Not collective
2431: Input Parameter:
2432: . prob - The PetscDS object
2434: Output Parameters:
2435: + basis - The basis function tabulation at quadrature points
2436: - basisDer - The basis function derivative tabulation at quadrature points
2438: Level: intermediate
2440: .seealso: PetscDSCreate()
2441: @*/
2442: PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer)
2443: {
2448: PetscDSSetUp(prob);
2451: return(0);
2452: }
2454: /*@C
2455: PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces
2457: Not collective
2459: Input Parameter:
2460: . prob - The PetscDS object
2462: Output Parameters:
2463: + basisFace - The basis function tabulation at quadrature points
2464: - basisDerFace - The basis function derivative tabulation at quadrature points
2466: Level: intermediate
2468: .seealso: PetscDSGetTabulation(), PetscDSCreate()
2469: @*/
2470: PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer)
2471: {
2476: PetscDSSetUp(prob);
2479: return(0);
2480: }
2482: PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x)
2483: {
2488: PetscDSSetUp(prob);
2492: return(0);
2493: }
2495: PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3)
2496: {
2501: PetscDSSetUp(prob);
2508: return(0);
2509: }
2511: PetscErrorCode PetscDSGetRefCoordArrays(PetscDS prob, PetscReal **x, PetscScalar **refSpaceDer)
2512: {
2517: PetscDSSetUp(prob);
2520: return(0);
2521: }
2523: /*@C
2524: PetscDSAddBoundary - Add a boundary condition to the model
2526: Input Parameters:
2527: + ds - The PetscDS object
2528: . type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2529: . name - The BC name
2530: . labelname - The label defining constrained points
2531: . field - The field to constrain
2532: . numcomps - The number of constrained field components
2533: . comps - An array of constrained component numbers
2534: . bcFunc - A pointwise function giving boundary values
2535: . numids - The number of DMLabel ids for constrained points
2536: . ids - An array of ids for constrained points
2537: - ctx - An optional user context for bcFunc
2539: Options Database Keys:
2540: + -bc_<boundary name> <num> - Overrides the boundary ids
2541: - -bc_<boundary name>_comp <num> - Overrides the boundary components
2543: Level: developer
2545: .seealso: PetscDSGetBoundary()
2546: @*/
2547: 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)
2548: {
2549: DSBoundary b;
2554: PetscNew(&b);
2555: PetscStrallocpy(name, (char **) &b->name);
2556: PetscStrallocpy(labelname, (char **) &b->labelname);
2557: PetscMalloc1(numcomps, &b->comps);
2558: if (numcomps) {PetscMemcpy(b->comps, comps, numcomps*sizeof(PetscInt));}
2559: PetscMalloc1(numids, &b->ids);
2560: if (numids) {PetscMemcpy(b->ids, ids, numids*sizeof(PetscInt));}
2561: b->type = type;
2562: b->field = field;
2563: b->numcomps = numcomps;
2564: b->func = bcFunc;
2565: b->numids = numids;
2566: b->ctx = ctx;
2567: b->next = ds->boundary;
2568: ds->boundary = b;
2569: return(0);
2570: }
2572: /*@
2573: PetscDSGetNumBoundary - Get the number of registered BC
2575: Input Parameters:
2576: . ds - The PetscDS object
2578: Output Parameters:
2579: . numBd - The number of BC
2581: Level: intermediate
2583: .seealso: PetscDSAddBoundary(), PetscDSGetBoundary()
2584: @*/
2585: PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd)
2586: {
2587: DSBoundary b = ds->boundary;
2592: *numBd = 0;
2593: while (b) {++(*numBd); b = b->next;}
2594: return(0);
2595: }
2597: /*@C
2598: PetscDSGetBoundary - Add a boundary condition to the model
2600: Input Parameters:
2601: + ds - The PetscDS object
2602: - bd - The BC number
2604: Output Parameters:
2605: + type - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2606: . name - The BC name
2607: . labelname - The label defining constrained points
2608: . field - The field to constrain
2609: . numcomps - The number of constrained field components
2610: . comps - An array of constrained component numbers
2611: . bcFunc - A pointwise function giving boundary values
2612: . numids - The number of DMLabel ids for constrained points
2613: . ids - An array of ids for constrained points
2614: - ctx - An optional user context for bcFunc
2616: Options Database Keys:
2617: + -bc_<boundary name> <num> - Overrides the boundary ids
2618: - -bc_<boundary name>_comp <num> - Overrides the boundary components
2620: Level: developer
2622: .seealso: PetscDSAddBoundary()
2623: @*/
2624: 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)
2625: {
2626: DSBoundary b = ds->boundary;
2627: PetscInt n = 0;
2631: while (b) {
2632: if (n == bd) break;
2633: b = b->next;
2634: ++n;
2635: }
2636: if (!b) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %d is not in [0, %d)", bd, n);
2637: if (type) {
2639: *type = b->type;
2640: }
2641: if (name) {
2643: *name = b->name;
2644: }
2645: if (labelname) {
2647: *labelname = b->labelname;
2648: }
2649: if (field) {
2651: *field = b->field;
2652: }
2653: if (numcomps) {
2655: *numcomps = b->numcomps;
2656: }
2657: if (comps) {
2659: *comps = b->comps;
2660: }
2661: if (func) {
2663: *func = b->func;
2664: }
2665: if (numids) {
2667: *numids = b->numids;
2668: }
2669: if (ids) {
2671: *ids = b->ids;
2672: }
2673: if (ctx) {
2675: *ctx = b->ctx;
2676: }
2677: return(0);
2678: }
2680: PetscErrorCode PetscDSCopyBoundary(PetscDS probA, PetscDS probB)
2681: {
2682: DSBoundary b, next, *lastnext;
2688: if (probA == probB) return(0);
2689: next = probB->boundary;
2690: while (next) {
2691: DSBoundary b = next;
2693: next = b->next;
2694: PetscFree(b->comps);
2695: PetscFree(b->ids);
2696: PetscFree(b->name);
2697: PetscFree(b->labelname);
2698: PetscFree(b);
2699: }
2700: lastnext = &(probB->boundary);
2701: for (b = probA->boundary; b; b = b->next) {
2702: DSBoundary bNew;
2704: PetscNew(&bNew);
2705: bNew->numcomps = b->numcomps;
2706: PetscMalloc1(bNew->numcomps, &bNew->comps);
2707: PetscMemcpy(bNew->comps, b->comps, bNew->numcomps*sizeof(PetscInt));
2708: bNew->numids = b->numids;
2709: PetscMalloc1(bNew->numids, &bNew->ids);
2710: PetscMemcpy(bNew->ids, b->ids, bNew->numids*sizeof(PetscInt));
2711: PetscStrallocpy(b->labelname,(char **) &(bNew->labelname));
2712: PetscStrallocpy(b->name,(char **) &(bNew->name));
2713: bNew->ctx = b->ctx;
2714: bNew->type = b->type;
2715: bNew->field = b->field;
2716: bNew->func = b->func;
2718: *lastnext = bNew;
2719: lastnext = &(bNew->next);
2720: }
2721: return(0);
2722: }
2724: /*@
2725: PetscDSCopyEquations - Copy all pointwise function pointers to the new problem
2727: Not collective
2729: Input Parameter:
2730: . prob - The PetscDS object
2732: Output Parameter:
2733: . newprob - The PetscDS copy
2735: Level: intermediate
2737: .seealso: PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
2738: @*/
2739: PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob)
2740: {
2741: PetscInt Nf, Ng, f, g;
2747: PetscDSGetNumFields(prob, &Nf);
2748: PetscDSGetNumFields(newprob, &Ng);
2749: if (Nf != Ng) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %D != %D", Nf, Ng);
2750: for (f = 0; f < Nf; ++f) {
2751: PetscPointFunc obj;
2752: PetscPointFunc f0, f1;
2753: PetscPointJac g0, g1, g2, g3;
2754: PetscBdPointFunc f0Bd, f1Bd;
2755: PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd;
2756: PetscRiemannFunc r;
2758: PetscDSGetObjective(prob, f, &obj);
2759: PetscDSGetResidual(prob, f, &f0, &f1);
2760: PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd);
2761: PetscDSGetRiemannSolver(prob, f, &r);
2762: PetscDSSetObjective(newprob, f, obj);
2763: PetscDSSetResidual(newprob, f, f0, f1);
2764: PetscDSSetBdResidual(newprob, f, f0Bd, f1Bd);
2765: PetscDSSetRiemannSolver(newprob, f, r);
2766: for (g = 0; g < Nf; ++g) {
2767: PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3);
2768: PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd);
2769: PetscDSSetJacobian(newprob, f, g, g0, g1, g2, g3);
2770: PetscDSSetBdJacobian(newprob, f, g, g0Bd, g1Bd, g2Bd, g3Bd);
2771: }
2772: }
2773: return(0);
2774: }
2776: static PetscErrorCode PetscDSDestroy_Basic(PetscDS prob)
2777: {
2778: PetscErrorCode ierr;
2781: PetscFree(prob->data);
2782: return(0);
2783: }
2785: static PetscErrorCode PetscDSInitialize_Basic(PetscDS prob)
2786: {
2788: prob->ops->setfromoptions = NULL;
2789: prob->ops->setup = NULL;
2790: prob->ops->view = NULL;
2791: prob->ops->destroy = PetscDSDestroy_Basic;
2792: return(0);
2793: }
2795: /*MC
2796: PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions
2798: Level: intermediate
2800: .seealso: PetscDSType, PetscDSCreate(), PetscDSSetType()
2801: M*/
2803: PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS prob)
2804: {
2805: PetscDS_Basic *b;
2806: PetscErrorCode ierr;
2810: PetscNewLog(prob, &b);
2811: prob->data = b;
2813: PetscDSInitialize_Basic(prob);
2814: return(0);
2815: }