Actual source code: plexproject.c
petsc-3.9.4 2018-09-11
1: #include <petsc/private/dmpleximpl.h>
3: #include <petsc/private/petscfeimpl.h>
5: static PetscErrorCode DMProjectPoint_Func_Private(DM dm, PetscDS prob, PetscReal time, PetscFECellGeom *fegeom, PetscFVCellGeom *fvgeom, PetscBool isFE[], PetscDualSpace sp[],
6: PetscErrorCode (**funcs)(PetscInt, PetscReal, const PetscReal [], PetscInt, PetscScalar *, void *), void **ctxs,
7: PetscScalar values[])
8: {
9: PetscInt Nf, *Nc, f, totDim, spDim, d, v;
13: PetscDSGetNumFields(prob, &Nf);
14: PetscDSGetComponents(prob, &Nc);
15: PetscDSGetTotalDimension(prob, &totDim);
16: /* Get values for closure */
17: for (f = 0, v = 0; f < Nf; ++f) {
18: void * const ctx = ctxs ? ctxs[f] : NULL;
20: if (!sp[f]) continue;
21: PetscDualSpaceGetDimension(sp[f], &spDim);
22: for (d = 0; d < spDim; ++d, ++v) {
23: if (funcs[f]) {
24: if (isFE[f]) {PetscDualSpaceApply(sp[f], d, time, fegeom, Nc[f], funcs[f], ctx, &values[v]);}
25: else {PetscDualSpaceApplyFVM(sp[f], d, time, fvgeom, Nc[f], funcs[f], ctx, &values[v]);}
26: } else {
27: values[v] = 0.0;
28: }
29: }
30: }
31: return(0);
32: }
34: static PetscErrorCode DMProjectPoint_Field_Private(DM dm, PetscDS prob, DM dmAux, PetscDS probAux, PetscReal time, Vec localU, Vec localA, PetscFECellGeom *fegeom, PetscDualSpace sp[], PetscInt p, PetscInt Ncc, const PetscInt comps[],
35: PetscReal **basisTab, PetscReal **basisDerTab, PetscReal **basisTabAux, PetscReal **basisDerTabAux,
36: void (**funcs)(PetscInt, PetscInt, PetscInt,
37: const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
38: const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
39: PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]), void **ctxs,
40: PetscScalar values[])
41: {
42: PetscSection section, sectionAux = NULL;
43: PetscScalar *u, *u_t = NULL, *u_x, *a = NULL, *a_t = NULL, *a_x = NULL, *refSpaceDer, *refSpaceDerAux = NULL, *bc;
44: PetscScalar *coefficients = NULL, *coefficientsAux = NULL;
45: PetscScalar *coefficients_t = NULL, *coefficientsAux_t = NULL;
46: const PetscScalar *constants;
47: PetscReal *x;
48: PetscInt *uOff, *uOff_x, *aOff = NULL, *aOff_x = NULL, *Nb, *Nc, *NbAux = NULL, *NcAux = NULL;
49: const PetscInt dim = fegeom->dim, dimEmbed = fegeom->dimEmbed;
50: PetscInt dimAux = 0, numConstants, Nf, NfAux = 0, f, spDim, d, c, v, tp = 0;
51: PetscErrorCode ierr;
54: PetscDSGetNumFields(prob, &Nf);
55: PetscDSGetDimensions(prob, &Nb);
56: PetscDSGetComponents(prob, &Nc);
57: PetscDSGetComponentOffsets(prob, &uOff);
58: PetscDSGetComponentDerivativeOffsets(prob, &uOff_x);
59: PetscDSGetEvaluationArrays(prob, &u, &bc /*&u_t*/, &u_x);
60: PetscDSGetRefCoordArrays(prob, &x, &refSpaceDer);
61: PetscDSGetConstants(prob, &numConstants, &constants);
62: DMGetDefaultSection(dm, §ion);
63: DMPlexVecGetClosure(dm, section, localU, p, NULL, &coefficients);
64: if (dmAux) {
65: DMLabel spmap;
66: PetscInt subp;
68: DMPlexGetSubpointMap(dmAux, &spmap);
69: if (spmap) {
70: IS subpointIS;
71: const PetscInt *subpoints;
72: PetscInt numSubpoints;
74: DMPlexCreateSubpointIS(dmAux, &subpointIS);
75: ISGetLocalSize(subpointIS, &numSubpoints);
76: ISGetIndices(subpointIS, &subpoints);
77: PetscFindInt(p, numSubpoints, subpoints, &subp);
78: if (subp < 0) SETERRQ1(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Point %d not found in submesh", p);
79: ISRestoreIndices(subpointIS, &subpoints);
80: ISDestroy(&subpointIS);
81: } else subp = p;
82: PetscDSGetSpatialDimension(probAux, &dimAux);
83: PetscDSGetNumFields(probAux, &NfAux);
84: PetscDSGetDimensions(probAux, &NbAux);
85: PetscDSGetComponents(probAux, &NcAux);
86: DMGetDefaultSection(dmAux, §ionAux);
87: PetscDSGetComponentOffsets(probAux, &aOff);
88: PetscDSGetComponentDerivativeOffsets(probAux, &aOff_x);
89: PetscDSGetEvaluationArrays(probAux, &a, NULL /*&a_t*/, &a_x);
90: PetscDSGetRefCoordArrays(probAux, NULL, &refSpaceDerAux);
91: DMPlexVecGetClosure(dmAux, sectionAux, localA, subp, NULL, &coefficientsAux);
92: }
93: /* Get values for closure */
94: for (f = 0, v = 0; f < Nf; ++f) {
95: if (!sp[f]) continue;
96: PetscDualSpaceGetDimension(sp[f], &spDim);
97: for (d = 0; d < spDim; ++d, ++v) {
98: if (funcs[f]) {
99: PetscQuadrature quad;
100: const PetscReal *points, *weights;
101: PetscInt numPoints, q;
103: PetscDualSpaceGetFunctional(sp[f], d, &quad);
104: PetscQuadratureGetData(quad, NULL, NULL, &numPoints, &points, &weights);
105: for (q = 0; q < numPoints; ++q, ++tp) {
106: CoordinatesRefToReal(dimEmbed, dim, fegeom->v0, fegeom->J, &points[q*dim], x);
107: EvaluateFieldJets(dim, Nf, Nb, Nc, tp, basisTab, basisDerTab, refSpaceDer, fegeom->invJ, coefficients, coefficients_t, u, u_x, u_t);
108: if (probAux) {EvaluateFieldJets(dimAux, NfAux, NbAux, NcAux, tp, basisTabAux, basisDerTabAux, refSpaceDerAux, fegeom->invJ, coefficientsAux, coefficientsAux_t, a, a_x, a_t);}
109: (*funcs[f])(dim, Nf, NfAux, uOff, uOff_x, u, u_t, u_x, aOff, aOff_x, a, a_t, a_x, time, x, numConstants, constants, bc);
110: if (Ncc) {
111: for (c = 0; c < Ncc; ++c) values[v] += bc[comps[c]]*weights[comps[c]];
112: } else {
113: for (c = 0; c < Nc[f]; ++c) values[v] += bc[c]*weights[c];
114: }
115: }
116: } else {
117: values[v] = 0.0;
118: }
119: }
120: }
121: DMPlexVecRestoreClosure(dm, section, localU, p, NULL, &coefficients);
122: if (dmAux) {DMPlexVecRestoreClosure(dmAux, sectionAux, localA, p, NULL, &coefficientsAux);}
123: return(0);
124: }
126: static PetscErrorCode DMProjectPoint_Private(DM dm, PetscDS prob, DM dmAux, PetscDS probAux, PetscInt effectiveHeight, PetscReal time, Vec localU, Vec localA, PetscBool hasFE, PetscBool hasFV, PetscBool isFE[],
127: PetscDualSpace sp[], PetscInt p, PetscInt Ncc, const PetscInt comps[],
128: PetscReal **basisTab, PetscReal **basisDerTab, PetscReal **basisTabAux, PetscReal **basisDerTabAux,
129: DMBoundaryConditionType type, void (**funcs)(void), void **ctxs, PetscBool fieldActive[], PetscScalar values[])
130: {
131: PetscFECellGeom fegeom;
132: PetscFVCellGeom fvgeom;
133: PetscInt dim, dimEmbed;
134: PetscErrorCode ierr;
137: DMGetDimension(dm, &dim);
138: DMGetCoordinateDim(dm, &dimEmbed);
139: if (hasFE) {
140: DMPlexComputeCellGeometryFEM(dm, p, NULL, fegeom.v0, fegeom.J, fegeom.invJ, &fegeom.detJ);
141: fegeom.dim = dim - effectiveHeight;
142: fegeom.dimEmbed = dimEmbed;
143: }
144: if (hasFV) {DMPlexComputeCellGeometryFVM(dm, p, &fvgeom.volume, fvgeom.centroid, NULL);}
145: switch (type) {
146: case DM_BC_ESSENTIAL:
147: case DM_BC_NATURAL:
148: DMProjectPoint_Func_Private(dm, prob, time, &fegeom, &fvgeom, isFE, sp, (PetscErrorCode (**)(PetscInt, PetscReal, const PetscReal [], PetscInt, PetscScalar *, void *)) funcs, ctxs, values);break;
149: case DM_BC_ESSENTIAL_FIELD:
150: case DM_BC_NATURAL_FIELD:
151: DMProjectPoint_Field_Private(dm, prob, dmAux, probAux, time, localU, localA, &fegeom, sp, p, Ncc, comps,
152: basisTab, basisDerTab, basisTabAux, basisDerTabAux,
153: (void (**)(PetscInt, PetscInt, PetscInt,
154: const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
155: const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
156: PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[])) funcs, ctxs, values);break;
157: default: SETERRQ1(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "Unknown boundary condition type: %d", (int) type);
158: }
159: return(0);
160: }
162: /*
163: This function iterates over a manifold, and interpolates the input function/field using the basis provided by the DS in our DM
165: There are several different scenarios:
167: 1) Volumetric mesh with volumetric auxiliary data
169: Here minHeight=0 since we loop over cells.
171: 2) Boundary mesh with boundary auxiliary data
173: Here minHeight=1 since we loop over faces. This normally happens since we hang cells off of our boundary meshes to facilitate computation.
175: 3) Volumetric mesh with boundary auxiliary data
177: Here minHeight=1 and auxbd=PETSC_TRUE since we loop over faces and use data only supported on those faces. This is common when imposing Dirichlet boundary conditions.
179: The maxHeight is used to support enforcement of constraints in DMForest.
181: If localU is given and not equal to localX, we call DMPlexInsertBoundaryValues() to complete it.
183: If we are using an input field (DM_BC_ESSENTIAL_FIELD or DM_BC_NATURAL_FIELD), we need to evaluate it at all the quadrature points of the dual basis functionals.
184: - We use effectiveHeight to mean the height above our incoming DS. For example, if the DS is for a submesh then the effective height is zero, whereas if the DS
185: is for the volumetric mesh, but we are iterating over a surface, then the effective height is nonzero. When th effective height is nonzero, we need to extract
186: dual spaces for the boundary from our input spaces.
187: - After extracting all quadrature points, we tabulate the input fields and auxiliary fields on them.
189: We check that the #dof(closure(p)) == #dual basis functionals(p) for a representative p in the iteration
191: If we have a label, we iterate over those points. This will probably break the maxHeight functionality since we do not check the height of those points.
192: */
193: static PetscErrorCode DMProjectLocal_Generic_Plex(DM dm, PetscReal time, Vec localU,
194: PetscInt Ncc, const PetscInt comps[], DMLabel label, PetscInt numIds, const PetscInt ids[],
195: DMBoundaryConditionType type, void (**funcs)(void), void **ctxs,
196: InsertMode mode, Vec localX)
197: {
198: DM dmAux = NULL;
199: PetscDS prob, probAux = NULL;
200: Vec localA = NULL;
201: PetscSection section;
202: PetscDualSpace *sp, *cellsp;
203: PetscReal **basisTab = NULL, **basisDerTab = NULL, **basisTabAux = NULL, **basisDerTabAux = NULL;
204: PetscInt *Nc;
205: PetscInt dim, dimEmbed, minHeight, maxHeight, h, Nf, NfAux = 0, f;
206: PetscBool *isFE, hasFE = PETSC_FALSE, hasFV = PETSC_FALSE, auxBd = PETSC_FALSE;
207: PetscErrorCode ierr;
210: PetscObjectQuery((PetscObject) dm, "dmAux", (PetscObject *) &dmAux);
211: PetscObjectQuery((PetscObject) dm, "A", (PetscObject *) &localA);
212: DMGetDimension(dm, &dim);
213: DMPlexGetVTKCellHeight(dm, &minHeight);
214: /* Auxiliary information can only be used with interpolation of field functions */
215: if (type == DM_BC_ESSENTIAL_FIELD || type == DM_BC_NATURAL_FIELD) {
216: if (!minHeight && dmAux) {
217: DMLabel spmap;
219: /* If dmAux is a surface, then force the projection to take place over a surface */
220: DMPlexGetSubpointMap(dmAux, &spmap);
221: if (spmap) {minHeight = 1; auxBd = PETSC_TRUE;}
222: }
223: }
224: DMPlexGetMaxProjectionHeight(dm, &maxHeight);
225: maxHeight = PetscMax(maxHeight, minHeight);
226: if (maxHeight < 0 || maxHeight > dim) {SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Maximum projection height %d not in [0, %d)\n", maxHeight, dim);}
227: DMGetDS(dm, &prob);
228: DMGetCoordinateDim(dm, &dimEmbed);
229: DMGetDefaultSection(dm, §ion);
230: PetscSectionGetNumFields(section, &Nf);
231: if (dmAux) {
232: DMGetDS(dmAux, &probAux);
233: PetscDSGetNumFields(probAux, &NfAux);
234: }
235: PetscDSGetComponents(prob, &Nc);
236: PetscMalloc2(Nf, &isFE, Nf, &sp);
237: if (maxHeight > 0) {PetscMalloc1(Nf, &cellsp);}
238: else {cellsp = sp;}
239: if (localU && localU != localX) {DMPlexInsertBoundaryValues(dm, PETSC_TRUE, localU, time, NULL, NULL, NULL);}
240: /* Get cell dual spaces */
241: for (f = 0; f < Nf; ++f) {
242: PetscObject obj;
243: PetscClassId id;
245: DMGetField(dm, f, &obj);
246: PetscObjectGetClassId(obj, &id);
247: if (id == PETSCFE_CLASSID) {
248: PetscFE fe = (PetscFE) obj;
250: hasFE = PETSC_TRUE;
251: isFE[f] = PETSC_TRUE;
252: PetscFEGetDualSpace(fe, &cellsp[f]);
253: } else if (id == PETSCFV_CLASSID) {
254: PetscFV fv = (PetscFV) obj;
256: hasFV = PETSC_TRUE;
257: isFE[f] = PETSC_FALSE;
258: PetscFVGetDualSpace(fv, &cellsp[f]);
259: } else SETERRQ1(PetscObjectComm((PetscObject) dm), PETSC_ERR_ARG_WRONG, "Unknown discretization type for field %d", f);
260: }
261: if (type == DM_BC_ESSENTIAL_FIELD || type == DM_BC_NATURAL_FIELD) {
262: PetscInt effectiveHeight = auxBd ? minHeight : 0;
263: PetscFE fem, subfem;
264: PetscReal *points;
265: PetscInt numPoints, spDim, qdim = 0, d;
267: if (maxHeight > minHeight) SETERRQ(PetscObjectComm((PetscObject) dm), PETSC_ERR_SUP, "Field projection not supported for face interpolation");
268: numPoints = 0;
269: for (f = 0; f < Nf; ++f) {
270: if (!effectiveHeight) {sp[f] = cellsp[f];}
271: else {PetscDualSpaceGetHeightSubspace(cellsp[f], effectiveHeight, &sp[f]);}
272: PetscDualSpaceGetDimension(sp[f], &spDim);
273: for (d = 0; d < spDim; ++d) {
274: if (funcs[f]) {
275: PetscQuadrature quad;
276: PetscInt Nq;
278: PetscDualSpaceGetFunctional(sp[f], d, &quad);
279: PetscQuadratureGetData(quad, &qdim, NULL, &Nq, NULL, NULL);
280: numPoints += Nq;
281: }
282: }
283: }
284: PetscMalloc1(numPoints*qdim, &points);
285: numPoints = 0;
286: for (f = 0; f < Nf; ++f) {
287: PetscDualSpaceGetDimension(sp[f], &spDim);
288: for (d = 0; d < spDim; ++d) {
289: if (funcs[f]) {
290: PetscQuadrature quad;
291: const PetscReal *qpoints;
292: PetscInt Nq, q;
294: PetscDualSpaceGetFunctional(sp[f], d, &quad);
295: PetscQuadratureGetData(quad, NULL, NULL, &Nq, &qpoints, NULL);
296: for (q = 0; q < Nq*qdim; ++q) points[numPoints*qdim+q] = qpoints[q];
297: numPoints += Nq;
298: }
299: }
300: }
301: PetscMalloc4(Nf, &basisTab, Nf, &basisDerTab, NfAux, &basisTabAux, NfAux, &basisDerTabAux);
302: for (f = 0; f < Nf; ++f) {
303: if (!isFE[f]) continue;
304: PetscDSGetDiscretization(prob, f, (PetscObject *) &fem);
305: if (!effectiveHeight) {subfem = fem;}
306: else {PetscFEGetHeightSubspace(fem, effectiveHeight, &subfem);}
307: PetscFEGetTabulation(subfem, numPoints, points, &basisTab[f], &basisDerTab[f], NULL);
308: }
309: for (f = 0; f < NfAux; ++f) {
310: PetscDSGetDiscretization(probAux, f, (PetscObject *) &fem);
311: if (!effectiveHeight || auxBd) {subfem = fem;}
312: else {PetscFEGetHeightSubspace(fem, effectiveHeight, &subfem);}
313: PetscFEGetTabulation(subfem, numPoints, points, &basisTabAux[f], &basisDerTabAux[f], NULL);
314: }
315: PetscFree(points);
316: }
317: /* Note: We make no attempt to optimize for height. Higher height things just overwrite the lower height results. */
318: for (h = minHeight; h <= maxHeight; h++) {
319: PetscInt effectiveHeight = h - (auxBd ? 0 : minHeight);
320: PetscDS probEff = prob;
321: PetscScalar *values;
322: PetscBool *fieldActive;
323: PetscInt pStart, pEnd, p, spDim, totDim, numValues;
325: if (effectiveHeight) {PetscDSGetHeightSubspace(prob, effectiveHeight, &probEff);}
326: DMPlexGetHeightStratum(dm, h, &pStart, &pEnd);
327: if (!h) {
328: PetscInt cEndInterior;
330: DMPlexGetHybridBounds(dm, &cEndInterior, NULL, NULL, NULL);
331: pEnd = cEndInterior < 0 ? pEnd : cEndInterior;
332: }
333: if (pEnd <= pStart) continue;
334: /* Compute totDim, the number of dofs in the closure of a point at this height */
335: totDim = 0;
336: for (f = 0; f < Nf; ++f) {
337: if (!effectiveHeight) {
338: sp[f] = cellsp[f];
339: } else {
340: PetscDualSpaceGetHeightSubspace(cellsp[f], effectiveHeight, &sp[f]);
341: if (!sp[f]) continue;
342: }
343: PetscDualSpaceGetDimension(sp[f], &spDim);
344: totDim += spDim;
345: }
346: DMPlexVecGetClosure(dm, section, localX, pStart, &numValues, NULL);
347: if (numValues != totDim) SETERRQ2(PETSC_COMM_SELF, PETSC_ERR_ARG_SIZ, "The section point closure size %d != dual space dimension %d", numValues, totDim);
348: if (!totDim) continue;
349: /* Loop over points at this height */
350: DMGetWorkArray(dm, numValues, MPIU_SCALAR, &values);
351: DMGetWorkArray(dm, Nf, MPI_INT, &fieldActive);
352: for (f = 0; f < Nf; ++f) fieldActive[f] = (funcs[f] && sp[f]) ? PETSC_TRUE : PETSC_FALSE;
353: if (label) {
354: PetscInt i;
356: for (i = 0; i < numIds; ++i) {
357: IS pointIS;
358: const PetscInt *points;
359: PetscInt n;
361: DMLabelGetStratumIS(label, ids[i], &pointIS);
362: if (!pointIS) continue; /* No points with that id on this process */
363: ISGetLocalSize(pointIS, &n);
364: ISGetIndices(pointIS, &points);
365: for (p = 0; p < n; ++p) {
366: const PetscInt point = points[p];
368: if ((point < pStart) || (point >= pEnd)) continue;
369: PetscMemzero(values, numValues * sizeof(PetscScalar));
370: DMProjectPoint_Private(dm, probEff, dmAux, probAux, effectiveHeight, time, localU, localA, hasFE, hasFV, isFE, sp, point, Ncc, comps, basisTab, basisDerTab, basisTabAux, basisDerTabAux, type, funcs, ctxs, fieldActive, values);
371: if (ierr) {
372: PetscErrorCode ierr2;
373: ierr2 = DMRestoreWorkArray(dm, numValues, MPIU_SCALAR, &values);CHKERRQ(ierr2);
374: ierr2 = DMRestoreWorkArray(dm, Nf, MPI_INT, &fieldActive);CHKERRQ(ierr2);
375:
376: }
377: DMPlexVecSetFieldClosure_Internal(dm, section, localX, fieldActive, point, Ncc, comps, values, mode);
378: }
379: ISRestoreIndices(pointIS, &points);
380: ISDestroy(&pointIS);
381: }
382: } else {
383: for (p = pStart; p < pEnd; ++p) {
384: PetscMemzero(values, numValues * sizeof(PetscScalar));
385: DMProjectPoint_Private(dm, probEff, dmAux, probAux, effectiveHeight, time, localU, localA, hasFE, hasFV, isFE, sp, p, Ncc, comps, basisTab, basisDerTab, basisTabAux, basisDerTabAux, type, funcs, ctxs, fieldActive, values);
386: if (ierr) {
387: PetscErrorCode ierr2;
388: ierr2 = DMRestoreWorkArray(dm, numValues, MPIU_SCALAR, &values);CHKERRQ(ierr2);
389: ierr2 = DMRestoreWorkArray(dm, Nf, MPI_INT, &fieldActive);CHKERRQ(ierr2);
390:
391: }
392: DMPlexVecSetFieldClosure_Internal(dm, section, localX, fieldActive, p, Ncc, comps, values, mode);
393: }
394: }
395: DMRestoreWorkArray(dm, numValues, MPIU_SCALAR, &values);
396: DMRestoreWorkArray(dm, Nf, MPI_INT, &fieldActive);
397: }
398: /* Cleanup */
399: if (type == DM_BC_ESSENTIAL_FIELD || type == DM_BC_NATURAL_FIELD) {
400: PetscInt effectiveHeight = auxBd ? minHeight : 0;
401: PetscFE fem, subfem;
403: for (f = 0; f < Nf; ++f) {
404: if (!isFE[f]) continue;
405: PetscDSGetDiscretization(prob, f, (PetscObject *) &fem);
406: if (!effectiveHeight) {subfem = fem;}
407: else {PetscFEGetHeightSubspace(fem, effectiveHeight, &subfem);}
408: PetscFERestoreTabulation(subfem, 0, NULL, &basisTab[f], &basisDerTab[f], NULL);
409: }
410: for (f = 0; f < NfAux; ++f) {
411: PetscDSGetDiscretization(probAux, f, (PetscObject *) &fem);
412: if (!effectiveHeight || auxBd) {subfem = fem;}
413: else {PetscFEGetHeightSubspace(fem, effectiveHeight, &subfem);}
414: PetscFERestoreTabulation(subfem, 0, NULL, &basisTabAux[f], &basisDerTabAux[f], NULL);
415: }
416: PetscFree4(basisTab, basisDerTab, basisTabAux, basisDerTabAux);
417: }
418: PetscFree2(isFE, sp);
419: if (maxHeight > 0) {PetscFree(cellsp);}
420: return(0);
421: }
423: PetscErrorCode DMProjectFunctionLocal_Plex(DM dm, PetscReal time, PetscErrorCode (**funcs)(PetscInt, PetscReal, const PetscReal [], PetscInt, PetscScalar *, void *), void **ctxs, InsertMode mode, Vec localX)
424: {
428: DMProjectLocal_Generic_Plex(dm, time, localX, 0, NULL, NULL, 0, NULL, DM_BC_ESSENTIAL, (void (**)(void)) funcs, ctxs, mode, localX);
429: return(0);
430: }
432: PetscErrorCode DMProjectFunctionLabelLocal_Plex(DM dm, PetscReal time, DMLabel label, PetscInt numIds, const PetscInt ids[], PetscInt Ncc, const PetscInt comps[], PetscErrorCode (**funcs)(PetscInt, PetscReal, const PetscReal [], PetscInt, PetscScalar *, void *), void **ctxs, InsertMode mode, Vec localX)
433: {
437: DMProjectLocal_Generic_Plex(dm, time, localX, Ncc, comps, label, numIds, ids, DM_BC_ESSENTIAL, (void (**)(void)) funcs, ctxs, mode, localX);
438: return(0);
439: }
441: PetscErrorCode DMProjectFieldLocal_Plex(DM dm, PetscReal time, Vec localU,
442: void (**funcs)(PetscInt, PetscInt, PetscInt,
443: const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
444: const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
445: PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
446: InsertMode mode, Vec localX)
447: {
451: DMProjectLocal_Generic_Plex(dm, time, localU, 0, NULL, NULL, 0, NULL, DM_BC_ESSENTIAL_FIELD, (void (**)(void)) funcs, NULL, mode, localX);
452: return(0);
453: }
455: PetscErrorCode DMProjectFieldLabelLocal_Plex(DM dm, PetscReal time, DMLabel label, PetscInt numIds, const PetscInt ids[], PetscInt Ncc, const PetscInt comps[], Vec localU,
456: void (**funcs)(PetscInt, PetscInt, PetscInt,
457: const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
458: const PetscInt[], const PetscInt[], const PetscScalar[], const PetscScalar[], const PetscScalar[],
459: PetscReal, const PetscReal[], PetscInt, const PetscScalar[], PetscScalar[]),
460: InsertMode mode, Vec localX)
461: {
465: DMProjectLocal_Generic_Plex(dm, time, localU, Ncc, comps, label, numIds, ids, DM_BC_ESSENTIAL_FIELD, (void (**)(void)) funcs, NULL, mode, localX);
466: return(0);
467: }