Actual source code: stagstencil.c

  1: /* Functions concerning getting and setting Vec and Mat values with DMStagStencil */
  2: #include <petsc/private/dmstagimpl.h>

  4: /* Strings corresponding to the types defined in $PETSC_DIR/include/petscdmstag.h */
  5: const char *const DMStagStencilTypes[] = {"NONE", "STAR", "BOX", "DMStagStencilType", "DM_STAG_STENCIL_", NULL};

  7: /* Strings corresponding the positions in $PETSC_DIR/include/petscdmstag.h */
  8: const char *const DMStagStencilLocations[] = {"NONE", "BACK_DOWN_LEFT", "BACK_DOWN", "BACK_DOWN_RIGHT", "BACK_LEFT", "BACK", "BACK_RIGHT", "BACK_UP_LEFT", "BACK_UP", "BACK_UP_RIGHT", "DOWN_LEFT", "DOWN", "DOWN_RIGHT", "LEFT", "ELEMENT", "RIGHT", "UP_LEFT", "UP", "UP_RIGHT", "FRONT_DOWN_LEFT", "FRONT_DOWN", "FRONT_DOWN_RIGHT", "FRONT_LEFT", "FRONT", "FRONT_RIGHT", "FRONT_UP_LEFT", "FRONT_UP", "FRONT_UP_RIGHT", "DMStagStencilLocation", "", NULL};

 10: /*@C
 11:   DMStagCreateISFromStencils - Create an `IS`, using global numberings, for a subset of DOF in a `DMSTAG` object

 13:   Collective

 15:   Input Parameters:
 16: + dm        - the `DMStag` object
 17: . n_stencil - the number of stencils provided
 18: - stencils  - an array of `DMStagStencil` objects (`i`, `j`, and `k` are ignored)

 20:   Output Parameter:
 21: . is - the global IS

 23:   Note:
 24:   Redundant entries in the stencils argument are ignored

 26:   Level: advanced

 28: .seealso: [](ch_stag), `DMSTAG`, `IS`, `DMStagStencil`, `DMCreateGlobalVector`
 29: @*/
 30: PetscErrorCode DMStagCreateISFromStencils(DM dm, PetscInt n_stencil, DMStagStencil *stencils, IS *is)
 31: {
 32:   PetscInt              *stencil_active;
 33:   DMStagStencil         *stencils_ordered_unique;
 34:   PetscInt              *idx, *idxLocal;
 35:   const PetscInt        *ltogidx;
 36:   PetscInt               n_stencil_unique, dim, count, nidx, nc_max;
 37:   ISLocalToGlobalMapping ltog;
 38:   PetscInt               start[DMSTAG_MAX_DIM], n[DMSTAG_MAX_DIM], extraPoint[DMSTAG_MAX_DIM];

 40:   PetscFunctionBegin;
 41:   PetscCall(DMGetDimension(dm, &dim));
 42:   PetscCheck(dim >= 1 && dim <= 3, PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Unsupported dimension %" PetscInt_FMT, dim);

 44:   /* To assert that the resulting IS has unique, sorted, entries, we perform
 45:      a bucket sort, taking advantage of the fact that DMStagStencilLocation
 46:      enum values are integers starting with 1, in canonical order */
 47:   nc_max           = 1; // maximum number of components to represent these stencils
 48:   n_stencil_unique = 0;
 49:   for (PetscInt p = 0; p < n_stencil; ++p) nc_max = PetscMax(nc_max, (stencils[p].c + 1));
 50:   PetscCall(PetscCalloc1(DMSTAG_NUMBER_LOCATIONS * nc_max, &stencil_active));
 51:   for (PetscInt p = 0; p < n_stencil; ++p) {
 52:     DMStagStencilLocation loc_canonical;
 53:     PetscInt              slot;

 55:     PetscCall(DMStagStencilLocationCanonicalize(stencils[p].loc, &loc_canonical));
 56:     slot = nc_max * ((PetscInt)loc_canonical) + stencils[p].c;
 57:     if (stencil_active[slot] == 0) {
 58:       stencil_active[slot] = 1;
 59:       ++n_stencil_unique;
 60:     }
 61:   }
 62:   PetscCall(PetscMalloc1(n_stencil_unique, &stencils_ordered_unique));
 63:   {
 64:     PetscInt p = 0;

 66:     for (PetscInt i = 1; i < DMSTAG_NUMBER_LOCATIONS; ++i) {
 67:       for (PetscInt c = 0; c < nc_max; ++c) {
 68:         if (stencil_active[nc_max * i + c] != 0) {
 69:           stencils_ordered_unique[p].loc = (DMStagStencilLocation)i;
 70:           stencils_ordered_unique[p].c   = c;
 71:           ++p;
 72:         }
 73:       }
 74:     }
 75:   }
 76:   PetscCall(PetscFree(stencil_active));

 78:   PetscCall(PetscMalloc1(n_stencil_unique, &idxLocal));
 79:   PetscCall(DMGetLocalToGlobalMapping(dm, &ltog));
 80:   PetscCall(ISLocalToGlobalMappingGetIndices(ltog, &ltogidx));
 81:   PetscCall(DMStagGetCorners(dm, &start[0], &start[1], &start[2], &n[0], &n[1], &n[2], &extraPoint[0], &extraPoint[1], &extraPoint[2]));
 82:   for (PetscInt d = dim; d < DMSTAG_MAX_DIM; ++d) {
 83:     start[d]      = 0;
 84:     n[d]          = 1; /* To allow for a single loop nest below */
 85:     extraPoint[d] = 0;
 86:   }
 87:   nidx = n_stencil_unique;
 88:   for (PetscInt d = 0; d < dim; ++d) nidx *= (n[d] + 1); /* Overestimate (always assumes extraPoint) */
 89:   PetscCall(PetscMalloc1(nidx, &idx));
 90:   count = 0;
 91:   /* Note that unused loop variables are not accessed, for lower dimensions */
 92:   for (PetscInt k = start[2]; k < start[2] + n[2] + extraPoint[2]; ++k) {
 93:     for (PetscInt j = start[1]; j < start[1] + n[1] + extraPoint[1]; ++j) {
 94:       for (PetscInt i = start[0]; i < start[0] + n[0] + extraPoint[0]; ++i) {
 95:         for (PetscInt p = 0; p < n_stencil_unique; ++p) {
 96:           stencils_ordered_unique[p].i = i;
 97:           stencils_ordered_unique[p].j = j;
 98:           stencils_ordered_unique[p].k = k;
 99:         }
100:         PetscCall(DMStagStencilToIndexLocal(dm, dim, n_stencil_unique, stencils_ordered_unique, idxLocal));
101:         for (PetscInt p = 0; p < n_stencil_unique; ++p) {
102:           const PetscInt gidx = ltogidx[idxLocal[p]];
103:           if (gidx >= 0) {
104:             idx[count] = gidx;
105:             ++count;
106:           }
107:         }
108:       }
109:     }
110:   }

112:   PetscCall(ISLocalToGlobalMappingRestoreIndices(ltog, &ltogidx));
113:   PetscCall(PetscFree(stencils_ordered_unique));
114:   PetscCall(PetscFree(idxLocal));

116:   PetscCall(ISCreateGeneral(PetscObjectComm((PetscObject)dm), count, idx, PETSC_OWN_POINTER, is));
117:   PetscCall(ISSetInfo(*is, IS_SORTED, IS_GLOBAL, PETSC_TRUE, PETSC_TRUE));
118:   PetscCall(ISSetInfo(*is, IS_UNIQUE, IS_GLOBAL, PETSC_TRUE, PETSC_TRUE));
119:   PetscFunctionReturn(PETSC_SUCCESS);
120: }

122: /*@C
123:   DMStagGetLocationDOF - Get number of DOF associated with a given point in a `DMSTAG` grid

125:   Not Collective

127:   Input Parameters:
128: + dm  - the `DMSTAG` object
129: - loc - grid point (see `DMStagStencilLocation`)

131:   Output Parameter:
132: . dof - the number of DOF (components) living at `loc` in `dm`

134:   Level: intermediate

136: .seealso: [](ch_stag), `DMSTAG`, `DMStagStencilLocation`, `DMStagStencil`, `DMDAGetDof()`
137: @*/
138: PetscErrorCode DMStagGetLocationDOF(DM dm, DMStagStencilLocation loc, PetscInt *dof)
139: {
140:   const DM_Stag *const stag = (DM_Stag *)dm->data;
141:   PetscInt             dim;

143:   PetscFunctionBegin;
145:   PetscCall(DMGetDimension(dm, &dim));
146:   switch (dim) {
147:   case 1:
148:     switch (loc) {
149:     case DMSTAG_LEFT:
150:     case DMSTAG_RIGHT:
151:       *dof = stag->dof[0];
152:       break;
153:     case DMSTAG_ELEMENT:
154:       *dof = stag->dof[1];
155:       break;
156:     default:
157:       SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Not implemented for location %s", DMStagStencilLocations[loc]);
158:     }
159:     break;
160:   case 2:
161:     switch (loc) {
162:     case DMSTAG_DOWN_LEFT:
163:     case DMSTAG_DOWN_RIGHT:
164:     case DMSTAG_UP_LEFT:
165:     case DMSTAG_UP_RIGHT:
166:       *dof = stag->dof[0];
167:       break;
168:     case DMSTAG_LEFT:
169:     case DMSTAG_RIGHT:
170:     case DMSTAG_UP:
171:     case DMSTAG_DOWN:
172:       *dof = stag->dof[1];
173:       break;
174:     case DMSTAG_ELEMENT:
175:       *dof = stag->dof[2];
176:       break;
177:     default:
178:       SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Not implemented for location %s", DMStagStencilLocations[loc]);
179:     }
180:     break;
181:   case 3:
182:     switch (loc) {
183:     case DMSTAG_BACK_DOWN_LEFT:
184:     case DMSTAG_BACK_DOWN_RIGHT:
185:     case DMSTAG_BACK_UP_LEFT:
186:     case DMSTAG_BACK_UP_RIGHT:
187:     case DMSTAG_FRONT_DOWN_LEFT:
188:     case DMSTAG_FRONT_DOWN_RIGHT:
189:     case DMSTAG_FRONT_UP_LEFT:
190:     case DMSTAG_FRONT_UP_RIGHT:
191:       *dof = stag->dof[0];
192:       break;
193:     case DMSTAG_BACK_DOWN:
194:     case DMSTAG_BACK_LEFT:
195:     case DMSTAG_BACK_RIGHT:
196:     case DMSTAG_BACK_UP:
197:     case DMSTAG_DOWN_LEFT:
198:     case DMSTAG_DOWN_RIGHT:
199:     case DMSTAG_UP_LEFT:
200:     case DMSTAG_UP_RIGHT:
201:     case DMSTAG_FRONT_DOWN:
202:     case DMSTAG_FRONT_LEFT:
203:     case DMSTAG_FRONT_RIGHT:
204:     case DMSTAG_FRONT_UP:
205:       *dof = stag->dof[1];
206:       break;
207:     case DMSTAG_LEFT:
208:     case DMSTAG_RIGHT:
209:     case DMSTAG_DOWN:
210:     case DMSTAG_UP:
211:     case DMSTAG_BACK:
212:     case DMSTAG_FRONT:
213:       *dof = stag->dof[2];
214:       break;
215:     case DMSTAG_ELEMENT:
216:       *dof = stag->dof[3];
217:       break;
218:     default:
219:       SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Not implemented for location %s", DMStagStencilLocations[loc]);
220:     }
221:     break;
222:   default:
223:     SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_SUP, "Unsupported dimension %" PetscInt_FMT, dim);
224:   }
225:   PetscFunctionReturn(PETSC_SUCCESS);
226: }

228: /*
229: Convert to a location value with only BACK, DOWN, LEFT, and ELEMENT involved
230: */
231: PETSC_INTERN PetscErrorCode DMStagStencilLocationCanonicalize(DMStagStencilLocation loc, DMStagStencilLocation *locCanonical)
232: {
233:   PetscFunctionBegin;
234:   switch (loc) {
235:   case DMSTAG_ELEMENT:
236:     *locCanonical = DMSTAG_ELEMENT;
237:     break;
238:   case DMSTAG_LEFT:
239:   case DMSTAG_RIGHT:
240:     *locCanonical = DMSTAG_LEFT;
241:     break;
242:   case DMSTAG_DOWN:
243:   case DMSTAG_UP:
244:     *locCanonical = DMSTAG_DOWN;
245:     break;
246:   case DMSTAG_BACK:
247:   case DMSTAG_FRONT:
248:     *locCanonical = DMSTAG_BACK;
249:     break;
250:   case DMSTAG_DOWN_LEFT:
251:   case DMSTAG_DOWN_RIGHT:
252:   case DMSTAG_UP_LEFT:
253:   case DMSTAG_UP_RIGHT:
254:     *locCanonical = DMSTAG_DOWN_LEFT;
255:     break;
256:   case DMSTAG_BACK_LEFT:
257:   case DMSTAG_BACK_RIGHT:
258:   case DMSTAG_FRONT_LEFT:
259:   case DMSTAG_FRONT_RIGHT:
260:     *locCanonical = DMSTAG_BACK_LEFT;
261:     break;
262:   case DMSTAG_BACK_DOWN:
263:   case DMSTAG_BACK_UP:
264:   case DMSTAG_FRONT_DOWN:
265:   case DMSTAG_FRONT_UP:
266:     *locCanonical = DMSTAG_BACK_DOWN;
267:     break;
268:   case DMSTAG_BACK_DOWN_LEFT:
269:   case DMSTAG_BACK_DOWN_RIGHT:
270:   case DMSTAG_BACK_UP_LEFT:
271:   case DMSTAG_BACK_UP_RIGHT:
272:   case DMSTAG_FRONT_DOWN_LEFT:
273:   case DMSTAG_FRONT_DOWN_RIGHT:
274:   case DMSTAG_FRONT_UP_LEFT:
275:   case DMSTAG_FRONT_UP_RIGHT:
276:     *locCanonical = DMSTAG_BACK_DOWN_LEFT;
277:     break;
278:   default:
279:     *locCanonical = DMSTAG_NULL_LOCATION;
280:     break;
281:   }
282:   PetscFunctionReturn(PETSC_SUCCESS);
283: }

285: /*@C
286:   DMStagMatGetValuesStencil - retrieve local matrix entries using grid indexing

288:   Not Collective

290:   Input Parameters:
291: + dm     - the `DMSTAG` object
292: . mat    - the matrix
293: . nRow   - number of rows
294: . posRow - grid locations (including components) of rows
295: . nCol   - number of columns
296: - posCol - grid locations (including components) of columns

298:   Output Parameter:
299: . val - logically two-dimensional array of values

301:   Level: advanced

303: .seealso: [](ch_stag), `DMSTAG`, `DMStagStencil`, `DMStagStencilLocation`, `DMStagVecGetValuesStencil()`, `DMStagVecSetValuesStencil()`, `DMStagMatSetValuesStencil()`, `MatSetValuesStencil()`, `MatAssemblyBegin()`, `MatAssemblyEnd()`, `DMCreateMatrix()`
304: @*/
305: PetscErrorCode DMStagMatGetValuesStencil(DM dm, Mat mat, PetscInt nRow, const DMStagStencil *posRow, PetscInt nCol, const DMStagStencil *posCol, PetscScalar *val)
306: {
307:   PetscInt  dim;
308:   PetscInt *ir, *ic;

310:   PetscFunctionBegin;
313:   PetscCall(DMGetDimension(dm, &dim));
314:   PetscCall(PetscMalloc2(nRow, &ir, nCol, &ic));
315:   PetscCall(DMStagStencilToIndexLocal(dm, dim, nRow, posRow, ir));
316:   PetscCall(DMStagStencilToIndexLocal(dm, dim, nCol, posCol, ic));
317:   PetscCall(MatGetValuesLocal(mat, nRow, ir, nCol, ic, val));
318:   PetscCall(PetscFree2(ir, ic));
319:   PetscFunctionReturn(PETSC_SUCCESS);
320: }

322: /*@C
323:   DMStagMatSetValuesStencil - insert or add matrix entries using grid indexing

325:   Not Collective

327:   Input Parameters:
328: + dm         - the `DMSTAG` object
329: . mat        - the matrix
330: . nRow       - number of rows
331: . posRow     - grid locations (including components) of rows
332: . nCol       - number of columns
333: . posCol     - grid locations (including components) of columns
334: . val        - logically two-dimensional array of values
335: - insertMode - `INSERT_VALUES` or `ADD_VALUES`

337:   Notes:
338:   See notes for `MatSetValuesStencil()`

340:   Level: intermediate

342: .seealso: [](ch_stag), `DMSTAG`, `DMStagStencil`, `DMStagStencilLocation`, `DMStagVecGetValuesStencil()`, `DMStagVecSetValuesStencil()`, `DMStagMatGetValuesStencil()`, `MatSetValuesStencil()`, `MatAssemblyBegin()`, `MatAssemblyEnd()`, `DMCreateMatrix()`
343: @*/
344: PetscErrorCode DMStagMatSetValuesStencil(DM dm, Mat mat, PetscInt nRow, const DMStagStencil *posRow, PetscInt nCol, const DMStagStencil *posCol, const PetscScalar *val, InsertMode insertMode)
345: {
346:   PetscInt *ir, *ic;

348:   PetscFunctionBegin;
351:   PetscCall(PetscMalloc2(nRow, &ir, nCol, &ic));
352:   PetscCall(DMStagStencilToIndexLocal(dm, dm->dim, nRow, posRow, ir));
353:   PetscCall(DMStagStencilToIndexLocal(dm, dm->dim, nCol, posCol, ic));
354:   PetscCall(MatSetValuesLocal(mat, nRow, ir, nCol, ic, val, insertMode));
355:   PetscCall(PetscFree2(ir, ic));
356:   PetscFunctionReturn(PETSC_SUCCESS);
357: }

359: /*@C
360:   DMStagStencilToIndexLocal - Convert an array of `DMStagStenci`l objects to an array of indices into a local vector.

362:   Not Collective

364:   Input Parameters:
365: + dm  - the `DMSTAG` object
366: . dim - the dimension of the `DMSTAG` object
367: . n   - the number of `DMStagStencil` objects
368: - pos - an array of `n` `DMStagStencil` objects

370:   Output Parameter:
371: . ix - output array of `n` indices

373:   Notes:
374:   The `DMStagStencil` objects in `pos` use global element indices.

376:   The `.c` fields in `pos` must always be set (even if to `0`).

378:   Developer Notes:
379:   This is a "hot" function, and accepts the dimension redundantly to avoid having to perform any error checking inside the function.

381:   Level: developer

383: .seealso: [](ch_stag), `DMSTAG`, `DMStagStencilLocation`, `DMStagStencil`, `DMGetLocalVector`, `DMCreateLocalVector`
384: @*/
385: PetscErrorCode DMStagStencilToIndexLocal(DM dm, PetscInt dim, PetscInt n, const DMStagStencil *pos, PetscInt *ix)
386: {
387:   const DM_Stag *const stag = (DM_Stag *)dm->data;
388:   const PetscInt       epe  = stag->entriesPerElement;

390:   PetscFunctionBeginHot;
391:   if (dim == 1) {
392:     for (PetscInt idx = 0; idx < n; ++idx) {
393:       const PetscInt eLocal = pos[idx].i - stag->startGhost[0];

395:       ix[idx] = eLocal * epe + stag->locationOffsets[pos[idx].loc] + pos[idx].c;
396:     }
397:   } else if (dim == 2) {
398:     const PetscInt epr = stag->nGhost[0];

400:     for (PetscInt idx = 0; idx < n; ++idx) {
401:       const PetscInt eLocalx = pos[idx].i - stag->startGhost[0];
402:       const PetscInt eLocaly = pos[idx].j - stag->startGhost[1];
403:       const PetscInt eLocal  = eLocalx + epr * eLocaly;

405:       ix[idx] = eLocal * epe + stag->locationOffsets[pos[idx].loc] + pos[idx].c;
406:     }
407:   } else if (dim == 3) {
408:     const PetscInt epr = stag->nGhost[0];
409:     const PetscInt epl = stag->nGhost[0] * stag->nGhost[1];

411:     for (PetscInt idx = 0; idx < n; ++idx) {
412:       const PetscInt eLocalx = pos[idx].i - stag->startGhost[0];
413:       const PetscInt eLocaly = pos[idx].j - stag->startGhost[1];
414:       const PetscInt eLocalz = pos[idx].k - stag->startGhost[2];
415:       const PetscInt eLocal  = epl * eLocalz + epr * eLocaly + eLocalx;

417:       ix[idx] = eLocal * epe + stag->locationOffsets[pos[idx].loc] + pos[idx].c;
418:     }
419:   } else SETERRQ(PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_OUTOFRANGE, "Unsupported dimension %" PetscInt_FMT, dim);
420:   PetscFunctionReturn(PETSC_SUCCESS);
421: }

423: /*@C
424:   DMStagVecGetValuesStencil - get vector values using grid indexing

426:   Not Collective

428:   Input Parameters:
429: + dm  - the `DMSTAG` object
430: . vec - the vector object
431: . n   - the number of values to obtain
432: - pos - locations to obtain values from (as an array of `DMStagStencil` values)

434:   Output Parameter:
435: . val - value at the point

437:   Notes:
438:   Accepts stencils which refer to global element numbers, but
439:   only allows access to entries in the local representation (including ghosts).

441:   This approach is not as efficient as getting values directly with `DMStagVecGetArray()`,
442:   which is recommended for matrix-free operators.

444:   Level: advanced

446: .seealso: [](ch_stag), `DMSTAG`, `DMStagStencil`, `DMStagStencilLocation`, `DMStagVecSetValuesStencil()`, `DMStagMatSetValuesStencil()`, `DMStagVecGetArray()`
447: @*/
448: PetscErrorCode DMStagVecGetValuesStencil(DM dm, Vec vec, PetscInt n, const DMStagStencil *pos, PetscScalar *val)
449: {
450:   DM_Stag *const     stag = (DM_Stag *)dm->data;
451:   PetscInt           nLocal, idx;
452:   PetscInt          *ix;
453:   PetscScalar const *arr;

455:   PetscFunctionBegin;
458:   PetscCall(VecGetLocalSize(vec, &nLocal));
459:   PetscCheck(nLocal == stag->entriesGhost, PETSC_COMM_SELF, PETSC_ERR_ARG_OUTOFRANGE, "Vector should be a local vector. Local size %" PetscInt_FMT " does not match expected %" PetscInt_FMT, nLocal, stag->entriesGhost);
460:   PetscCall(PetscMalloc1(n, &ix));
461:   PetscCall(DMStagStencilToIndexLocal(dm, dm->dim, n, pos, ix));
462:   PetscCall(VecGetArrayRead(vec, &arr));
463:   for (idx = 0; idx < n; ++idx) val[idx] = arr[ix[idx]];
464:   PetscCall(VecRestoreArrayRead(vec, &arr));
465:   PetscCall(PetscFree(ix));
466:   PetscFunctionReturn(PETSC_SUCCESS);
467: }

469: /*@C
470:   DMStagVecSetValuesStencil - Set `Vec` values using global grid indexing

472:   Not Collective

474:   Input Parameters:
475: + dm         - the `DMSTAG` object
476: . vec        - the `Vec`
477: . n          - the number of values to set
478: . pos        - the locations to set values, as an array of `DMStagStencil` structs
479: . val        - the values to set
480: - insertMode - `INSERT_VALUES` or `ADD_VALUES`

482:   Notes:
483:   The vector is expected to be a global vector compatible with the DM (usually obtained by `DMGetGlobalVector()` or `DMCreateGlobalVector()`).

485:   This approach is not as efficient as setting values directly with `DMStagVecGetArray()`, which is recommended for matrix-free operators.
486:   For assembling systems, where overhead may be less important than convenience, this routine could be helpful in assembling a righthand side and a matrix (using `DMStagMatSetValuesStencil()`).

488:   Level: advanced

490: .seealso: [](ch_stag), `DMSTAG`, `Vec`, `DMStagStencil`, `DMStagStencilLocation`, `DMStagVecGetValuesStencil()`, `DMStagMatSetValuesStencil()`, `DMCreateGlobalVector()`, `DMGetLocalVector()`, `DMStagVecGetArray()`
491: @*/
492: PetscErrorCode DMStagVecSetValuesStencil(DM dm, Vec vec, PetscInt n, const DMStagStencil *pos, const PetscScalar *val, InsertMode insertMode)
493: {
494:   DM_Stag *const stag = (DM_Stag *)dm->data;
495:   PetscInt       nLocal;
496:   PetscInt      *ix;

498:   PetscFunctionBegin;
501:   PetscCall(VecGetLocalSize(vec, &nLocal));
502:   PetscCheck(nLocal == stag->entries, PetscObjectComm((PetscObject)dm), PETSC_ERR_ARG_WRONG, "Provided vec has a different number of local entries (%" PetscInt_FMT ") than expected (%" PetscInt_FMT "). It should be a global vector", nLocal, stag->entries);
503:   PetscCall(PetscMalloc1(n, &ix));
504:   PetscCall(DMStagStencilToIndexLocal(dm, dm->dim, n, pos, ix));
505:   PetscCall(VecSetValuesLocal(vec, n, ix, val, insertMode));
506:   PetscCall(PetscFree(ix));
507:   PetscFunctionReturn(PETSC_SUCCESS);
508: }