Actual source code: dtds.c

petsc-3.12.1 2019-10-22
Report Typos and Errors
  1:  #include <petsc/private/petscdsimpl.h>

  3: PetscClassId PETSCDS_CLASSID = 0;

  5: PetscFunctionList PetscDSList              = NULL;
  6: PetscBool         PetscDSRegisterAllCalled = PETSC_FALSE;

  8: /* A PetscDS (Discrete System) encodes a set of equations posed in a discrete space, which represents a set of
  9:    nonlinear continuum equations. The equations can have multiple fields, each field having a different
 10:    discretization. In addition, different pieces of the domain can have different field combinations and equations.

 12:    The DS provides the user a description of the approximation space on any given cell. It also gives pointwise
 13:    functions representing the equations.

 15:    Each field is associated with a label, marking the cells on which it is supported. Note that a field can be
 16:    supported on the closure of a cell not in the label due to overlap of the boundary of neighboring cells. The DM
 17:    then creates a DS for each set of cells with identical approximation spaces. When assembling, the user asks for
 18:    the space associated with a given cell. DMPlex uses the labels associated with each DS in the default integration loop.
 19: */

 21: /*@C
 22:   PetscDSRegister - Adds a new PetscDS implementation

 24:   Not Collective

 26:   Input Parameters:
 27: + name        - The name of a new user-defined creation routine
 28: - create_func - The creation routine itself

 30:   Notes:
 31:   PetscDSRegister() may be called multiple times to add several user-defined PetscDSs

 33:   Sample usage:
 34: .vb
 35:     PetscDSRegister("my_ds", MyPetscDSCreate);
 36: .ve

 38:   Then, your PetscDS type can be chosen with the procedural interface via
 39: .vb
 40:     PetscDSCreate(MPI_Comm, PetscDS *);
 41:     PetscDSSetType(PetscDS, "my_ds");
 42: .ve
 43:    or at runtime via the option
 44: .vb
 45:     -petscds_type my_ds
 46: .ve

 48:   Level: advanced

 50:    Not available from Fortran

 52: .seealso: PetscDSRegisterAll(), PetscDSRegisterDestroy()

 54: @*/
 55: PetscErrorCode PetscDSRegister(const char sname[], PetscErrorCode (*function)(PetscDS))
 56: {

 60:   PetscFunctionListAdd(&PetscDSList, sname, function);
 61:   return(0);
 62: }

 64: /*@C
 65:   PetscDSSetType - Builds a particular PetscDS

 67:   Collective on prob

 69:   Input Parameters:
 70: + prob - The PetscDS object
 71: - name - The kind of system

 73:   Options Database Key:
 74: . -petscds_type <type> - Sets the PetscDS type; use -help for a list of available types

 76:   Level: intermediate

 78:    Not available from Fortran

 80: .seealso: PetscDSGetType(), PetscDSCreate()
 81: @*/
 82: PetscErrorCode PetscDSSetType(PetscDS prob, PetscDSType name)
 83: {
 84:   PetscErrorCode (*r)(PetscDS);
 85:   PetscBool      match;

 90:   PetscObjectTypeCompare((PetscObject) prob, name, &match);
 91:   if (match) return(0);

 93:   PetscDSRegisterAll();
 94:   PetscFunctionListFind(PetscDSList, name, &r);
 95:   if (!r) SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_UNKNOWN_TYPE, "Unknown PetscDS type: %s", name);

 97:   if (prob->ops->destroy) {
 98:     (*prob->ops->destroy)(prob);
 99:     prob->ops->destroy = NULL;
100:   }
101:   (*r)(prob);
102:   PetscObjectChangeTypeName((PetscObject) prob, name);
103:   return(0);
104: }

106: /*@C
107:   PetscDSGetType - Gets the PetscDS type name (as a string) from the object.

109:   Not Collective

111:   Input Parameter:
112: . prob  - The PetscDS

114:   Output Parameter:
115: . name - The PetscDS type name

117:   Level: intermediate

119:    Not available from Fortran

121: .seealso: PetscDSSetType(), PetscDSCreate()
122: @*/
123: PetscErrorCode PetscDSGetType(PetscDS prob, PetscDSType *name)
124: {

130:   PetscDSRegisterAll();
131:   *name = ((PetscObject) prob)->type_name;
132:   return(0);
133: }

135: static PetscErrorCode PetscDSView_Ascii(PetscDS prob, PetscViewer viewer)
136: {
137:   PetscViewerFormat  format;
138:   const PetscScalar *constants;
139:   PetscInt           numConstants, f;
140:   PetscErrorCode     ierr;

143:   PetscViewerGetFormat(viewer, &format);
144:   PetscViewerASCIIPrintf(viewer, "Discrete System with %d fields\n", prob->Nf);
145:   PetscViewerASCIIPushTab(viewer);
146:   PetscViewerASCIIPrintf(viewer, "  cell total dim %D total comp %D\n", prob->totDim, prob->totComp);
147:   if (prob->isHybrid) {PetscViewerASCIIPrintf(viewer, "  hybrid cell\n");}
148:   for (f = 0; f < prob->Nf; ++f) {
149:     DSBoundary      b;
150:     PetscObject     obj;
151:     PetscClassId    id;
152:     PetscQuadrature q;
153:     const char     *name;
154:     PetscInt        Nc, Nq, Nqc;

156:     PetscDSGetDiscretization(prob, f, &obj);
157:     PetscObjectGetClassId(obj, &id);
158:     PetscObjectGetName(obj, &name);
159:     PetscViewerASCIIPrintf(viewer, "Field %s", name ? name : "<unknown>");
160:     PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
161:     if (id == PETSCFE_CLASSID)      {
162:       PetscFEGetNumComponents((PetscFE) obj, &Nc);
163:       PetscFEGetQuadrature((PetscFE) obj, &q);
164:       PetscViewerASCIIPrintf(viewer, " FEM");
165:     } else if (id == PETSCFV_CLASSID) {
166:       PetscFVGetNumComponents((PetscFV) obj, &Nc);
167:       PetscFVGetQuadrature((PetscFV) obj, &q);
168:       PetscViewerASCIIPrintf(viewer, " FVM");
169:     }
170:     else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %D", f);
171:     if (Nc > 1) {PetscViewerASCIIPrintf(viewer, "%D components", Nc);}
172:     else        {PetscViewerASCIIPrintf(viewer, "%D component ", Nc);}
173:     if (prob->implicit[f]) {PetscViewerASCIIPrintf(viewer, " (implicit)");}
174:     else                   {PetscViewerASCIIPrintf(viewer, " (explicit)");}
175:     if (q) {
176:       PetscQuadratureGetData(q, NULL, &Nqc, &Nq, NULL, NULL);
177:       PetscViewerASCIIPrintf(viewer, " (Nq %D Nqc %D)", Nq, Nqc);
178:     }
179:     PetscViewerASCIIPrintf(viewer, "\n");
180:     PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
181:     PetscViewerASCIIPushTab(viewer);
182:     if (id == PETSCFE_CLASSID)      {PetscFEView((PetscFE) obj, viewer);}
183:     else if (id == PETSCFV_CLASSID) {PetscFVView((PetscFV) obj, viewer);}
184:     PetscViewerASCIIPopTab(viewer);

186:     for (b = prob->boundary; b; b = b->next) {
187:       PetscInt c, i;

189:       if (b->field != f) continue;
190:       PetscViewerASCIIPushTab(viewer);
191:       PetscViewerASCIIPrintf(viewer, "Boundary %s (%s) %s\n", b->name, b->labelname, DMBoundaryConditionTypes[b->type]);
192:       if (!b->numcomps) {
193:         PetscViewerASCIIPrintf(viewer, "  all components\n");
194:       } else {
195:         PetscViewerASCIIPrintf(viewer, "  components: ");
196:         PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
197:         for (c = 0; c < b->numcomps; ++c) {
198:           if (c > 0) {PetscViewerASCIIPrintf(viewer, ", ");}
199:           PetscViewerASCIIPrintf(viewer, "%D", b->comps[c]);
200:         }
201:         PetscViewerASCIIPrintf(viewer, "\n");
202:         PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
203:       }
204:       PetscViewerASCIIPrintf(viewer, "  ids: ");
205:       PetscViewerASCIIUseTabs(viewer, PETSC_FALSE);
206:       for (i = 0; i < b->numids; ++i) {
207:         if (i > 0) {PetscViewerASCIIPrintf(viewer, ", ");}
208:         PetscViewerASCIIPrintf(viewer, "%D", b->ids[i]);
209:       }
210:       PetscViewerASCIIPrintf(viewer, "\n");
211:       PetscViewerASCIIUseTabs(viewer, PETSC_TRUE);
212:       PetscViewerASCIIPopTab(viewer);
213:     }
214:   }
215:   PetscDSGetConstants(prob, &numConstants, &constants);
216:   if (numConstants) {
217:     PetscViewerASCIIPrintf(viewer, "%D constants\n", numConstants);
218:     PetscViewerASCIIPushTab(viewer);
219:     for (f = 0; f < numConstants; ++f) {PetscViewerASCIIPrintf(viewer, "%g\n", (double) PetscRealPart(constants[f]));}
220:     PetscViewerASCIIPopTab(viewer);
221:   }
222:   PetscViewerASCIIPopTab(viewer);
223:   return(0);
224: }

226: /*@C
227:   PetscDSView - Views a PetscDS

229:   Collective on prob

231:   Input Parameter:
232: + prob - the PetscDS object to view
233: - v  - the viewer

235:   Level: developer

237: .seealso PetscDSDestroy()
238: @*/
239: PetscErrorCode PetscDSView(PetscDS prob, PetscViewer v)
240: {
241:   PetscBool      iascii;

246:   if (!v) {PetscViewerASCIIGetStdout(PetscObjectComm((PetscObject) prob), &v);}
248:   PetscObjectTypeCompare((PetscObject) v, PETSCVIEWERASCII, &iascii);
249:   if (iascii) {PetscDSView_Ascii(prob, v);}
250:   if (prob->ops->view) {(*prob->ops->view)(prob, v);}
251:   return(0);
252: }

254: /*@
255:   PetscDSSetFromOptions - sets parameters in a PetscDS from the options database

257:   Collective on prob

259:   Input Parameter:
260: . prob - the PetscDS object to set options for

262:   Options Database:
263: + -petscds_type <type>     : Set the DS type
264: . -petscds_view <view opt> : View the DS
265: . -petscds_jac_pre         : Turn formation of a separate Jacobian preconditioner on and off
266: . -bc_<name> <ids>         : Specify a list of label ids for a boundary condition
267: - -bc_<name>_comp <comps>  : Specify a list of field components to constrain for a boundary condition

269:   Level: developer

271: .seealso PetscDSView()
272: @*/
273: PetscErrorCode PetscDSSetFromOptions(PetscDS prob)
274: {
275:   DSBoundary     b;
276:   const char    *defaultType;
277:   char           name[256];
278:   PetscBool      flg;

283:   if (!((PetscObject) prob)->type_name) {
284:     defaultType = PETSCDSBASIC;
285:   } else {
286:     defaultType = ((PetscObject) prob)->type_name;
287:   }
288:   PetscDSRegisterAll();

290:   PetscObjectOptionsBegin((PetscObject) prob);
291:   for (b = prob->boundary; b; b = b->next) {
292:     char       optname[1024];
293:     PetscInt   ids[1024], len = 1024;
294:     PetscBool  flg;

296:     PetscSNPrintf(optname, sizeof(optname), "-bc_%s", b->name);
297:     PetscMemzero(ids, sizeof(ids));
298:     PetscOptionsIntArray(optname, "List of boundary IDs", "", ids, &len, &flg);
299:     if (flg) {
300:       b->numids = len;
301:       PetscFree(b->ids);
302:       PetscMalloc1(len, &b->ids);
303:       PetscArraycpy(b->ids, ids, len);
304:     }
305:     len = 1024;
306:     PetscSNPrintf(optname, sizeof(optname), "-bc_%s_comp", b->name);
307:     PetscMemzero(ids, sizeof(ids));
308:     PetscOptionsIntArray(optname, "List of boundary field components", "", ids, &len, &flg);
309:     if (flg) {
310:       b->numcomps = len;
311:       PetscFree(b->comps);
312:       PetscMalloc1(len, &b->comps);
313:       PetscArraycpy(b->comps, ids, len);
314:     }
315:   }
316:   PetscOptionsFList("-petscds_type", "Discrete System", "PetscDSSetType", PetscDSList, defaultType, name, 256, &flg);
317:   if (flg) {
318:     PetscDSSetType(prob, name);
319:   } else if (!((PetscObject) prob)->type_name) {
320:     PetscDSSetType(prob, defaultType);
321:   }
322:   PetscOptionsBool("-petscds_jac_pre", "Discrete System", "PetscDSUseJacobianPreconditioner", prob->useJacPre, &prob->useJacPre, &flg);
323:   if (prob->ops->setfromoptions) {(*prob->ops->setfromoptions)(prob);}
324:   /* process any options handlers added with PetscObjectAddOptionsHandler() */
325:   PetscObjectProcessOptionsHandlers(PetscOptionsObject,(PetscObject) prob);
326:   PetscOptionsEnd();
327:   if (prob->Nf) {PetscDSViewFromOptions(prob, NULL, "-petscds_view");}
328:   return(0);
329: }

331: /*@C
332:   PetscDSSetUp - Construct data structures for the PetscDS

334:   Collective on prob

336:   Input Parameter:
337: . prob - the PetscDS object to setup

339:   Level: developer

341: .seealso PetscDSView(), PetscDSDestroy()
342: @*/
343: PetscErrorCode PetscDSSetUp(PetscDS prob)
344: {
345:   const PetscInt Nf = prob->Nf;
346:   PetscInt       dim, dimEmbed, NbMax = 0, NcMax = 0, NqMax = 0, NsMax = 1, f;

351:   if (prob->setup) return(0);
352:   /* Calculate sizes */
353:   PetscDSGetSpatialDimension(prob, &dim);
354:   PetscDSGetCoordinateDimension(prob, &dimEmbed);
355:   prob->totDim = prob->totComp = 0;
356:   PetscMalloc2(Nf,&prob->Nc,Nf,&prob->Nb);
357:   PetscCalloc2(Nf+1,&prob->off,Nf+1,&prob->offDer);
358:   PetscMalloc4(Nf,&prob->basis,Nf,&prob->basisDer,Nf,&prob->basisFace,Nf,&prob->basisDerFace);
359:   for (f = 0; f < Nf; ++f) {
360:     PetscObject     obj;
361:     PetscClassId    id;
362:     PetscQuadrature q;
363:     PetscInt        Nq = 0, Nb, Nc;

365:     PetscDSGetDiscretization(prob, f, &obj);
366:     PetscObjectGetClassId(obj, &id);
367:     if (id == PETSCFE_CLASSID)      {
368:       PetscFE fe = (PetscFE) obj;

370:       PetscFEGetQuadrature(fe, &q);
371:       PetscFEGetDimension(fe, &Nb);
372:       PetscFEGetNumComponents(fe, &Nc);
373:       PetscFEGetDefaultTabulation(fe, &prob->basis[f], &prob->basisDer[f], NULL);
374:       PetscFEGetFaceTabulation(fe, &prob->basisFace[f], &prob->basisDerFace[f], NULL);
375:     } else if (id == PETSCFV_CLASSID) {
376:       PetscFV fv = (PetscFV) obj;

378:       PetscFVGetQuadrature(fv, &q);
379:       PetscFVGetNumComponents(fv, &Nc);
380:       Nb   = Nc;
381:       PetscFVGetDefaultTabulation(fv, &prob->basis[f], &prob->basisDer[f], NULL);
382:       /* TODO: should PetscFV also have face tabulation? Otherwise there will be a null pointer in prob->basisFace */
383:     } else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", f);
384:     prob->Nc[f]       = Nc;
385:     prob->Nb[f]       = Nb;
386:     prob->off[f+1]    = Nc     + prob->off[f];
387:     prob->offDer[f+1] = Nc*dim + prob->offDer[f];
388:     if (q) {PetscQuadratureGetData(q, NULL, NULL, &Nq, NULL, NULL);}
389:     NqMax          = PetscMax(NqMax, Nq);
390:     NbMax          = PetscMax(NbMax, Nb);
391:     NcMax          = PetscMax(NcMax, Nc);
392:     prob->totDim  += Nb;
393:     prob->totComp += Nc;
394:     /* There are two faces for all fields but the cohesive field on a hybrid cell */
395:     if (prob->isHybrid && (f < Nf-1)) prob->totDim += Nb;
396:   }
397:   /* Allocate works space */
398:   if (prob->isHybrid) NsMax = 2;
399:   PetscMalloc3(NsMax*prob->totComp,&prob->u,NsMax*prob->totComp,&prob->u_t,NsMax*prob->totComp*dimEmbed,&prob->u_x);
400:   PetscMalloc5(dimEmbed,&prob->x,NbMax*NcMax,&prob->basisReal,NbMax*NcMax*dimEmbed,&prob->basisDerReal,NbMax*NcMax,&prob->testReal,NbMax*NcMax*dimEmbed,&prob->testDerReal);
401:   PetscMalloc6(NsMax*NqMax*NcMax,&prob->f0,NsMax*NqMax*NcMax*dim,&prob->f1,
402:                       NsMax*NsMax*NqMax*NcMax*NcMax,&prob->g0,NsMax*NsMax*NqMax*NcMax*NcMax*dim,&prob->g1,
403:                       NsMax*NsMax*NqMax*NcMax*NcMax*dim,&prob->g2,NsMax*NsMax*NqMax*NcMax*NcMax*dim*dim,&prob->g3);
404:   if (prob->ops->setup) {(*prob->ops->setup)(prob);}
405:   prob->setup = PETSC_TRUE;
406:   return(0);
407: }

409: static PetscErrorCode PetscDSDestroyStructs_Static(PetscDS prob)
410: {

414:   PetscFree2(prob->Nc,prob->Nb);
415:   PetscFree2(prob->off,prob->offDer);
416:   PetscFree4(prob->basis,prob->basisDer,prob->basisFace,prob->basisDerFace);
417:   PetscFree3(prob->u,prob->u_t,prob->u_x);
418:   PetscFree5(prob->x,prob->basisReal, prob->basisDerReal,prob->testReal,prob->testDerReal);
419:   PetscFree6(prob->f0,prob->f1,prob->g0,prob->g1,prob->g2,prob->g3);
420:   return(0);
421: }

423: static PetscErrorCode PetscDSEnlarge_Static(PetscDS prob, PetscInt NfNew)
424: {
425:   PetscObject      *tmpd;
426:   PetscBool        *tmpi;
427:   PetscPointFunc   *tmpobj, *tmpf, *tmpup;
428:   PetscPointJac    *tmpg, *tmpgp, *tmpgt;
429:   PetscBdPointFunc *tmpfbd;
430:   PetscBdPointJac  *tmpgbd;
431:   PetscRiemannFunc *tmpr;
432:   PetscSimplePointFunc *tmpexactSol;
433:   void                **tmpexactCtx;
434:   void            **tmpctx;
435:   PetscInt          Nf = prob->Nf, f;
436:   PetscErrorCode    ierr;

439:   if (Nf >= NfNew) return(0);
440:   prob->setup = PETSC_FALSE;
441:   PetscDSDestroyStructs_Static(prob);
442:   PetscMalloc2(NfNew, &tmpd, NfNew, &tmpi);
443:   for (f = 0; f < Nf; ++f) {tmpd[f] = prob->disc[f]; tmpi[f] = prob->implicit[f];}
444:   for (f = Nf; f < NfNew; ++f) {tmpd[f] = NULL; tmpi[f] = PETSC_TRUE;}
445:   PetscFree2(prob->disc, prob->implicit);
446:   prob->Nf        = NfNew;
447:   prob->disc      = tmpd;
448:   prob->implicit  = tmpi;
449:   PetscCalloc7(NfNew, &tmpobj, NfNew*2, &tmpf, NfNew*NfNew*4, &tmpg, NfNew*NfNew*4, &tmpgp, NfNew*NfNew*4, &tmpgt, NfNew, &tmpr, NfNew, &tmpctx);
450:   PetscCalloc1(NfNew, &tmpup);
451:   for (f = 0; f < Nf; ++f) tmpobj[f] = prob->obj[f];
452:   for (f = 0; f < Nf*2; ++f) tmpf[f] = prob->f[f];
453:   for (f = 0; f < Nf*Nf*4; ++f) tmpg[f] = prob->g[f];
454:   for (f = 0; f < Nf*Nf*4; ++f) tmpgp[f] = prob->gp[f];
455:   for (f = 0; f < Nf; ++f) tmpr[f] = prob->r[f];
456:   for (f = 0; f < Nf; ++f) tmpup[f] = prob->update[f];
457:   for (f = 0; f < Nf; ++f) tmpctx[f] = prob->ctx[f];
458:   for (f = Nf; f < NfNew; ++f) tmpobj[f] = NULL;
459:   for (f = Nf*2; f < NfNew*2; ++f) tmpf[f] = NULL;
460:   for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpg[f] = NULL;
461:   for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgp[f] = NULL;
462:   for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgt[f] = NULL;
463:   for (f = Nf; f < NfNew; ++f) tmpr[f] = NULL;
464:   for (f = Nf; f < NfNew; ++f) tmpup[f] = NULL;
465:   for (f = Nf; f < NfNew; ++f) tmpctx[f] = NULL;
466:   PetscFree7(prob->obj, prob->f, prob->g, prob->gp, prob->gt, prob->r, prob->ctx);
467:   PetscFree(prob->update);
468:   prob->obj = tmpobj;
469:   prob->f   = tmpf;
470:   prob->g   = tmpg;
471:   prob->gp  = tmpgp;
472:   prob->gt  = tmpgt;
473:   prob->r   = tmpr;
474:   prob->update = tmpup;
475:   prob->ctx = tmpctx;
476:   PetscCalloc4(NfNew*2, &tmpfbd, NfNew*NfNew*4, &tmpgbd, NfNew, &tmpexactSol, NfNew, &tmpexactCtx);
477:   for (f = 0; f < Nf*2; ++f) tmpfbd[f] = prob->fBd[f];
478:   for (f = 0; f < Nf*Nf*4; ++f) tmpgbd[f] = prob->gBd[f];
479:   for (f = 0; f < Nf; ++f) tmpexactSol[f] = prob->exactSol[f];
480:   for (f = 0; f < Nf; ++f) tmpexactCtx[f] = prob->exactCtx[f];
481:   for (f = Nf*2; f < NfNew*2; ++f) tmpfbd[f] = NULL;
482:   for (f = Nf*Nf*4; f < NfNew*NfNew*4; ++f) tmpgbd[f] = NULL;
483:   for (f = Nf; f < NfNew; ++f) tmpexactSol[f] = NULL;
484:   for (f = Nf; f < NfNew; ++f) tmpexactCtx[f] = NULL;
485:   PetscFree4(prob->fBd, prob->gBd, prob->exactSol, prob->exactCtx);
486:   prob->fBd = tmpfbd;
487:   prob->gBd = tmpgbd;
488:   prob->exactSol = tmpexactSol;
489:   prob->exactCtx = tmpexactCtx;
490:   return(0);
491: }

493: /*@
494:   PetscDSDestroy - Destroys a PetscDS object

496:   Collective on prob

498:   Input Parameter:
499: . prob - the PetscDS object to destroy

501:   Level: developer

503: .seealso PetscDSView()
504: @*/
505: PetscErrorCode PetscDSDestroy(PetscDS *prob)
506: {
507:   PetscInt       f;
508:   DSBoundary     next;

512:   if (!*prob) return(0);

515:   if (--((PetscObject)(*prob))->refct > 0) {*prob = 0; return(0);}
516:   ((PetscObject) (*prob))->refct = 0;
517:   if ((*prob)->subprobs) {
518:     PetscInt dim, d;

520:     PetscDSGetSpatialDimension(*prob, &dim);
521:     for (d = 0; d < dim; ++d) {PetscDSDestroy(&(*prob)->subprobs[d]);}
522:   }
523:   PetscFree((*prob)->subprobs);
524:   PetscDSDestroyStructs_Static(*prob);
525:   for (f = 0; f < (*prob)->Nf; ++f) {
526:     PetscObjectDereference((*prob)->disc[f]);
527:   }
528:   PetscFree2((*prob)->disc, (*prob)->implicit);
529:   PetscFree7((*prob)->obj,(*prob)->f,(*prob)->g,(*prob)->gp,(*prob)->gt,(*prob)->r,(*prob)->ctx);
530:   PetscFree((*prob)->update);
531:   PetscFree4((*prob)->fBd,(*prob)->gBd,(*prob)->exactSol,(*prob)->exactCtx);
532:   if ((*prob)->ops->destroy) {(*(*prob)->ops->destroy)(*prob);}
533:   next = (*prob)->boundary;
534:   while (next) {
535:     DSBoundary b = next;

537:     next = b->next;
538:     PetscFree(b->comps);
539:     PetscFree(b->ids);
540:     PetscFree(b->name);
541:     PetscFree(b->labelname);
542:     PetscFree(b);
543:   }
544:   PetscFree((*prob)->constants);
545:   PetscHeaderDestroy(prob);
546:   return(0);
547: }

549: /*@
550:   PetscDSCreate - Creates an empty PetscDS object. The type can then be set with PetscDSSetType().

552:   Collective

554:   Input Parameter:
555: . comm - The communicator for the PetscDS object

557:   Output Parameter:
558: . prob - The PetscDS object

560:   Level: beginner

562: .seealso: PetscDSSetType(), PETSCDSBASIC
563: @*/
564: PetscErrorCode PetscDSCreate(MPI_Comm comm, PetscDS *prob)
565: {
566:   PetscDS   p;

571:   *prob  = NULL;
572:   PetscDSInitializePackage();

574:   PetscHeaderCreate(p, PETSCDS_CLASSID, "PetscDS", "Discrete System", "PetscDS", comm, PetscDSDestroy, PetscDSView);

576:   p->Nf           = 0;
577:   p->setup        = PETSC_FALSE;
578:   p->numConstants = 0;
579:   p->constants    = NULL;
580:   p->dimEmbed     = -1;
581:   p->useJacPre    = PETSC_TRUE;

583:   *prob = p;
584:   return(0);
585: }

587: /*@
588:   PetscDSGetNumFields - Returns the number of fields in the DS

590:   Not collective

592:   Input Parameter:
593: . prob - The PetscDS object

595:   Output Parameter:
596: . Nf - The number of fields

598:   Level: beginner

600: .seealso: PetscDSGetSpatialDimension(), PetscDSCreate()
601: @*/
602: PetscErrorCode PetscDSGetNumFields(PetscDS prob, PetscInt *Nf)
603: {
607:   *Nf = prob->Nf;
608:   return(0);
609: }

611: /*@
612:   PetscDSGetSpatialDimension - Returns the spatial dimension of the DS, meaning the topological dimension of the discretizations

614:   Not collective

616:   Input Parameter:
617: . prob - The PetscDS object

619:   Output Parameter:
620: . dim - The spatial dimension

622:   Level: beginner

624: .seealso: PetscDSGetCoordinateDimension(), PetscDSGetNumFields(), PetscDSCreate()
625: @*/
626: PetscErrorCode PetscDSGetSpatialDimension(PetscDS prob, PetscInt *dim)
627: {

633:   *dim = 0;
634:   if (prob->Nf) {
635:     PetscObject  obj;
636:     PetscClassId id;

638:     PetscDSGetDiscretization(prob, 0, &obj);
639:     PetscObjectGetClassId(obj, &id);
640:     if (id == PETSCFE_CLASSID)      {PetscFEGetSpatialDimension((PetscFE) obj, dim);}
641:     else if (id == PETSCFV_CLASSID) {PetscFVGetSpatialDimension((PetscFV) obj, dim);}
642:     else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", 0);
643:   }
644:   return(0);
645: }

647: /*@
648:   PetscDSGetCoordinateDimension - Returns the coordinate dimension of the DS, meaning the dimension of the space into which the discretiaztions are embedded

650:   Not collective

652:   Input Parameter:
653: . prob - The PetscDS object

655:   Output Parameter:
656: . dimEmbed - The coordinate dimension

658:   Level: beginner

660: .seealso: PetscDSSetCoordinateDimension(), PetscDSGetSpatialDimension(), PetscDSGetNumFields(), PetscDSCreate()
661: @*/
662: PetscErrorCode PetscDSGetCoordinateDimension(PetscDS prob, PetscInt *dimEmbed)
663: {
667:   if (prob->dimEmbed < 0) SETERRQ(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONGSTATE, "No coordinate dimension set for this DS");
668:   *dimEmbed = prob->dimEmbed;
669:   return(0);
670: }

672: /*@
673:   PetscDSSetCoordinateDimension - Set the coordinate dimension of the DS, meaning the dimension of the space into which the discretiaztions are embedded

675:   Logically collective on prob

677:   Input Parameters:
678: + prob - The PetscDS object
679: - dimEmbed - The coordinate dimension

681:   Level: beginner

683: .seealso: PetscDSGetCoordinateDimension(), PetscDSGetSpatialDimension(), PetscDSGetNumFields(), PetscDSCreate()
684: @*/
685: PetscErrorCode PetscDSSetCoordinateDimension(PetscDS prob, PetscInt dimEmbed)
686: {
689:   if (dimEmbed < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Coordinate dimension must be non-negative, not %D", dimEmbed);
690:   prob->dimEmbed = dimEmbed;
691:   return(0);
692: }

694: /*@
695:   PetscDSGetHybrid - Returns the flag for a hybrid (cohesive) cell

697:   Not collective

699:   Input Parameter:
700: . prob - The PetscDS object

702:   Output Parameter:
703: . isHybrid - The flag

705:   Level: developer

707: .seealso: PetscDSSetHybrid(), PetscDSCreate()
708: @*/
709: PetscErrorCode PetscDSGetHybrid(PetscDS prob, PetscBool *isHybrid)
710: {
714:   *isHybrid = prob->isHybrid;
715:   return(0);
716: }

718: /*@
719:   PetscDSSetHybrid - Set the flag for a hybrid (cohesive) cell

721:   Not collective

723:   Input Parameters:
724: + prob - The PetscDS object
725: - isHybrid - The flag

727:   Level: developer

729: .seealso: PetscDSGetHybrid(), PetscDSCreate()
730: @*/
731: PetscErrorCode PetscDSSetHybrid(PetscDS prob, PetscBool isHybrid)
732: {
735:   prob->isHybrid = isHybrid;
736:   return(0);
737: }

739: /*@
740:   PetscDSGetTotalDimension - Returns the total size of the approximation space for this system

742:   Not collective

744:   Input Parameter:
745: . prob - The PetscDS object

747:   Output Parameter:
748: . dim - The total problem dimension

750:   Level: beginner

752: .seealso: PetscDSGetNumFields(), PetscDSCreate()
753: @*/
754: PetscErrorCode PetscDSGetTotalDimension(PetscDS prob, PetscInt *dim)
755: {

760:   PetscDSSetUp(prob);
762:   *dim = prob->totDim;
763:   return(0);
764: }

766: /*@
767:   PetscDSGetTotalComponents - Returns the total number of components in this system

769:   Not collective

771:   Input Parameter:
772: . prob - The PetscDS object

774:   Output Parameter:
775: . dim - The total number of components

777:   Level: beginner

779: .seealso: PetscDSGetNumFields(), PetscDSCreate()
780: @*/
781: PetscErrorCode PetscDSGetTotalComponents(PetscDS prob, PetscInt *Nc)
782: {

787:   PetscDSSetUp(prob);
789:   *Nc = prob->totComp;
790:   return(0);
791: }

793: /*@
794:   PetscDSGetDiscretization - Returns the discretization object for the given field

796:   Not collective

798:   Input Parameters:
799: + prob - The PetscDS object
800: - f - The field number

802:   Output Parameter:
803: . disc - The discretization object

805:   Level: beginner

807: .seealso: PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
808: @*/
809: PetscErrorCode PetscDSGetDiscretization(PetscDS prob, PetscInt f, PetscObject *disc)
810: {
814:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
815:   *disc = prob->disc[f];
816:   return(0);
817: }

819: /*@
820:   PetscDSSetDiscretization - Sets the discretization object for the given field

822:   Not collective

824:   Input Parameters:
825: + prob - The PetscDS object
826: . f - The field number
827: - disc - The discretization object

829:   Level: beginner

831: .seealso: PetscDSGetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
832: @*/
833: PetscErrorCode PetscDSSetDiscretization(PetscDS prob, PetscInt f, PetscObject disc)
834: {

840:   if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
841:   PetscDSEnlarge_Static(prob, f+1);
842:   PetscObjectDereference(prob->disc[f]);
843:   prob->disc[f] = disc;
844:   PetscObjectReference(disc);
845:   {
846:     PetscClassId id;

848:     PetscObjectGetClassId(disc, &id);
849:     if (id == PETSCFE_CLASSID) {
850:       PetscDSSetImplicit(prob, f, PETSC_TRUE);
851:     } else if (id == PETSCFV_CLASSID) {
852:       PetscDSSetImplicit(prob, f, PETSC_FALSE);
853:     }
854:   }
855:   return(0);
856: }

858: /*@
859:   PetscDSAddDiscretization - Adds a discretization object

861:   Not collective

863:   Input Parameters:
864: + prob - The PetscDS object
865: - disc - The boundary discretization object

867:   Level: beginner

869: .seealso: PetscDSGetDiscretization(), PetscDSSetDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
870: @*/
871: PetscErrorCode PetscDSAddDiscretization(PetscDS prob, PetscObject disc)
872: {

876:   PetscDSSetDiscretization(prob, prob->Nf, disc);
877:   return(0);
878: }

880: /*@
881:   PetscDSGetImplicit - Returns the flag for implicit solve for this field. This is just a guide for IMEX

883:   Not collective

885:   Input Parameters:
886: + prob - The PetscDS object
887: - f - The field number

889:   Output Parameter:
890: . implicit - The flag indicating what kind of solve to use for this field

892:   Level: developer

894: .seealso: PetscDSSetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
895: @*/
896: PetscErrorCode PetscDSGetImplicit(PetscDS prob, PetscInt f, PetscBool *implicit)
897: {
901:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
902:   *implicit = prob->implicit[f];
903:   return(0);
904: }

906: /*@
907:   PetscDSSetImplicit - Set the flag for implicit solve for this field. This is just a guide for IMEX

909:   Not collective

911:   Input Parameters:
912: + prob - The PetscDS object
913: . f - The field number
914: - implicit - The flag indicating what kind of solve to use for this field

916:   Level: developer

918: .seealso: PetscDSGetImplicit(), PetscDSSetDiscretization(), PetscDSAddDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
919: @*/
920: PetscErrorCode PetscDSSetImplicit(PetscDS prob, PetscInt f, PetscBool implicit)
921: {
924:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
925:   prob->implicit[f] = implicit;
926:   return(0);
927: }

929: PetscErrorCode PetscDSGetObjective(PetscDS prob, PetscInt f,
930:                                    void (**obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
931:                                                 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
932:                                                 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
933:                                                 PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
934: {
938:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
939:   *obj = prob->obj[f];
940:   return(0);
941: }

943: PetscErrorCode PetscDSSetObjective(PetscDS prob, PetscInt f,
944:                                    void (*obj)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
945:                                                const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
946:                                                const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
947:                                                PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar obj[]))
948: {

954:   if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
955:   PetscDSEnlarge_Static(prob, f+1);
956:   prob->obj[f] = obj;
957:   return(0);
958: }

960: /*@C
961:   PetscDSGetResidual - Get the pointwise residual function for a given test field

963:   Not collective

965:   Input Parameters:
966: + prob - The PetscDS
967: - f    - The test field number

969:   Output Parameters:
970: + f0 - integrand for the test function term
971: - f1 - integrand for the test function gradient term

973:   Note: We are using a first order FEM model for the weak form:

975:   \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t)

977: The calling sequence for the callbacks f0 and f1 is given by:

979: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
980: $    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
981: $    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
982: $    PetscReal t, const PetscReal x[], PetscScalar f0[])

984: + dim - the spatial dimension
985: . Nf - the number of fields
986: . uOff - the offset into u[] and u_t[] for each field
987: . uOff_x - the offset into u_x[] for each field
988: . u - each field evaluated at the current point
989: . u_t - the time derivative of each field evaluated at the current point
990: . u_x - the gradient of each field evaluated at the current point
991: . aOff - the offset into a[] and a_t[] for each auxiliary field
992: . aOff_x - the offset into a_x[] for each auxiliary field
993: . a - each auxiliary field evaluated at the current point
994: . a_t - the time derivative of each auxiliary field evaluated at the current point
995: . a_x - the gradient of auxiliary each field evaluated at the current point
996: . t - current time
997: . x - coordinates of the current point
998: . numConstants - number of constant parameters
999: . constants - constant parameters
1000: - f0 - output values at the current point

1002:   Level: intermediate

1004: .seealso: PetscDSSetResidual()
1005: @*/
1006: PetscErrorCode PetscDSGetResidual(PetscDS prob, PetscInt f,
1007:                                   void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1008:                                               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1009:                                               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1010:                                               PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1011:                                   void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1012:                                               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1013:                                               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1014:                                               PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1015: {
1018:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1021:   return(0);
1022: }

1024: /*@C
1025:   PetscDSSetResidual - Set the pointwise residual function for a given test field

1027:   Not collective

1029:   Input Parameters:
1030: + prob - The PetscDS
1031: . f    - The test field number
1032: . f0 - integrand for the test function term
1033: - f1 - integrand for the test function gradient term

1035:   Note: We are using a first order FEM model for the weak form:

1037:   \int_\Omega \phi f_0(u, u_t, \nabla u, x, t) + \nabla\phi \cdot {\vec f}_1(u, u_t, \nabla u, x, t)

1039: The calling sequence for the callbacks f0 and f1 is given by:

1041: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1042: $    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1043: $    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1044: $    PetscReal t, const PetscReal x[], PetscScalar f0[])

1046: + dim - the spatial dimension
1047: . Nf - the number of fields
1048: . uOff - the offset into u[] and u_t[] for each field
1049: . uOff_x - the offset into u_x[] for each field
1050: . u - each field evaluated at the current point
1051: . u_t - the time derivative of each field evaluated at the current point
1052: . u_x - the gradient of each field evaluated at the current point
1053: . aOff - the offset into a[] and a_t[] for each auxiliary field
1054: . aOff_x - the offset into a_x[] for each auxiliary field
1055: . a - each auxiliary field evaluated at the current point
1056: . a_t - the time derivative of each auxiliary field evaluated at the current point
1057: . a_x - the gradient of auxiliary each field evaluated at the current point
1058: . t - current time
1059: . x - coordinates of the current point
1060: . numConstants - number of constant parameters
1061: . constants - constant parameters
1062: - f0 - output values at the current point

1064:   Level: intermediate

1066: .seealso: PetscDSGetResidual()
1067: @*/
1068: PetscErrorCode PetscDSSetResidual(PetscDS prob, PetscInt f,
1069:                                   void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1070:                                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1071:                                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1072:                                              PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1073:                                   void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1074:                                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1075:                                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1076:                                              PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1077: {

1084:   if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1085:   PetscDSEnlarge_Static(prob, f+1);
1086:   prob->f[f*2+0] = f0;
1087:   prob->f[f*2+1] = f1;
1088:   return(0);
1089: }

1091: /*@C
1092:   PetscDSHasJacobian - Signals that Jacobian functions have been set

1094:   Not collective

1096:   Input Parameter:
1097: . prob - The PetscDS

1099:   Output Parameter:
1100: . hasJac - flag that pointwise function for the Jacobian has been set

1102:   Level: intermediate

1104: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1105: @*/
1106: PetscErrorCode PetscDSHasJacobian(PetscDS prob, PetscBool *hasJac)
1107: {
1108:   PetscInt f, g, h;

1112:   *hasJac = PETSC_FALSE;
1113:   for (f = 0; f < prob->Nf; ++f) {
1114:     for (g = 0; g < prob->Nf; ++g) {
1115:       for (h = 0; h < 4; ++h) {
1116:         if (prob->g[(f*prob->Nf + g)*4+h]) *hasJac = PETSC_TRUE;
1117:       }
1118:     }
1119:   }
1120:   return(0);
1121: }

1123: /*@C
1124:   PetscDSGetJacobian - Get the pointwise Jacobian function for given test and basis field

1126:   Not collective

1128:   Input Parameters:
1129: + prob - The PetscDS
1130: . f    - The test field number
1131: - g    - The field number

1133:   Output Parameters:
1134: + g0 - integrand for the test and basis function term
1135: . g1 - integrand for the test function and basis function gradient term
1136: . g2 - integrand for the test function gradient and basis function term
1137: - g3 - integrand for the test function gradient and basis function gradient term

1139:   Note: We are using a first order FEM model for the weak form:

1141:   \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi

1143: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:

1145: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1146: $    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1147: $    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1148: $    PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])

1150: + dim - the spatial dimension
1151: . Nf - the number of fields
1152: . uOff - the offset into u[] and u_t[] for each field
1153: . uOff_x - the offset into u_x[] for each field
1154: . u - each field evaluated at the current point
1155: . u_t - the time derivative of each field evaluated at the current point
1156: . u_x - the gradient of each field evaluated at the current point
1157: . aOff - the offset into a[] and a_t[] for each auxiliary field
1158: . aOff_x - the offset into a_x[] for each auxiliary field
1159: . a - each auxiliary field evaluated at the current point
1160: . a_t - the time derivative of each auxiliary field evaluated at the current point
1161: . a_x - the gradient of auxiliary each field evaluated at the current point
1162: . t - current time
1163: . u_tShift - the multiplier a for dF/dU_t
1164: . x - coordinates of the current point
1165: . numConstants - number of constant parameters
1166: . constants - constant parameters
1167: - g0 - output values at the current point

1169:   Level: intermediate

1171: .seealso: PetscDSSetJacobian()
1172: @*/
1173: PetscErrorCode PetscDSGetJacobian(PetscDS prob, PetscInt f, PetscInt g,
1174:                                   void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1175:                                               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1176:                                               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1177:                                               PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1178:                                   void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1179:                                               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1180:                                               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1181:                                               PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1182:                                   void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1183:                                               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1184:                                               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1185:                                               PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1186:                                   void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1187:                                               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1188:                                               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1189:                                               PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1190: {
1193:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1194:   if ((g < 0) || (g >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", g, prob->Nf);
1199:   return(0);
1200: }

1202: /*@C
1203:   PetscDSSetJacobian - Set the pointwise Jacobian function for given test and basis fields

1205:   Not collective

1207:   Input Parameters:
1208: + prob - The PetscDS
1209: . f    - The test field number
1210: . g    - The field number
1211: . g0 - integrand for the test and basis function term
1212: . g1 - integrand for the test function and basis function gradient term
1213: . g2 - integrand for the test function gradient and basis function term
1214: - g3 - integrand for the test function gradient and basis function gradient term

1216:   Note: We are using a first order FEM model for the weak form:

1218:   \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi

1220: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:

1222: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1223: $    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1224: $    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1225: $    PetscReal t, const PetscReal x[], PetscScalar g0[])

1227: + dim - the spatial dimension
1228: . Nf - the number of fields
1229: . uOff - the offset into u[] and u_t[] for each field
1230: . uOff_x - the offset into u_x[] for each field
1231: . u - each field evaluated at the current point
1232: . u_t - the time derivative of each field evaluated at the current point
1233: . u_x - the gradient of each field evaluated at the current point
1234: . aOff - the offset into a[] and a_t[] for each auxiliary field
1235: . aOff_x - the offset into a_x[] for each auxiliary field
1236: . a - each auxiliary field evaluated at the current point
1237: . a_t - the time derivative of each auxiliary field evaluated at the current point
1238: . a_x - the gradient of auxiliary each field evaluated at the current point
1239: . t - current time
1240: . u_tShift - the multiplier a for dF/dU_t
1241: . x - coordinates of the current point
1242: . numConstants - number of constant parameters
1243: . constants - constant parameters
1244: - g0 - output values at the current point

1246:   Level: intermediate

1248: .seealso: PetscDSGetJacobian()
1249: @*/
1250: PetscErrorCode PetscDSSetJacobian(PetscDS prob, PetscInt f, PetscInt g,
1251:                                   void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1252:                                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1253:                                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1254:                                              PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1255:                                   void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1256:                                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1257:                                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1258:                                              PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1259:                                   void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1260:                                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1261:                                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1262:                                              PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1263:                                   void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1264:                                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1265:                                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1266:                                              PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1267: {

1276:   if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1277:   if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1278:   PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1279:   prob->g[(f*prob->Nf + g)*4+0] = g0;
1280:   prob->g[(f*prob->Nf + g)*4+1] = g1;
1281:   prob->g[(f*prob->Nf + g)*4+2] = g2;
1282:   prob->g[(f*prob->Nf + g)*4+3] = g3;
1283:   return(0);
1284: }

1286: /*@C
1287:   PetscDSUseJacobianPreconditioner - Whether to construct a Jacobian preconditioner

1289:   Not collective

1291:   Input Parameters:
1292: + prob - The PetscDS
1293: - useJacPre - flag that enables construction of a Jacobian preconditioner

1295:   Level: intermediate

1297: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1298: @*/
1299: PetscErrorCode PetscDSUseJacobianPreconditioner(PetscDS prob, PetscBool useJacPre)
1300: {
1303:   prob->useJacPre = useJacPre;
1304:   return(0);
1305: }

1307: /*@C
1308:   PetscDSHasJacobianPreconditioner - Signals that a Jacobian preconditioner matrix has been set

1310:   Not collective

1312:   Input Parameter:
1313: . prob - The PetscDS

1315:   Output Parameter:
1316: . hasJacPre - flag that pointwise function for Jacobian preconditioner matrix has been set

1318:   Level: intermediate

1320: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1321: @*/
1322: PetscErrorCode PetscDSHasJacobianPreconditioner(PetscDS prob, PetscBool *hasJacPre)
1323: {
1324:   PetscInt f, g, h;

1328:   *hasJacPre = PETSC_FALSE;
1329:   if (!prob->useJacPre) return(0);
1330:   for (f = 0; f < prob->Nf; ++f) {
1331:     for (g = 0; g < prob->Nf; ++g) {
1332:       for (h = 0; h < 4; ++h) {
1333:         if (prob->gp[(f*prob->Nf + g)*4+h]) *hasJacPre = PETSC_TRUE;
1334:       }
1335:     }
1336:   }
1337:   return(0);
1338: }

1340: /*@C
1341:   PetscDSGetJacobianPreconditioner - Get the pointwise Jacobian preconditioner function for given test and basis field. If this is missing, the system matrix is used to build the preconditioner.

1343:   Not collective

1345:   Input Parameters:
1346: + prob - The PetscDS
1347: . f    - The test field number
1348: - g    - The field number

1350:   Output Parameters:
1351: + g0 - integrand for the test and basis function term
1352: . g1 - integrand for the test function and basis function gradient term
1353: . g2 - integrand for the test function gradient and basis function term
1354: - g3 - integrand for the test function gradient and basis function gradient term

1356:   Note: We are using a first order FEM model for the weak form:

1358:   \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi

1360: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:

1362: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1363: $    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1364: $    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1365: $    PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])

1367: + dim - the spatial dimension
1368: . Nf - the number of fields
1369: . uOff - the offset into u[] and u_t[] for each field
1370: . uOff_x - the offset into u_x[] for each field
1371: . u - each field evaluated at the current point
1372: . u_t - the time derivative of each field evaluated at the current point
1373: . u_x - the gradient of each field evaluated at the current point
1374: . aOff - the offset into a[] and a_t[] for each auxiliary field
1375: . aOff_x - the offset into a_x[] for each auxiliary field
1376: . a - each auxiliary field evaluated at the current point
1377: . a_t - the time derivative of each auxiliary field evaluated at the current point
1378: . a_x - the gradient of auxiliary each field evaluated at the current point
1379: . t - current time
1380: . u_tShift - the multiplier a for dF/dU_t
1381: . x - coordinates of the current point
1382: . numConstants - number of constant parameters
1383: . constants - constant parameters
1384: - g0 - output values at the current point

1386:   Level: intermediate

1388: .seealso: PetscDSSetJacobianPreconditioner(), PetscDSGetJacobian()
1389: @*/
1390: PetscErrorCode PetscDSGetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g,
1391:                                   void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1392:                                               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1393:                                               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1394:                                               PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1395:                                   void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1396:                                               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1397:                                               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1398:                                               PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1399:                                   void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1400:                                               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1401:                                               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1402:                                               PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1403:                                   void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1404:                                               const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1405:                                               const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1406:                                               PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1407: {
1410:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1411:   if ((g < 0) || (g >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", g, prob->Nf);
1416:   return(0);
1417: }

1419: /*@C
1420:   PetscDSSetJacobianPreconditioner - Set the pointwise Jacobian preconditioner function for given test and basis fields. If this is missing, the system matrix is used to build the preconditioner.

1422:   Not collective

1424:   Input Parameters:
1425: + prob - The PetscDS
1426: . f    - The test field number
1427: . g    - The field number
1428: . g0 - integrand for the test and basis function term
1429: . g1 - integrand for the test function and basis function gradient term
1430: . g2 - integrand for the test function gradient and basis function term
1431: - g3 - integrand for the test function gradient and basis function gradient term

1433:   Note: We are using a first order FEM model for the weak form:

1435:   \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi

1437: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:

1439: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1440: $    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1441: $    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1442: $    PetscReal t, const PetscReal x[], PetscScalar g0[])

1444: + dim - the spatial dimension
1445: . Nf - the number of fields
1446: . uOff - the offset into u[] and u_t[] for each field
1447: . uOff_x - the offset into u_x[] for each field
1448: . u - each field evaluated at the current point
1449: . u_t - the time derivative of each field evaluated at the current point
1450: . u_x - the gradient of each field evaluated at the current point
1451: . aOff - the offset into a[] and a_t[] for each auxiliary field
1452: . aOff_x - the offset into a_x[] for each auxiliary field
1453: . a - each auxiliary field evaluated at the current point
1454: . a_t - the time derivative of each auxiliary field evaluated at the current point
1455: . a_x - the gradient of auxiliary each field evaluated at the current point
1456: . t - current time
1457: . u_tShift - the multiplier a for dF/dU_t
1458: . x - coordinates of the current point
1459: . numConstants - number of constant parameters
1460: . constants - constant parameters
1461: - g0 - output values at the current point

1463:   Level: intermediate

1465: .seealso: PetscDSGetJacobianPreconditioner(), PetscDSSetJacobian()
1466: @*/
1467: PetscErrorCode PetscDSSetJacobianPreconditioner(PetscDS prob, PetscInt f, PetscInt g,
1468:                                   void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1469:                                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1470:                                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1471:                                              PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1472:                                   void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1473:                                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1474:                                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1475:                                              PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1476:                                   void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1477:                                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1478:                                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1479:                                              PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1480:                                   void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1481:                                              const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1482:                                              const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1483:                                              PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1484: {

1493:   if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1494:   if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1495:   PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1496:   prob->gp[(f*prob->Nf + g)*4+0] = g0;
1497:   prob->gp[(f*prob->Nf + g)*4+1] = g1;
1498:   prob->gp[(f*prob->Nf + g)*4+2] = g2;
1499:   prob->gp[(f*prob->Nf + g)*4+3] = g3;
1500:   return(0);
1501: }

1503: /*@C
1504:   PetscDSHasDynamicJacobian - Signals that a dynamic Jacobian, dF/du_t, has been set

1506:   Not collective

1508:   Input Parameter:
1509: . prob - The PetscDS

1511:   Output Parameter:
1512: . hasDynJac - flag that pointwise function for dynamic Jacobian has been set

1514:   Level: intermediate

1516: .seealso: PetscDSGetDynamicJacobian(), PetscDSSetDynamicJacobian(), PetscDSGetJacobian()
1517: @*/
1518: PetscErrorCode PetscDSHasDynamicJacobian(PetscDS prob, PetscBool *hasDynJac)
1519: {
1520:   PetscInt f, g, h;

1524:   *hasDynJac = PETSC_FALSE;
1525:   for (f = 0; f < prob->Nf; ++f) {
1526:     for (g = 0; g < prob->Nf; ++g) {
1527:       for (h = 0; h < 4; ++h) {
1528:         if (prob->gt[(f*prob->Nf + g)*4+h]) *hasDynJac = PETSC_TRUE;
1529:       }
1530:     }
1531:   }
1532:   return(0);
1533: }

1535: /*@C
1536:   PetscDSGetDynamicJacobian - Get the pointwise dynamic Jacobian, dF/du_t, function for given test and basis field

1538:   Not collective

1540:   Input Parameters:
1541: + prob - The PetscDS
1542: . f    - The test field number
1543: - g    - The field number

1545:   Output Parameters:
1546: + g0 - integrand for the test and basis function term
1547: . g1 - integrand for the test function and basis function gradient term
1548: . g2 - integrand for the test function gradient and basis function term
1549: - g3 - integrand for the test function gradient and basis function gradient term

1551:   Note: We are using a first order FEM model for the weak form:

1553:   \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi

1555: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:

1557: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1558: $    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1559: $    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1560: $    PetscReal t, const PetscReal u_tShift, const PetscReal x[], PetscScalar g0[])

1562: + dim - the spatial dimension
1563: . Nf - the number of fields
1564: . uOff - the offset into u[] and u_t[] for each field
1565: . uOff_x - the offset into u_x[] for each field
1566: . u - each field evaluated at the current point
1567: . u_t - the time derivative of each field evaluated at the current point
1568: . u_x - the gradient of each field evaluated at the current point
1569: . aOff - the offset into a[] and a_t[] for each auxiliary field
1570: . aOff_x - the offset into a_x[] for each auxiliary field
1571: . a - each auxiliary field evaluated at the current point
1572: . a_t - the time derivative of each auxiliary field evaluated at the current point
1573: . a_x - the gradient of auxiliary each field evaluated at the current point
1574: . t - current time
1575: . u_tShift - the multiplier a for dF/dU_t
1576: . x - coordinates of the current point
1577: . numConstants - number of constant parameters
1578: . constants - constant parameters
1579: - g0 - output values at the current point

1581:   Level: intermediate

1583: .seealso: PetscDSSetJacobian()
1584: @*/
1585: PetscErrorCode PetscDSGetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g,
1586:                                          void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1587:                                                      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1588:                                                      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1589:                                                      PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1590:                                          void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1591:                                                      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1592:                                                      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1593:                                                      PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1594:                                          void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1595:                                                      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1596:                                                      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1597:                                                      PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1598:                                          void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1599:                                                      const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1600:                                                      const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1601:                                                      PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1602: {
1605:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1606:   if ((g < 0) || (g >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", g, prob->Nf);
1611:   return(0);
1612: }

1614: /*@C
1615:   PetscDSSetDynamicJacobian - Set the pointwise dynamic Jacobian, dF/du_t, function for given test and basis fields

1617:   Not collective

1619:   Input Parameters:
1620: + prob - The PetscDS
1621: . f    - The test field number
1622: . g    - The field number
1623: . g0 - integrand for the test and basis function term
1624: . g1 - integrand for the test function and basis function gradient term
1625: . g2 - integrand for the test function gradient and basis function term
1626: - g3 - integrand for the test function gradient and basis function gradient term

1628:   Note: We are using a first order FEM model for the weak form:

1630:   \int_\Omega \phi g_0(u, u_t, \nabla u, x, t) \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \nabla \psi

1632: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:

1634: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1635: $    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1636: $    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1637: $    PetscReal t, const PetscReal x[], PetscScalar g0[])

1639: + dim - the spatial dimension
1640: . Nf - the number of fields
1641: . uOff - the offset into u[] and u_t[] for each field
1642: . uOff_x - the offset into u_x[] for each field
1643: . u - each field evaluated at the current point
1644: . u_t - the time derivative of each field evaluated at the current point
1645: . u_x - the gradient of each field evaluated at the current point
1646: . aOff - the offset into a[] and a_t[] for each auxiliary field
1647: . aOff_x - the offset into a_x[] for each auxiliary field
1648: . a - each auxiliary field evaluated at the current point
1649: . a_t - the time derivative of each auxiliary field evaluated at the current point
1650: . a_x - the gradient of auxiliary each field evaluated at the current point
1651: . t - current time
1652: . u_tShift - the multiplier a for dF/dU_t
1653: . x - coordinates of the current point
1654: . numConstants - number of constant parameters
1655: . constants - constant parameters
1656: - g0 - output values at the current point

1658:   Level: intermediate

1660: .seealso: PetscDSGetJacobian()
1661: @*/
1662: PetscErrorCode PetscDSSetDynamicJacobian(PetscDS prob, PetscInt f, PetscInt g,
1663:                                          void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1664:                                                     const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1665:                                                     const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1666:                                                     PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
1667:                                          void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1668:                                                     const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1669:                                                     const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1670:                                                     PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
1671:                                          void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1672:                                                     const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1673:                                                     const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1674:                                                     PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
1675:                                          void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1676:                                                     const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1677:                                                     const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1678:                                                     PetscReal t, PetscReal u_tShift, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
1679: {

1688:   if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1689:   if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
1690:   PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
1691:   prob->gt[(f*prob->Nf + g)*4+0] = g0;
1692:   prob->gt[(f*prob->Nf + g)*4+1] = g1;
1693:   prob->gt[(f*prob->Nf + g)*4+2] = g2;
1694:   prob->gt[(f*prob->Nf + g)*4+3] = g3;
1695:   return(0);
1696: }

1698: /*@C
1699:   PetscDSGetRiemannSolver - Returns the Riemann solver for the given field

1701:   Not collective

1703:   Input Arguments:
1704: + prob - The PetscDS object
1705: - f    - The field number

1707:   Output Argument:
1708: . r    - Riemann solver

1710:   Calling sequence for r:

1712: $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)

1714: + dim  - The spatial dimension
1715: . Nf   - The number of fields
1716: . x    - The coordinates at a point on the interface
1717: . n    - The normal vector to the interface
1718: . uL   - The state vector to the left of the interface
1719: . uR   - The state vector to the right of the interface
1720: . flux - output array of flux through the interface
1721: . numConstants - number of constant parameters
1722: . constants - constant parameters
1723: - ctx  - optional user context

1725:   Level: intermediate

1727: .seealso: PetscDSSetRiemannSolver()
1728: @*/
1729: PetscErrorCode PetscDSGetRiemannSolver(PetscDS prob, PetscInt f,
1730:                                        void (**r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscInt numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx))
1731: {
1734:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1736:   *r = prob->r[f];
1737:   return(0);
1738: }

1740: /*@C
1741:   PetscDSSetRiemannSolver - Sets the Riemann solver for the given field

1743:   Not collective

1745:   Input Arguments:
1746: + prob - The PetscDS object
1747: . f    - The field number
1748: - r    - Riemann solver

1750:   Calling sequence for r:

1752: $ r(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscScalar flux[], void *ctx)

1754: + dim  - The spatial dimension
1755: . Nf   - The number of fields
1756: . x    - The coordinates at a point on the interface
1757: . n    - The normal vector to the interface
1758: . uL   - The state vector to the left of the interface
1759: . uR   - The state vector to the right of the interface
1760: . flux - output array of flux through the interface
1761: . numConstants - number of constant parameters
1762: . constants - constant parameters
1763: - ctx  - optional user context

1765:   Level: intermediate

1767: .seealso: PetscDSGetRiemannSolver()
1768: @*/
1769: PetscErrorCode PetscDSSetRiemannSolver(PetscDS prob, PetscInt f,
1770:                                        void (*r)(PetscInt dim, PetscInt Nf, const PetscReal x[], const PetscReal n[], const PetscScalar uL[], const PetscScalar uR[], PetscInt numConstants, const PetscScalar constants[], PetscScalar flux[], void *ctx))
1771: {

1777:   if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1778:   PetscDSEnlarge_Static(prob, f+1);
1779:   prob->r[f] = r;
1780:   return(0);
1781: }

1783: /*@C
1784:   PetscDSGetUpdate - Get the pointwise update function for a given field

1786:   Not collective

1788:   Input Parameters:
1789: + prob - The PetscDS
1790: - f    - The field number

1792:   Output Parameters:
1793: . update - update function

1795:   Note: The calling sequence for the callback update is given by:

1797: $ update(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1798: $        const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1799: $        const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1800: $        PetscReal t, const PetscReal x[], PetscScalar uNew[])

1802: + dim - the spatial dimension
1803: . Nf - the number of fields
1804: . uOff - the offset into u[] and u_t[] for each field
1805: . uOff_x - the offset into u_x[] for each field
1806: . u - each field evaluated at the current point
1807: . u_t - the time derivative of each field evaluated at the current point
1808: . u_x - the gradient of each field evaluated at the current point
1809: . aOff - the offset into a[] and a_t[] for each auxiliary field
1810: . aOff_x - the offset into a_x[] for each auxiliary field
1811: . a - each auxiliary field evaluated at the current point
1812: . a_t - the time derivative of each auxiliary field evaluated at the current point
1813: . a_x - the gradient of auxiliary each field evaluated at the current point
1814: . t - current time
1815: . x - coordinates of the current point
1816: - uNew - new value for field at the current point

1818:   Level: intermediate

1820: .seealso: PetscDSSetUpdate(), PetscDSSetResidual()
1821: @*/
1822: PetscErrorCode PetscDSGetUpdate(PetscDS prob, PetscInt f,
1823:                                   void (**update)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1824:                                                   const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1825:                                                   const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1826:                                                   PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
1827: {
1830:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1832:   return(0);
1833: }

1835: /*@C
1836:   PetscDSSetUpdate - Set the pointwise update function for a given field

1838:   Not collective

1840:   Input Parameters:
1841: + prob   - The PetscDS
1842: . f      - The field number
1843: - update - update function

1845:   Note: The calling sequence for the callback update is given by:

1847: $ update(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1848: $        const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1849: $        const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1850: $        PetscReal t, const PetscReal x[], PetscScalar uNew[])

1852: + dim - the spatial dimension
1853: . Nf - the number of fields
1854: . uOff - the offset into u[] and u_t[] for each field
1855: . uOff_x - the offset into u_x[] for each field
1856: . u - each field evaluated at the current point
1857: . u_t - the time derivative of each field evaluated at the current point
1858: . u_x - the gradient of each field evaluated at the current point
1859: . aOff - the offset into a[] and a_t[] for each auxiliary field
1860: . aOff_x - the offset into a_x[] for each auxiliary field
1861: . a - each auxiliary field evaluated at the current point
1862: . a_t - the time derivative of each auxiliary field evaluated at the current point
1863: . a_x - the gradient of auxiliary each field evaluated at the current point
1864: . t - current time
1865: . x - coordinates of the current point
1866: - uNew - new field values at the current point

1868:   Level: intermediate

1870: .seealso: PetscDSGetResidual()
1871: @*/
1872: PetscErrorCode PetscDSSetUpdate(PetscDS prob, PetscInt f,
1873:                                 void (*update)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1874:                                                const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1875:                                                const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1876:                                                PetscReal t, const PetscReal x[], PetscInt numConstants, const PetscScalar constants[], PetscScalar uNew[]))
1877: {

1883:   if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1884:   PetscDSEnlarge_Static(prob, f+1);
1885:   prob->update[f] = update;
1886:   return(0);
1887: }

1889: PetscErrorCode PetscDSGetContext(PetscDS prob, PetscInt f, void **ctx)
1890: {
1893:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1895:   *ctx = prob->ctx[f];
1896:   return(0);
1897: }

1899: PetscErrorCode PetscDSSetContext(PetscDS prob, PetscInt f, void *ctx)
1900: {

1905:   if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
1906:   PetscDSEnlarge_Static(prob, f+1);
1907:   prob->ctx[f] = ctx;
1908:   return(0);
1909: }

1911: /*@C
1912:   PetscDSGetBdResidual - Get the pointwise boundary residual function for a given test field

1914:   Not collective

1916:   Input Parameters:
1917: + prob - The PetscDS
1918: - f    - The test field number

1920:   Output Parameters:
1921: + f0 - boundary integrand for the test function term
1922: - f1 - boundary integrand for the test function gradient term

1924:   Note: We are using a first order FEM model for the weak form:

1926:   \int_\Gamma \phi {\vec f}_0(u, u_t, \nabla u, x, t) \cdot \hat n + \nabla\phi \cdot {\overleftrightarrow f}_1(u, u_t, \nabla u, x, t) \cdot \hat n

1928: The calling sequence for the callbacks f0 and f1 is given by:

1930: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1931: $    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1932: $    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1933: $    PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[])

1935: + dim - the spatial dimension
1936: . Nf - the number of fields
1937: . uOff - the offset into u[] and u_t[] for each field
1938: . uOff_x - the offset into u_x[] for each field
1939: . u - each field evaluated at the current point
1940: . u_t - the time derivative of each field evaluated at the current point
1941: . u_x - the gradient of each field evaluated at the current point
1942: . aOff - the offset into a[] and a_t[] for each auxiliary field
1943: . aOff_x - the offset into a_x[] for each auxiliary field
1944: . a - each auxiliary field evaluated at the current point
1945: . a_t - the time derivative of each auxiliary field evaluated at the current point
1946: . a_x - the gradient of auxiliary each field evaluated at the current point
1947: . t - current time
1948: . x - coordinates of the current point
1949: . n - unit normal at the current point
1950: . numConstants - number of constant parameters
1951: . constants - constant parameters
1952: - f0 - output values at the current point

1954:   Level: intermediate

1956: .seealso: PetscDSSetBdResidual()
1957: @*/
1958: PetscErrorCode PetscDSGetBdResidual(PetscDS prob, PetscInt f,
1959:                                     void (**f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1960:                                                 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1961:                                                 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1962:                                                 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
1963:                                     void (**f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1964:                                                 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1965:                                                 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1966:                                                 PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
1967: {
1970:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
1973:   return(0);
1974: }

1976: /*@C
1977:   PetscDSSetBdResidual - Get the pointwise boundary residual function for a given test field

1979:   Not collective

1981:   Input Parameters:
1982: + prob - The PetscDS
1983: . f    - The test field number
1984: . f0 - boundary integrand for the test function term
1985: - f1 - boundary integrand for the test function gradient term

1987:   Note: We are using a first order FEM model for the weak form:

1989:   \int_\Gamma \phi {\vec f}_0(u, u_t, \nabla u, x, t) \cdot \hat n + \nabla\phi \cdot {\overleftrightarrow f}_1(u, u_t, \nabla u, x, t) \cdot \hat n

1991: The calling sequence for the callbacks f0 and f1 is given by:

1993: $ f0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
1994: $    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
1995: $    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
1996: $    PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar f0[])

1998: + dim - the spatial dimension
1999: . Nf - the number of fields
2000: . uOff - the offset into u[] and u_t[] for each field
2001: . uOff_x - the offset into u_x[] for each field
2002: . u - each field evaluated at the current point
2003: . u_t - the time derivative of each field evaluated at the current point
2004: . u_x - the gradient of each field evaluated at the current point
2005: . aOff - the offset into a[] and a_t[] for each auxiliary field
2006: . aOff_x - the offset into a_x[] for each auxiliary field
2007: . a - each auxiliary field evaluated at the current point
2008: . a_t - the time derivative of each auxiliary field evaluated at the current point
2009: . a_x - the gradient of auxiliary each field evaluated at the current point
2010: . t - current time
2011: . x - coordinates of the current point
2012: . n - unit normal at the current point
2013: . numConstants - number of constant parameters
2014: . constants - constant parameters
2015: - f0 - output values at the current point

2017:   Level: intermediate

2019: .seealso: PetscDSGetBdResidual()
2020: @*/
2021: PetscErrorCode PetscDSSetBdResidual(PetscDS prob, PetscInt f,
2022:                                     void (*f0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2023:                                                const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2024:                                                const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2025:                                                PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f0[]),
2026:                                     void (*f1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2027:                                                const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2028:                                                const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2029:                                                PetscReal t, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar f1[]))
2030: {

2035:   if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2036:   PetscDSEnlarge_Static(prob, f+1);
2039:   return(0);
2040: }

2042: /*@C
2043:   PetscDSGetBdJacobian - Get the pointwise boundary Jacobian function for given test and basis field

2045:   Not collective

2047:   Input Parameters:
2048: + prob - The PetscDS
2049: . f    - The test field number
2050: - g    - The field number

2052:   Output Parameters:
2053: + g0 - integrand for the test and basis function term
2054: . g1 - integrand for the test function and basis function gradient term
2055: . g2 - integrand for the test function gradient and basis function term
2056: - g3 - integrand for the test function gradient and basis function gradient term

2058:   Note: We are using a first order FEM model for the weak form:

2060:   \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi

2062: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:

2064: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2065: $    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2066: $    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2067: $    PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[])

2069: + dim - the spatial dimension
2070: . Nf - the number of fields
2071: . uOff - the offset into u[] and u_t[] for each field
2072: . uOff_x - the offset into u_x[] for each field
2073: . u - each field evaluated at the current point
2074: . u_t - the time derivative of each field evaluated at the current point
2075: . u_x - the gradient of each field evaluated at the current point
2076: . aOff - the offset into a[] and a_t[] for each auxiliary field
2077: . aOff_x - the offset into a_x[] for each auxiliary field
2078: . a - each auxiliary field evaluated at the current point
2079: . a_t - the time derivative of each auxiliary field evaluated at the current point
2080: . a_x - the gradient of auxiliary each field evaluated at the current point
2081: . t - current time
2082: . u_tShift - the multiplier a for dF/dU_t
2083: . x - coordinates of the current point
2084: . n - normal at the current point
2085: . numConstants - number of constant parameters
2086: . constants - constant parameters
2087: - g0 - output values at the current point

2089:   Level: intermediate

2091: .seealso: PetscDSSetBdJacobian()
2092: @*/
2093: PetscErrorCode PetscDSGetBdJacobian(PetscDS prob, PetscInt f, PetscInt g,
2094:                                     void (**g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2095:                                                 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2096:                                                 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2097:                                                 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
2098:                                     void (**g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2099:                                                 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2100:                                                 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2101:                                                 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
2102:                                     void (**g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2103:                                                 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2104:                                                 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2105:                                                 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
2106:                                     void (**g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2107:                                                 const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2108:                                                 const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2109:                                                 PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2110: {
2113:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
2114:   if ((g < 0) || (g >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", g, prob->Nf);
2119:   return(0);
2120: }

2122: /*@C
2123:   PetscDSSetBdJacobian - Set the pointwise boundary Jacobian function for given test and basis field

2125:   Not collective

2127:   Input Parameters:
2128: + prob - The PetscDS
2129: . f    - The test field number
2130: . g    - The field number
2131: . g0 - integrand for the test and basis function term
2132: . g1 - integrand for the test function and basis function gradient term
2133: . g2 - integrand for the test function gradient and basis function term
2134: - g3 - integrand for the test function gradient and basis function gradient term

2136:   Note: We are using a first order FEM model for the weak form:

2138:   \int_\Gamma \phi {\vec g}_0(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \phi {\vec g}_1(u, u_t, \nabla u, x, t) \cdot \hat n \nabla \psi + \nabla\phi \cdot {\vec g}_2(u, u_t, \nabla u, x, t) \cdot \hat n \psi + \nabla\phi \cdot {\overleftrightarrow g}_3(u, u_t, \nabla u, x, t) \cdot \hat n \cdot \nabla \psi

2140: The calling sequence for the callbacks g0, g1, g2 and g3 is given by:

2142: $ g0(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2143: $    const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2144: $    const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2145: $    PetscReal t, const PetscReal x[], const PetscReal n[], PetscScalar g0[])

2147: + dim - the spatial dimension
2148: . Nf - the number of fields
2149: . uOff - the offset into u[] and u_t[] for each field
2150: . uOff_x - the offset into u_x[] for each field
2151: . u - each field evaluated at the current point
2152: . u_t - the time derivative of each field evaluated at the current point
2153: . u_x - the gradient of each field evaluated at the current point
2154: . aOff - the offset into a[] and a_t[] for each auxiliary field
2155: . aOff_x - the offset into a_x[] for each auxiliary field
2156: . a - each auxiliary field evaluated at the current point
2157: . a_t - the time derivative of each auxiliary field evaluated at the current point
2158: . a_x - the gradient of auxiliary each field evaluated at the current point
2159: . t - current time
2160: . u_tShift - the multiplier a for dF/dU_t
2161: . x - coordinates of the current point
2162: . n - normal at the current point
2163: . numConstants - number of constant parameters
2164: . constants - constant parameters
2165: - g0 - output values at the current point

2167:   Level: intermediate

2169: .seealso: PetscDSGetBdJacobian()
2170: @*/
2171: PetscErrorCode PetscDSSetBdJacobian(PetscDS prob, PetscInt f, PetscInt g,
2172:                                     void (*g0)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2173:                                                const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2174:                                                const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2175:                                                PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g0[]),
2176:                                     void (*g1)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2177:                                                const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2178:                                                const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2179:                                                PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g1[]),
2180:                                     void (*g2)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2181:                                                const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2182:                                                const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2183:                                                PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g2[]),
2184:                                     void (*g3)(PetscInt dim, PetscInt Nf, PetscInt NfAux,
2185:                                                const PetscInt uOff[], const PetscInt uOff_x[], const PetscScalar u[], const PetscScalar u_t[], const PetscScalar u_x[],
2186:                                                const PetscInt aOff[], const PetscInt aOff_x[], const PetscScalar a[], const PetscScalar a_t[], const PetscScalar a_x[],
2187:                                                PetscReal t, PetscReal u_tShift, const PetscReal x[], const PetscReal n[], PetscInt numConstants, const PetscScalar constants[], PetscScalar g3[]))
2188: {

2197:   if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2198:   if (g < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", g);
2199:   PetscDSEnlarge_Static(prob, PetscMax(f, g)+1);
2200:   prob->gBd[(f*prob->Nf + g)*4+0] = g0;
2201:   prob->gBd[(f*prob->Nf + g)*4+1] = g1;
2202:   prob->gBd[(f*prob->Nf + g)*4+2] = g2;
2203:   prob->gBd[(f*prob->Nf + g)*4+3] = g3;
2204:   return(0);
2205: }

2207: /*@C
2208:   PetscDSGetExactSolution - Get the pointwise exact solution function for a given test field

2210:   Not collective

2212:   Input Parameters:
2213: + prob - The PetscDS
2214: - f    - The test field number

2216:   Output Parameter:
2217: + exactSol - exact solution for the test field
2218: - exactCtx - exact solution context

2220:   Note: The calling sequence for the solution functions is given by:

2222: $ sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)

2224: + dim - the spatial dimension
2225: . t - current time
2226: . x - coordinates of the current point
2227: . Nc - the number of field components
2228: . u - the solution field evaluated at the current point
2229: - ctx - a user context

2231:   Level: intermediate

2233: .seealso: PetscDSSetExactSolution()
2234: @*/
2235: PetscErrorCode PetscDSGetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (**sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void **ctx)
2236: {
2239:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
2242:   return(0);
2243: }

2245: /*@C
2246:   PetscDSSetExactSolution - Set the pointwise exact solution function for a given test field

2248:   Not collective

2250:   Input Parameters:
2251: + prob - The PetscDS
2252: . f    - The test field number
2253: . sol  - solution function for the test fields
2254: - ctx  - solution context or NULL

2256:   Note: The calling sequence for solution functions is given by:

2258: $ sol(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx)

2260: + dim - the spatial dimension
2261: . t - current time
2262: . x - coordinates of the current point
2263: . Nc - the number of field components
2264: . u - the solution field evaluated at the current point
2265: - ctx - a user context

2267:   Level: intermediate

2269: .seealso: PetscDSGetExactSolution()
2270: @*/
2271: PetscErrorCode PetscDSSetExactSolution(PetscDS prob, PetscInt f, PetscErrorCode (*sol)(PetscInt dim, PetscReal t, const PetscReal x[], PetscInt Nc, PetscScalar u[], void *ctx), void *ctx)
2272: {

2277:   if (f < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be non-negative", f);
2278:   PetscDSEnlarge_Static(prob, f+1);
2281:   return(0);
2282: }

2284: /*@C
2285:   PetscDSGetConstants - Returns the array of constants passed to point functions

2287:   Not collective

2289:   Input Parameter:
2290: . prob - The PetscDS object

2292:   Output Parameters:
2293: + numConstants - The number of constants
2294: - constants    - The array of constants, NULL if there are none

2296:   Level: intermediate

2298: .seealso: PetscDSSetConstants(), PetscDSCreate()
2299: @*/
2300: PetscErrorCode PetscDSGetConstants(PetscDS prob, PetscInt *numConstants, const PetscScalar *constants[])
2301: {
2306:   return(0);
2307: }

2309: /*@C
2310:   PetscDSSetConstants - Set the array of constants passed to point functions

2312:   Not collective

2314:   Input Parameters:
2315: + prob         - The PetscDS object
2316: . numConstants - The number of constants
2317: - constants    - The array of constants, NULL if there are none

2319:   Level: intermediate

2321: .seealso: PetscDSGetConstants(), PetscDSCreate()
2322: @*/
2323: PetscErrorCode PetscDSSetConstants(PetscDS prob, PetscInt numConstants, PetscScalar constants[])
2324: {

2329:   if (numConstants != prob->numConstants) {
2330:     PetscFree(prob->constants);
2331:     prob->numConstants = numConstants;
2332:     if (prob->numConstants) {
2333:       PetscMalloc1(prob->numConstants, &prob->constants);
2334:     } else {
2335:       prob->constants = NULL;
2336:     }
2337:   }
2338:   if (prob->numConstants) {
2340:     PetscArraycpy(prob->constants, constants, prob->numConstants);
2341:   }
2342:   return(0);
2343: }

2345: /*@
2346:   PetscDSGetFieldIndex - Returns the index of the given field

2348:   Not collective

2350:   Input Parameters:
2351: + prob - The PetscDS object
2352: - disc - The discretization object

2354:   Output Parameter:
2355: . f - The field number

2357:   Level: beginner

2359: .seealso: PetscGetDiscretization(), PetscDSGetNumFields(), PetscDSCreate()
2360: @*/
2361: PetscErrorCode PetscDSGetFieldIndex(PetscDS prob, PetscObject disc, PetscInt *f)
2362: {
2363:   PetscInt g;

2368:   *f = -1;
2369:   for (g = 0; g < prob->Nf; ++g) {if (disc == prob->disc[g]) break;}
2370:   if (g == prob->Nf) SETERRQ(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Field not found in PetscDS.");
2371:   *f = g;
2372:   return(0);
2373: }

2375: /*@
2376:   PetscDSGetFieldSize - Returns the size of the given field in the full space basis

2378:   Not collective

2380:   Input Parameters:
2381: + prob - The PetscDS object
2382: - f - The field number

2384:   Output Parameter:
2385: . size - The size

2387:   Level: beginner

2389: .seealso: PetscDSGetFieldOffset(), PetscDSGetNumFields(), PetscDSCreate()
2390: @*/
2391: PetscErrorCode PetscDSGetFieldSize(PetscDS prob, PetscInt f, PetscInt *size)
2392: {

2398:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
2399:   PetscDSSetUp(prob);
2400:   *size = prob->Nb[f];
2401:   return(0);
2402: }

2404: /*@
2405:   PetscDSGetFieldOffset - Returns the offset of the given field in the full space basis

2407:   Not collective

2409:   Input Parameters:
2410: + prob - The PetscDS object
2411: - f - The field number

2413:   Output Parameter:
2414: . off - The offset

2416:   Level: beginner

2418: .seealso: PetscDSGetFieldSize(), PetscDSGetNumFields(), PetscDSCreate()
2419: @*/
2420: PetscErrorCode PetscDSGetFieldOffset(PetscDS prob, PetscInt f, PetscInt *off)
2421: {
2422:   PetscInt       size, g;

2428:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
2429:   *off = 0;
2430:   for (g = 0; g < f; ++g) {
2431:     PetscDSGetFieldSize(prob, g, &size);
2432:     *off += size;
2433:   }
2434:   return(0);
2435: }

2437: /*@
2438:   PetscDSGetDimensions - Returns the size of the approximation space for each field on an evaluation point

2440:   Not collective

2442:   Input Parameter:
2443: . prob - The PetscDS object

2445:   Output Parameter:
2446: . dimensions - The number of dimensions

2448:   Level: beginner

2450: .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate()
2451: @*/
2452: PetscErrorCode PetscDSGetDimensions(PetscDS prob, PetscInt *dimensions[])
2453: {

2458:   PetscDSSetUp(prob);
2460:   *dimensions = prob->Nb;
2461:   return(0);
2462: }

2464: /*@
2465:   PetscDSGetComponents - Returns the number of components for each field on an evaluation point

2467:   Not collective

2469:   Input Parameter:
2470: . prob - The PetscDS object

2472:   Output Parameter:
2473: . components - The number of components

2475:   Level: beginner

2477: .seealso: PetscDSGetComponentOffsets(), PetscDSGetNumFields(), PetscDSCreate()
2478: @*/
2479: PetscErrorCode PetscDSGetComponents(PetscDS prob, PetscInt *components[])
2480: {

2485:   PetscDSSetUp(prob);
2487:   *components = prob->Nc;
2488:   return(0);
2489: }

2491: /*@
2492:   PetscDSGetComponentOffset - Returns the offset of the given field on an evaluation point

2494:   Not collective

2496:   Input Parameters:
2497: + prob - The PetscDS object
2498: - f - The field number

2500:   Output Parameter:
2501: . off - The offset

2503:   Level: beginner

2505: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2506: @*/
2507: PetscErrorCode PetscDSGetComponentOffset(PetscDS prob, PetscInt f, PetscInt *off)
2508: {
2512:   if ((f < 0) || (f >= prob->Nf)) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Field number %d must be in [0, %d)", f, prob->Nf);
2513:   *off = prob->off[f];
2514:   return(0);
2515: }

2517: /*@
2518:   PetscDSGetComponentOffsets - Returns the offset of each field on an evaluation point

2520:   Not collective

2522:   Input Parameter:
2523: . prob - The PetscDS object

2525:   Output Parameter:
2526: . offsets - The offsets

2528:   Level: beginner

2530: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2531: @*/
2532: PetscErrorCode PetscDSGetComponentOffsets(PetscDS prob, PetscInt *offsets[])
2533: {
2537:   *offsets = prob->off;
2538:   return(0);
2539: }

2541: /*@
2542:   PetscDSGetComponentDerivativeOffsets - Returns the offset of each field derivative on an evaluation point

2544:   Not collective

2546:   Input Parameter:
2547: . prob - The PetscDS object

2549:   Output Parameter:
2550: . offsets - The offsets

2552:   Level: beginner

2554: .seealso: PetscDSGetNumFields(), PetscDSCreate()
2555: @*/
2556: PetscErrorCode PetscDSGetComponentDerivativeOffsets(PetscDS prob, PetscInt *offsets[])
2557: {
2561:   *offsets = prob->offDer;
2562:   return(0);
2563: }

2565: /*@C
2566:   PetscDSGetTabulation - Return the basis tabulation at quadrature points for the volume discretization

2568:   Not collective

2570:   Input Parameter:
2571: . prob - The PetscDS object

2573:   Output Parameters:
2574: + basis - The basis function tabulation at quadrature points
2575: - basisDer - The basis function derivative tabulation at quadrature points

2577:   Level: intermediate

2579: .seealso: PetscDSCreate()
2580: @*/
2581: PetscErrorCode PetscDSGetTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer)
2582: {

2587:   PetscDSSetUp(prob);
2590:   return(0);
2591: }

2593: /*@C
2594:   PetscDSGetFaceTabulation - Return the basis tabulation at quadrature points on the faces

2596:   Not collective

2598:   Input Parameter:
2599: . prob - The PetscDS object

2601:   Output Parameters:
2602: + basisFace - The basis function tabulation at quadrature points
2603: - basisDerFace - The basis function derivative tabulation at quadrature points

2605:   Level: intermediate

2607: .seealso: PetscDSGetTabulation(), PetscDSCreate()
2608: @*/
2609: PetscErrorCode PetscDSGetFaceTabulation(PetscDS prob, PetscReal ***basis, PetscReal ***basisDer)
2610: {

2615:   PetscDSSetUp(prob);
2618:   return(0);
2619: }

2621: PetscErrorCode PetscDSGetEvaluationArrays(PetscDS prob, PetscScalar **u, PetscScalar **u_t, PetscScalar **u_x)
2622: {

2627:   PetscDSSetUp(prob);
2631:   return(0);
2632: }

2634: PetscErrorCode PetscDSGetWeakFormArrays(PetscDS prob, PetscScalar **f0, PetscScalar **f1, PetscScalar **g0, PetscScalar **g1, PetscScalar **g2, PetscScalar **g3)
2635: {

2640:   PetscDSSetUp(prob);
2647:   return(0);
2648: }

2650: PetscErrorCode PetscDSGetWorkspace(PetscDS prob, PetscReal **x, PetscScalar **basisReal, PetscScalar **basisDerReal, PetscScalar **testReal, PetscScalar **testDerReal)
2651: {

2656:   PetscDSSetUp(prob);
2662:   return(0);
2663: }

2665: /*@C
2666:   PetscDSAddBoundary - Add a boundary condition to the model

2668:   Input Parameters:
2669: + ds          - The PetscDS object
2670: . type        - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2671: . name        - The BC name
2672: . labelname   - The label defining constrained points
2673: . field       - The field to constrain
2674: . numcomps    - The number of constrained field components (0 will constrain all fields)
2675: . comps       - An array of constrained component numbers
2676: . bcFunc      - A pointwise function giving boundary values
2677: . numids      - The number of DMLabel ids for constrained points
2678: . ids         - An array of ids for constrained points
2679: - ctx         - An optional user context for bcFunc

2681:   Options Database Keys:
2682: + -bc_<boundary name> <num> - Overrides the boundary ids
2683: - -bc_<boundary name>_comp <num> - Overrides the boundary components

2685:   Level: developer

2687: .seealso: PetscDSGetBoundary()
2688: @*/
2689: PetscErrorCode PetscDSAddBoundary(PetscDS ds, DMBoundaryConditionType type, const char name[], const char labelname[], PetscInt field, PetscInt numcomps, const PetscInt *comps, void (*bcFunc)(void), PetscInt numids, const PetscInt *ids, void *ctx)
2690: {
2691:   DSBoundary     b;

2696:   PetscNew(&b);
2697:   PetscStrallocpy(name, (char **) &b->name);
2698:   PetscStrallocpy(labelname, (char **) &b->labelname);
2699:   PetscMalloc1(numcomps, &b->comps);
2700:   if (numcomps) {PetscArraycpy(b->comps, comps, numcomps);}
2701:   PetscMalloc1(numids, &b->ids);
2702:   if (numids) {PetscArraycpy(b->ids, ids, numids);}
2703:   b->type            = type;
2704:   b->field           = field;
2705:   b->numcomps        = numcomps;
2706:   b->func            = bcFunc;
2707:   b->numids          = numids;
2708:   b->ctx             = ctx;
2709:   b->next            = ds->boundary;
2710:   ds->boundary       = b;
2711:   return(0);
2712: }

2714: /*@C
2715:   PetscDSUpdateBoundary - Change a boundary condition for the model

2717:   Input Parameters:
2718: + ds          - The PetscDS object
2719: . bd          - The boundary condition number
2720: . type        - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2721: . name        - The BC name
2722: . labelname   - The label defining constrained points
2723: . field       - The field to constrain
2724: . numcomps    - The number of constrained field components
2725: . comps       - An array of constrained component numbers
2726: . bcFunc      - A pointwise function giving boundary values
2727: . numids      - The number of DMLabel ids for constrained points
2728: . ids         - An array of ids for constrained points
2729: - ctx         - An optional user context for bcFunc

2731:   Note: The boundary condition number is the order in which it was registered. The user can get the number of boundary conditions from PetscDSGetNumBoundary().

2733:   Level: developer

2735: .seealso: PetscDSAddBoundary(), PetscDSGetBoundary(), PetscDSGetNumBoundary()
2736: @*/
2737: PetscErrorCode PetscDSUpdateBoundary(PetscDS ds, PetscInt bd, DMBoundaryConditionType type, const char name[], const char labelname[], PetscInt field, PetscInt numcomps, const PetscInt *comps, void (*bcFunc)(void), PetscInt numids, const PetscInt *ids, void *ctx)
2738: {
2739:   DSBoundary     b = ds->boundary;
2740:   PetscInt       n = 0;

2745:   while (b) {
2746:     if (n == bd) break;
2747:     b = b->next;
2748:     ++n;
2749:   }
2750:   if (!b) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %d is not in [0, %d)", bd, n);
2751:   if (name) {
2752:     PetscFree(b->name);
2753:     PetscStrallocpy(name, (char **) &b->name);
2754:   }
2755:   if (labelname) {
2756:     PetscFree(b->labelname);
2757:     PetscStrallocpy(labelname, (char **) &b->labelname);
2758:   }
2759:   if (numcomps >= 0 && numcomps != b->numcomps) {
2760:     b->numcomps = numcomps;
2761:     PetscFree(b->comps);
2762:     PetscMalloc1(numcomps, &b->comps);
2763:     if (numcomps) {PetscArraycpy(b->comps, comps, numcomps);}
2764:   }
2765:   if (numids >= 0 && numids != b->numids) {
2766:     b->numids = numids;
2767:     PetscFree(b->ids);
2768:     PetscMalloc1(numids, &b->ids);
2769:     if (numids) {PetscArraycpy(b->ids, ids, numids);}
2770:   }
2771:   b->type = type;
2772:   if (field >= 0) {b->field  = field;}
2773:   if (bcFunc)     {b->func   = bcFunc;}
2774:   if (ctx)        {b->ctx    = ctx;}
2775:   return(0);
2776: }

2778: /*@
2779:   PetscDSGetNumBoundary - Get the number of registered BC

2781:   Input Parameters:
2782: . ds - The PetscDS object

2784:   Output Parameters:
2785: . numBd - The number of BC

2787:   Level: intermediate

2789: .seealso: PetscDSAddBoundary(), PetscDSGetBoundary()
2790: @*/
2791: PetscErrorCode PetscDSGetNumBoundary(PetscDS ds, PetscInt *numBd)
2792: {
2793:   DSBoundary b = ds->boundary;

2798:   *numBd = 0;
2799:   while (b) {++(*numBd); b = b->next;}
2800:   return(0);
2801: }

2803: /*@C
2804:   PetscDSGetBoundary - Gets a boundary condition to the model

2806:   Input Parameters:
2807: + ds          - The PetscDS object
2808: - bd          - The BC number

2810:   Output Parameters:
2811: + type        - The type of condition, e.g. DM_BC_ESSENTIAL/DM_BC_ESSENTIAL_FIELD (Dirichlet), or DM_BC_NATURAL (Neumann)
2812: . name        - The BC name
2813: . labelname   - The label defining constrained points
2814: . field       - The field to constrain
2815: . numcomps    - The number of constrained field components
2816: . comps       - An array of constrained component numbers
2817: . bcFunc      - A pointwise function giving boundary values
2818: . numids      - The number of DMLabel ids for constrained points
2819: . ids         - An array of ids for constrained points
2820: - ctx         - An optional user context for bcFunc

2822:   Options Database Keys:
2823: + -bc_<boundary name> <num> - Overrides the boundary ids
2824: - -bc_<boundary name>_comp <num> - Overrides the boundary components

2826:   Level: developer

2828: .seealso: PetscDSAddBoundary()
2829: @*/
2830: PetscErrorCode PetscDSGetBoundary(PetscDS ds, PetscInt bd, DMBoundaryConditionType *type, const char **name, const char **labelname, PetscInt *field, PetscInt *numcomps, const PetscInt **comps, void (**func)(void), PetscInt *numids, const PetscInt **ids, void **ctx)
2831: {
2832:   DSBoundary b    = ds->boundary;
2833:   PetscInt   n    = 0;

2837:   while (b) {
2838:     if (n == bd) break;
2839:     b = b->next;
2840:     ++n;
2841:   }
2842:   if (!b) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Boundary %d is not in [0, %d)", bd, n);
2843:   if (type) {
2845:     *type = b->type;
2846:   }
2847:   if (name) {
2849:     *name = b->name;
2850:   }
2851:   if (labelname) {
2853:     *labelname = b->labelname;
2854:   }
2855:   if (field) {
2857:     *field = b->field;
2858:   }
2859:   if (numcomps) {
2861:     *numcomps = b->numcomps;
2862:   }
2863:   if (comps) {
2865:     *comps = b->comps;
2866:   }
2867:   if (func) {
2869:     *func = b->func;
2870:   }
2871:   if (numids) {
2873:     *numids = b->numids;
2874:   }
2875:   if (ids) {
2877:     *ids = b->ids;
2878:   }
2879:   if (ctx) {
2881:     *ctx = b->ctx;
2882:   }
2883:   return(0);
2884: }

2886: /*@
2887:   PetscDSCopyBoundary - Copy all boundary condition objects to the new problem

2889:   Not collective

2891:   Input Parameter:
2892: . prob - The PetscDS object

2894:   Output Parameter:
2895: . newprob - The PetscDS copy

2897:   Level: intermediate

2899: .seealso: PetscDSCopyEquations(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
2900: @*/
2901: PetscErrorCode PetscDSCopyBoundary(PetscDS probA, PetscDS probB)
2902: {
2903:   DSBoundary     b, next, *lastnext;

2909:   if (probA == probB) return(0);
2910:   next = probB->boundary;
2911:   while (next) {
2912:     DSBoundary b = next;

2914:     next = b->next;
2915:     PetscFree(b->comps);
2916:     PetscFree(b->ids);
2917:     PetscFree(b->name);
2918:     PetscFree(b->labelname);
2919:     PetscFree(b);
2920:   }
2921:   lastnext = &(probB->boundary);
2922:   for (b = probA->boundary; b; b = b->next) {
2923:     DSBoundary bNew;

2925:     PetscNew(&bNew);
2926:     bNew->numcomps = b->numcomps;
2927:     PetscMalloc1(bNew->numcomps, &bNew->comps);
2928:     PetscArraycpy(bNew->comps, b->comps, bNew->numcomps);
2929:     bNew->numids = b->numids;
2930:     PetscMalloc1(bNew->numids, &bNew->ids);
2931:     PetscArraycpy(bNew->ids, b->ids, bNew->numids);
2932:     PetscStrallocpy(b->labelname,(char **) &(bNew->labelname));
2933:     PetscStrallocpy(b->name,(char **) &(bNew->name));
2934:     bNew->ctx   = b->ctx;
2935:     bNew->type  = b->type;
2936:     bNew->field = b->field;
2937:     bNew->func  = b->func;

2939:     *lastnext = bNew;
2940:     lastnext = &(bNew->next);
2941:   }
2942:   return(0);
2943: }

2945: /*@C
2946:   PetscDSSelectEquations - Copy pointwise function pointers to the new problem with different field layout

2948:   Not collective

2950:   Input Parameter:
2951: + prob - The PetscDS object
2952: . numFields - Number of new fields
2953: - fields - Old field number for each new field

2955:   Output Parameter:
2956: . newprob - The PetscDS copy

2958:   Level: intermediate

2960: .seealso: PetscDSCopyBoundary(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
2961: @*/
2962: PetscErrorCode PetscDSSelectEquations(PetscDS prob, PetscInt numFields, const PetscInt fields[], PetscDS newprob)
2963: {
2964:   PetscInt       Nf, Nfn, fn, gn;

2971:   PetscDSGetNumFields(prob, &Nf);
2972:   PetscDSGetNumFields(newprob, &Nfn);
2973:   if (numFields > Nfn) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Number of fields %D to transfer must not be greater then the total number of fields %D", numFields, Nfn);
2974:   for (fn = 0; fn < numFields; ++fn) {
2975:     const PetscInt   f = fields ? fields[fn] : fn;
2976:     PetscPointFunc   obj;
2977:     PetscPointFunc   f0, f1;
2978:     PetscBdPointFunc f0Bd, f1Bd;
2979:     PetscRiemannFunc r;

2981:     if (f >= Nf) continue;
2982:     PetscDSGetObjective(prob, f, &obj);
2983:     PetscDSGetResidual(prob, f, &f0, &f1);
2984:     PetscDSGetBdResidual(prob, f, &f0Bd, &f1Bd);
2985:     PetscDSGetRiemannSolver(prob, f, &r);
2986:     PetscDSSetObjective(newprob, fn, obj);
2987:     PetscDSSetResidual(newprob, fn, f0, f1);
2988:     PetscDSSetBdResidual(newprob, fn, f0Bd, f1Bd);
2989:     PetscDSSetRiemannSolver(newprob, fn, r);
2990:     for (gn = 0; gn < numFields; ++gn) {
2991:       const PetscInt  g = fields ? fields[gn] : gn;
2992:       PetscPointJac   g0, g1, g2, g3;
2993:       PetscPointJac   g0p, g1p, g2p, g3p;
2994:       PetscBdPointJac g0Bd, g1Bd, g2Bd, g3Bd;

2996:       if (g >= Nf) continue;
2997:       PetscDSGetJacobian(prob, f, g, &g0, &g1, &g2, &g3);
2998:       PetscDSGetJacobianPreconditioner(prob, f, g, &g0p, &g1p, &g2p, &g3p);
2999:       PetscDSGetBdJacobian(prob, f, g, &g0Bd, &g1Bd, &g2Bd, &g3Bd);
3000:       PetscDSSetJacobian(newprob, fn, gn, g0, g1, g2, g3);
3001:       PetscDSSetJacobianPreconditioner(prob, fn, gn, g0p, g1p, g2p, g3p);
3002:       PetscDSSetBdJacobian(newprob, fn, gn, g0Bd, g1Bd, g2Bd, g3Bd);
3003:     }
3004:   }
3005:   return(0);
3006: }

3008: /*@
3009:   PetscDSCopyEquations - Copy all pointwise function pointers to the new problem

3011:   Not collective

3013:   Input Parameter:
3014: . prob - The PetscDS object

3016:   Output Parameter:
3017: . newprob - The PetscDS copy

3019:   Level: intermediate

3021: .seealso: PetscDSCopyBoundary(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
3022: @*/
3023: PetscErrorCode PetscDSCopyEquations(PetscDS prob, PetscDS newprob)
3024: {
3025:   PetscInt       Nf, Ng;

3031:   PetscDSGetNumFields(prob, &Nf);
3032:   PetscDSGetNumFields(newprob, &Ng);
3033:   if (Nf != Ng) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_SIZ, "Number of fields must match %D != %D", Nf, Ng);
3034:   PetscDSSelectEquations(prob, Nf, NULL, newprob);
3035:   return(0);
3036: }
3037: /*@
3038:   PetscDSCopyConstants - Copy all constants to the new problem

3040:   Not collective

3042:   Input Parameter:
3043: . prob - The PetscDS object

3045:   Output Parameter:
3046: . newprob - The PetscDS copy

3048:   Level: intermediate

3050: .seealso: PetscDSCopyBoundary(), PetscDSCopyEquations(), PetscDSSetResidual(), PetscDSSetJacobian(), PetscDSSetRiemannSolver(), PetscDSSetBdResidual(), PetscDSSetBdJacobian(), PetscDSCreate()
3051: @*/
3052: PetscErrorCode PetscDSCopyConstants(PetscDS prob, PetscDS newprob)
3053: {
3054:   PetscInt           Nc;
3055:   const PetscScalar *constants;
3056:   PetscErrorCode     ierr;

3061:   PetscDSGetConstants(prob, &Nc, &constants);
3062:   PetscDSSetConstants(newprob, Nc, (PetscScalar *) constants);
3063:   return(0);
3064: }

3066: PetscErrorCode PetscDSGetHeightSubspace(PetscDS prob, PetscInt height, PetscDS *subprob)
3067: {
3068:   PetscInt       dim, Nf, f;

3074:   if (height == 0) {*subprob = prob; return(0);}
3075:   PetscDSGetNumFields(prob, &Nf);
3076:   PetscDSGetSpatialDimension(prob, &dim);
3077:   if (height > dim) SETERRQ2(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_OUTOFRANGE, "DS can only handle height in [0, %D], not %D", dim, height);
3078:   if (!prob->subprobs) {PetscCalloc1(dim, &prob->subprobs);}
3079:   if (!prob->subprobs[height-1]) {
3080:     PetscInt cdim;

3082:     PetscDSCreate(PetscObjectComm((PetscObject) prob), &prob->subprobs[height-1]);
3083:     PetscDSGetCoordinateDimension(prob, &cdim);
3084:     PetscDSSetCoordinateDimension(prob->subprobs[height-1], cdim);
3085:     for (f = 0; f < Nf; ++f) {
3086:       PetscFE      subfe;
3087:       PetscObject  obj;
3088:       PetscClassId id;

3090:       PetscDSGetDiscretization(prob, f, &obj);
3091:       PetscObjectGetClassId(obj, &id);
3092:       if (id == PETSCFE_CLASSID) {PetscFEGetHeightSubspace((PetscFE) obj, height, &subfe);}
3093:       else SETERRQ1(PetscObjectComm((PetscObject) prob), PETSC_ERR_ARG_WRONG, "Unsupported discretization type for field %d", f);
3094:       PetscDSSetDiscretization(prob->subprobs[height-1], f, (PetscObject) subfe);
3095:     }
3096:   }
3097:   *subprob = prob->subprobs[height-1];
3098:   return(0);
3099: }

3101: static PetscErrorCode PetscDSDestroy_Basic(PetscDS prob)
3102: {
3103:   PetscErrorCode      ierr;

3106:   PetscFree(prob->data);
3107:   return(0);
3108: }

3110: static PetscErrorCode PetscDSInitialize_Basic(PetscDS prob)
3111: {
3113:   prob->ops->setfromoptions = NULL;
3114:   prob->ops->setup          = NULL;
3115:   prob->ops->view           = NULL;
3116:   prob->ops->destroy        = PetscDSDestroy_Basic;
3117:   return(0);
3118: }

3120: /*MC
3121:   PETSCDSBASIC = "basic" - A discrete system with pointwise residual and boundary residual functions

3123:   Level: intermediate

3125: .seealso: PetscDSType, PetscDSCreate(), PetscDSSetType()
3126: M*/

3128: PETSC_EXTERN PetscErrorCode PetscDSCreate_Basic(PetscDS prob)
3129: {
3130:   PetscDS_Basic *b;
3131:   PetscErrorCode      ierr;

3135:   PetscNewLog(prob, &b);
3136:   prob->data = b;

3138:   PetscDSInitialize_Basic(prob);
3139:   return(0);
3140: }