Actual source code: dagetarray.c


  2: #include <petsc/private/dmdaimpl.h>

  4: /*MC
  5:   DMDAVecGetArrayF90 - check Fortran Notes at `DMDAVecGetArray()`

  7:   Level: intermediate
  8: M*/

 10: /*@C
 11:    DMDAVecGetArray - Returns a multiple dimension array that shares data with
 12:       the underlying vector and is indexed using the global dimensions.

 14:    Logically Collective

 16:    Input Parameters:
 17: +  da - the distributed array
 18: -  vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`

 20:    Output Parameter:
 21: .  array - the array

 23:   Level: intermediate

 25:    Notes:
 26:     Call `DMDAVecRestoreArray()` once you have finished accessing the vector entries.

 28:     In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]!

 30:     If `vec` is a local vector (obtained with DMCreateLocalVector() etc) then the ghost point locations are accessible. If it is
 31:     a global vector then the ghost points are not accessible. Of course, with a local vector you will have had to do the
 32:     appropriate `DMGlobalToLocalBegin()` and `DMGlobalToLocalEnd()` to have correct values in the ghost locations.

 34:     The accessible indices are `array[zs:zs+zm-1][ys:ys+ym-1][xs:xs+xm-1]` where the values are obtained from
 35:    `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.

 37:   Fortran Notes:
 38:   Use `DMDAVecGetArrayF90()` and pass for the array type `PetscScalar`,pointer :: array(:,...,:) of the appropriate
 39:   dimension. For a `DMDA` created with a dof of 1 use the dimension of the `DMDA`, for a `DMDA` created with a dof greater than 1 use one more than the
 40:   dimension of the `DMDA`.

 42:   The order of the indices is `array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` (when dof is 1) otherwise
 43:   `array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` where the values are obtained from
 44:   `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.

 46: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArray()`, `DMDAVecRestoreArrayDOF()`
 47:           `DMDAVecGetArrayDOF()`, `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`,
 48:           `DMStagVecGetArray()`
 49: @*/
 50: PetscErrorCode DMDAVecGetArray(DM da, Vec vec, void *array)
 51: {
 52:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

 54:   PetscFunctionBegin;
 58:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
 59:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
 60:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

 62:   /* Handle case where user passes in global vector as opposed to local */
 63:   PetscCall(VecGetLocalSize(vec, &N));
 64:   if (N == xm * ym * zm * dof) {
 65:     gxm = xm;
 66:     gym = ym;
 67:     gzm = zm;
 68:     gxs = xs;
 69:     gys = ys;
 70:     gzs = zs;
 71:   } else PetscCheck(N == gxm * gym * gzm * dof, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Vector local size %" PetscInt_FMT " is not compatible with DMDA local sizes %" PetscInt_FMT " %" PetscInt_FMT, N, xm * ym * zm * dof, gxm * gym * gzm * dof);

 73:   if (dim == 1) {
 74:     PetscCall(VecGetArray1d(vec, gxm * dof, gxs * dof, (PetscScalar **)array));
 75:   } else if (dim == 2) {
 76:     PetscCall(VecGetArray2d(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array));
 77:   } else if (dim == 3) {
 78:     PetscCall(VecGetArray3d(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array));
 79:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
 80:   PetscFunctionReturn(PETSC_SUCCESS);
 81: }

 83: /*MC
 84:   DMDAVecRestoreArrayF90 - check Fortran Notes at `DMDAVecRestoreArray()`

 86:   Level: intermediate
 87: M*/

 89: /*@
 90:    DMDAVecRestoreArray - Restores a multiple dimension array obtained with `DMDAVecGetArray()`

 92:    Logically Collective

 94:    Input Parameters:
 95: +  da - the distributed array
 96: .  vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
 97: -  array - the `array` pointer is zeroed

 99:   Level: intermediate

101:   Fortran Note:
102:   Use `DMDAVecRestoreArayF90()`

104: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `DMDAVecRestoreArayF90()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArray()`,
105:           `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`,
106:           `DMDStagVecRestoreArray()`
107: @*/
108: PetscErrorCode DMDAVecRestoreArray(DM da, Vec vec, void *array)
109: {
110:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

112:   PetscFunctionBegin;
116:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
117:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
118:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

120:   /* Handle case where user passes in global vector as opposed to local */
121:   PetscCall(VecGetLocalSize(vec, &N));
122:   if (N == xm * ym * zm * dof) {
123:     gxm = xm;
124:     gym = ym;
125:     gzm = zm;
126:     gxs = xs;
127:     gys = ys;
128:     gzs = zs;
129:   } else PetscCheck(N == gxm * gym * gzm * dof, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Vector local size %" PetscInt_FMT " is not compatible with DMDA local sizes %" PetscInt_FMT " %" PetscInt_FMT, N, xm * ym * zm * dof, gxm * gym * gzm * dof);

131:   if (dim == 1) {
132:     PetscCall(VecRestoreArray1d(vec, gxm * dof, gxs * dof, (PetscScalar **)array));
133:   } else if (dim == 2) {
134:     PetscCall(VecRestoreArray2d(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array));
135:   } else if (dim == 3) {
136:     PetscCall(VecRestoreArray3d(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array));
137:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
138:   PetscFunctionReturn(PETSC_SUCCESS);
139: }

141: /*MC
142:   DMDAVecGetArrayWriteF90 - check Fortran Notes at `DMDAVecGetArrayWrite()`

144:   Level: intermediate
145: M*/

147: /*@C
148:    DMDAVecGetArrayWrite - Returns a multiple dimension array that shares data with
149:       the underlying vector and is indexed using the global dimensions.

151:    Logically Collective

153:    Input Parameters:
154: +  da - the distributed array
155: -  vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`

157:    Output Parameter:
158: .  array - the array

160:   Level: intermediate

162:    Notes:
163:     Call `DMDAVecRestoreArray()` once you have finished accessing the vector entries.

165:     In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]!

167:     if `vec` is a local vector (obtained with `DMCreateLocalVector()` etc) then the ghost point locations are accessible. If it is
168:     a global vector then the ghost points are not accessible. Of course with the local vector you will have had to do the
169:     appropriate `DMGlobalToLocalBegin()` and `DMGlobalToLocalEnd()` to have correct values in the ghost locations.

171:      The accessible indices are `array[zs:zs+zm-1][ys:ys+ym-1][xs:xs+xm-1]` where the values are obtained from
172:    `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.

174:    Fortran Notes:
175:    From Fortran use `DMDAVecGetArrayWriteF90()` and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
176:    dimension. For a `DMDA` created with a dof of 1 use the dimension of the `DMDA`, for a `DMDA` created with a dof greater than 1 use one more than the
177:    dimension of the `DMDA`.

179:    The order of the indices is `array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` (when dof is 1) otherwise
180:    `array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` where the values are obtained from
181:    `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.

183:   Developer Note:
184:   This has code duplication with `DMDAVecGetArray()` and `DMDAVecGetArrayRead()`

186: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecRestoreArrayDOF()`
187:           `DMDAVecGetArrayDOF()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`
188: @*/
189: PetscErrorCode DMDAVecGetArrayWrite(DM da, Vec vec, void *array)
190: {
191:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

193:   PetscFunctionBegin;
197:   if (da->localSection) {
198:     PetscCall(VecGetArrayWrite(vec, (PetscScalar **)array));
199:     PetscFunctionReturn(PETSC_SUCCESS);
200:   }
201:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
202:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
203:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

205:   /* Handle case where user passes in global vector as opposed to local */
206:   PetscCall(VecGetLocalSize(vec, &N));
207:   if (N == xm * ym * zm * dof) {
208:     gxm = xm;
209:     gym = ym;
210:     gzm = zm;
211:     gxs = xs;
212:     gys = ys;
213:     gzs = zs;
214:   } else PetscCheck(N == gxm * gym * gzm * dof, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Vector local size %" PetscInt_FMT " is not compatible with DMDA local sizes %" PetscInt_FMT " %" PetscInt_FMT, N, xm * ym * zm * dof, gxm * gym * gzm * dof);

216:   if (dim == 1) {
217:     PetscCall(VecGetArray1dWrite(vec, gxm * dof, gxs * dof, (PetscScalar **)array));
218:   } else if (dim == 2) {
219:     PetscCall(VecGetArray2dWrite(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array));
220:   } else if (dim == 3) {
221:     PetscCall(VecGetArray3dWrite(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array));
222:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
223:   PetscFunctionReturn(PETSC_SUCCESS);
224: }

226: /*MC
227:   DMDAVecRestoreArrayWriteF90 - check Fortran Notes at `DMDAVecRestoreArrayWrite()`

229:   Level: intermediate
230: M*/

232: /*@
233:    DMDAVecRestoreArrayWrite - Restores a multiple dimension array obtained with `DMDAVecGetArrayWrite()`

235:    Logically Collective

237:    Input Parameters:
238: +  da - the distributed array
239: .  vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
240: -  array - the `array` pointer is zeroed

242:   Level: intermediate

244:   Fortran Note:
245:   Use `DMDAVecRestoreArayWriteF90()`

247: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArrayWrite()`,
248:           `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`
249: @*/
250: PetscErrorCode DMDAVecRestoreArrayWrite(DM da, Vec vec, void *array)
251: {
252:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

254:   PetscFunctionBegin;
258:   if (da->localSection) {
259:     PetscCall(VecRestoreArray(vec, (PetscScalar **)array));
260:     PetscFunctionReturn(PETSC_SUCCESS);
261:   }
262:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
263:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
264:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

266:   /* Handle case where user passes in global vector as opposed to local */
267:   PetscCall(VecGetLocalSize(vec, &N));
268:   if (N == xm * ym * zm * dof) {
269:     gxm = xm;
270:     gym = ym;
271:     gzm = zm;
272:     gxs = xs;
273:     gys = ys;
274:     gzs = zs;
275:   } else PetscCheck(N == gxm * gym * gzm * dof, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Vector local size %" PetscInt_FMT " is not compatible with DMDA local sizes %" PetscInt_FMT " %" PetscInt_FMT, N, xm * ym * zm * dof, gxm * gym * gzm * dof);

277:   if (dim == 1) {
278:     PetscCall(VecRestoreArray1dWrite(vec, gxm * dof, gxs * dof, (PetscScalar **)array));
279:   } else if (dim == 2) {
280:     PetscCall(VecRestoreArray2dWrite(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array));
281:   } else if (dim == 3) {
282:     PetscCall(VecRestoreArray3dWrite(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array));
283:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
284:   PetscFunctionReturn(PETSC_SUCCESS);
285: }

287: /*@C
288:    DMDAVecGetArrayDOF - Returns a multiple dimension array that shares data with
289:       the underlying vector and is indexed using the global dimensions.

291:    Logically Collective

293:    Input Parameters:
294: +  da - the distributed array
295: -  vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`

297:    Output Parameter:
298: .  array - the `array` pointer is zeroed

300:   Level: intermediate

302:    Notes:
303:     Call `DMDAVecRestoreArrayDOF()` once you have finished accessing the vector entries.

305:     In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]

307:     The accessible indices are `array[zs:zs+zm-1][ys:ys+ym-1][xs:xs+xm-1][0:ndof-1]` where the values are obtained from
308:    `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.

310:    Fortran Notes:
311:     Use `DMDAVecGetArrayF90()` and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
312:    dimension. For a `DMDA` created with a dof of 1 use the dimension of the `DMDA`, for a `DMDA` created with a dof greater than 1 use one more than the
313:    dimension of the `DMDA`.

315:    The order of the indices is `array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` (when ndof is 1) otherwise
316:    `array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` where the values are obtained from
317:    `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.

319: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecRestoreArrayDOF()`,
320:           `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`, `DMDAVecGetArrayDOFRead()`
321: @*/
322: PetscErrorCode DMDAVecGetArrayDOF(DM da, Vec vec, void *array)
323: {
324:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

326:   PetscFunctionBegin;
327:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
328:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
329:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

331:   /* Handle case where user passes in global vector as opposed to local */
332:   PetscCall(VecGetLocalSize(vec, &N));
333:   if (N == xm * ym * zm * dof) {
334:     gxm = xm;
335:     gym = ym;
336:     gzm = zm;
337:     gxs = xs;
338:     gys = ys;
339:     gzs = zs;
340:   } else PetscCheck(N == gxm * gym * gzm * dof, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Vector local size %" PetscInt_FMT " is not compatible with DMDA local sizes %" PetscInt_FMT " %" PetscInt_FMT, N, xm * ym * zm * dof, gxm * gym * gzm * dof);

342:   if (dim == 1) {
343:     PetscCall(VecGetArray2d(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
344:   } else if (dim == 2) {
345:     PetscCall(VecGetArray3d(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
346:   } else if (dim == 3) {
347:     PetscCall(VecGetArray4d(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
348:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
349:   PetscFunctionReturn(PETSC_SUCCESS);
350: }

352: /*@
353:    DMDAVecRestoreArrayDOF - Restores a multiple dimension array obtained with `DMDAVecGetArrayDOF()`

355:    Logically Collective

357:    Input Parameters:
358: +  da - the distributed array
359: .  vec - vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
360: -  array - the `array` point is zeroed

362:   Level: intermediate

364:   Fortran Note:
365:   Use `DMDAVecRestoreArayF90()`

367: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`, `DMDAVecRestoreArrayDOF()`,
368:           `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`
369: @*/
370: PetscErrorCode DMDAVecRestoreArrayDOF(DM da, Vec vec, void *array)
371: {
372:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

374:   PetscFunctionBegin;
375:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
376:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
377:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

379:   /* Handle case where user passes in global vector as opposed to local */
380:   PetscCall(VecGetLocalSize(vec, &N));
381:   if (N == xm * ym * zm * dof) {
382:     gxm = xm;
383:     gym = ym;
384:     gzm = zm;
385:     gxs = xs;
386:     gys = ys;
387:     gzs = zs;
388:   }

390:   if (dim == 1) {
391:     PetscCall(VecRestoreArray2d(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
392:   } else if (dim == 2) {
393:     PetscCall(VecRestoreArray3d(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
394:   } else if (dim == 3) {
395:     PetscCall(VecRestoreArray4d(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
396:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
397:   PetscFunctionReturn(PETSC_SUCCESS);
398: }

400: /*MC
401:   DMDAVecGetArrayReadF90 - check Fortran Notes at `DMDAVecGetArrayRead()`

403:   Level: intermediate
404: M*/

406: /*@C
407:    DMDAVecGetArrayRead - Returns a multiple dimension array that shares data with
408:       the underlying vector and is indexed using the global dimensions.

410:    Not Collective

412:    Input Parameters:
413: +  da - the distributed array
414: -  vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`

416:    Output Parameter:
417: .  array - the array

419:   Level: intermediate

421:    Notes:
422:     Call `DMDAVecRestoreArrayRead()` once you have finished accessing the vector entries.

424:     In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]!

426:     If `vec` is a local vector (obtained with `DMCreateLocalVector()` etc) then the ghost point locations are accessible. If it is
427:     a global vector then the ghost points are not accessible. Of course with the local vector you will have had to do the
428:     appropriate `DMGlobalToLocalBegin()` and `DMGlobalToLocalEnd()` to have correct values in the ghost locations.

430:     The accessible indices are `array[zs:zs+zm-1][ys:ys+ym-1][xs:xs+xm-1]` where the values are obtained from
431:     `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.

433:   Fortran Notes:
434:   Use `DMDAVecGetArrayReadF90()` and pass for the array type `PetscScalar`,pointer :: array(:,...,:) of the appropriate
435:   dimension. For a `DMDA` created with a dof of 1 use the dimension of the `DMDA`, for a `DMDA` created with a dof greater than 1 use one more than the
436:   dimension of the `DMDA`.

438:   The order of the indices is `array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` (when dof is 1) otherwise
439:   `array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` where the values are obtained from
440:   `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.

442: .seealso: `DM`, `DMDA`, `DMDAVecGetArrayReadF90()`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArrayRead()`, `DMDAVecRestoreArrayDOF()`
443:           `DMDAVecGetArrayDOF()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`,
444:           `DMStagVecGetArrayRead()`
445: @*/
446: PetscErrorCode DMDAVecGetArrayRead(DM da, Vec vec, void *array)
447: {
448:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

450:   PetscFunctionBegin;
454:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
455:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
456:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

458:   /* Handle case where user passes in global vector as opposed to local */
459:   PetscCall(VecGetLocalSize(vec, &N));
460:   if (N == xm * ym * zm * dof) {
461:     gxm = xm;
462:     gym = ym;
463:     gzm = zm;
464:     gxs = xs;
465:     gys = ys;
466:     gzs = zs;
467:   } else PetscCheck(N == gxm * gym * gzm * dof, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Vector local size %" PetscInt_FMT " is not compatible with DMDA local sizes %" PetscInt_FMT " %" PetscInt_FMT, N, xm * ym * zm * dof, gxm * gym * gzm * dof);

469:   if (dim == 1) {
470:     PetscCall(VecGetArray1dRead(vec, gxm * dof, gxs * dof, (PetscScalar **)array));
471:   } else if (dim == 2) {
472:     PetscCall(VecGetArray2dRead(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array));
473:   } else if (dim == 3) {
474:     PetscCall(VecGetArray3dRead(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array));
475:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
476:   PetscFunctionReturn(PETSC_SUCCESS);
477: }

479: /*MC
480:   DMDAVecRestoreArrayReadF90 - check Fortran Notes at `DMDAVecRestoreArrayRead()`

482:   Level: intermediate
483: M*/

485: /*@
486:    DMDAVecRestoreArrayRead - Restores a multiple dimension array obtained with `DMDAVecGetArrayRead()`

488:    Not Collective

490:    Input Parameters:
491: +  da - the distributed array
492: .  vec - vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
493: -  array - the `array` pointer is zeroed

495:   Level: intermediate

497:   Fortran Note:
498:   Use `DMDAVecRestoreArrayReadF90()`

500: .seealso: `DM`, `DMDA`, `DMDAVecRestoreArrayReadF90()`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArrayRead()`,
501:           `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`,
502:           `DMStagVecRestoreArrayRead()`
503: @*/
504: PetscErrorCode DMDAVecRestoreArrayRead(DM da, Vec vec, void *array)
505: {
506:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

508:   PetscFunctionBegin;
512:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
513:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
514:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

516:   /* Handle case where user passes in global vector as opposed to local */
517:   PetscCall(VecGetLocalSize(vec, &N));
518:   if (N == xm * ym * zm * dof) {
519:     gxm = xm;
520:     gym = ym;
521:     gzm = zm;
522:     gxs = xs;
523:     gys = ys;
524:     gzs = zs;
525:   } else PetscCheck(N == gxm * gym * gzm * dof, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Vector local size %" PetscInt_FMT " is not compatible with DMDA local sizes %" PetscInt_FMT " %" PetscInt_FMT, N, xm * ym * zm * dof, gxm * gym * gzm * dof);

527:   if (dim == 1) {
528:     PetscCall(VecRestoreArray1dRead(vec, gxm * dof, gxs * dof, (PetscScalar **)array));
529:   } else if (dim == 2) {
530:     PetscCall(VecRestoreArray2dRead(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array));
531:   } else if (dim == 3) {
532:     PetscCall(VecRestoreArray3dRead(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array));
533:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
534:   PetscFunctionReturn(PETSC_SUCCESS);
535: }

537: /*@C
538:    DMDAVecGetArrayDOFRead - Returns a multiple dimension array that shares data with
539:       the underlying vector and is indexed using the global dimensions.

541:    Not Collective

543:    Input Parameters:
544: +  da - the distributed array
545: -  vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`

547:    Output Parameter:
548: .  array - the array

550:   Level: intermediate

552:    Notes:
553:    Call `DMDAVecRestoreArrayDOFRead()` once you have finished accessing the vector entries.

555:    In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]!

557:    The accessible indices are `array[zs:zs+zm-1][ys:ys+ym-1][xs:xs+xm-1]` where the values are obtained from
558:    `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.

560:    Fortran Note:
561:    Use  `DMDAVecGetArrayReadF90()` and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
562:    dimension. For a `DMDA` created with a dof of 1 use the dimension of the `DMDA`, for a `DMDA` created with a dof greater than 1 use one more than the
563:    dimension of the `DMDA`.

565:    The order of the indices is `array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` (when dof is 1) otherwise
566:    `array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` where the values are obtained from
567:    `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.

569: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`,
570:           `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`
571: @*/
572: PetscErrorCode DMDAVecGetArrayDOFRead(DM da, Vec vec, void *array)
573: {
574:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

576:   PetscFunctionBegin;
577:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
578:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
579:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

581:   /* Handle case where user passes in global vector as opposed to local */
582:   PetscCall(VecGetLocalSize(vec, &N));
583:   if (N == xm * ym * zm * dof) {
584:     gxm = xm;
585:     gym = ym;
586:     gzm = zm;
587:     gxs = xs;
588:     gys = ys;
589:     gzs = zs;
590:   } else PetscCheck(N == gxm * gym * gzm * dof, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Vector local size %" PetscInt_FMT " is not compatible with DMDA local sizes %" PetscInt_FMT " %" PetscInt_FMT, N, xm * ym * zm * dof, gxm * gym * gzm * dof);

592:   if (dim == 1) {
593:     PetscCall(VecGetArray2dRead(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
594:   } else if (dim == 2) {
595:     PetscCall(VecGetArray3dRead(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
596:   } else if (dim == 3) {
597:     PetscCall(VecGetArray4dRead(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
598:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
599:   PetscFunctionReturn(PETSC_SUCCESS);
600: }

602: /*@
603:    DMDAVecRestoreArrayDOFRead - Restores a multiple dimension array obtained with `DMDAVecGetArrayDOFRead()`

605:    Not Collective

607:    Input Parameters:
608: +  da - the distributed array
609: .  vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
610: -  array - the `array` pointer is zeroed

612:   Level: intermediate

614:   Fortran Note:
615:   Use `DMDAVecRestoreArrayReadF90()`

617: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`, `DMDAVecRestoreArrayDOF()`,
618:           `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`
619: @*/
620: PetscErrorCode DMDAVecRestoreArrayDOFRead(DM da, Vec vec, void *array)
621: {
622:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

624:   PetscFunctionBegin;
625:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
626:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
627:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

629:   /* Handle case where user passes in global vector as opposed to local */
630:   PetscCall(VecGetLocalSize(vec, &N));
631:   if (N == xm * ym * zm * dof) {
632:     gxm = xm;
633:     gym = ym;
634:     gzm = zm;
635:     gxs = xs;
636:     gys = ys;
637:     gzs = zs;
638:   }

640:   if (dim == 1) {
641:     PetscCall(VecRestoreArray2dRead(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
642:   } else if (dim == 2) {
643:     PetscCall(VecRestoreArray3dRead(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
644:   } else if (dim == 3) {
645:     PetscCall(VecRestoreArray4dRead(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
646:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
647:   PetscFunctionReturn(PETSC_SUCCESS);
648: }

650: /*@C
651:    DMDAVecGetArrayDOFWrite - Returns a multiple dimension array that shares data with
652:       the underlying vector and is indexed using the global dimensions.

654:    Not Collective

656:    Input Parameters:
657: +  da - the distributed array
658: -  vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`

660:    Output Parameter:
661: .  array - the array

663:    Level: intermediate

665:    Notes:
666:     Call `DMDAVecRestoreArrayDOFWrite()` once you have finished accessing the vector entries.

668:     In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]!

670:    The accessible indices are `array[zs:zs+zm-1][ys:ys+ym-1][xs:xs+xm-1][0:dof-1]` where the values are obtained from
671:    `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.

673:    Fortran Note:
674:    Use  `DMDAVecGetArrayWriteF90()` and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
675:    dimension. For a `DMDA` created with a dof of 1 use the dimension of the `DMDA`, for a `DMDA` created with a dof greater than 1 use one more than the
676:    dimension of the `DMDA`.

678:    The order of the indices is `array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` (when dof is 1) otherwise
679:    `array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` where the values are obtained from
680:    `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.

682: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `DMDAVecRestoreArrayWriteF90()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`,
683:           `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`
684: @*/
685: PetscErrorCode DMDAVecGetArrayDOFWrite(DM da, Vec vec, void *array)
686: {
687:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

689:   PetscFunctionBegin;
690:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
691:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
692:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

694:   /* Handle case where user passes in global vector as opposed to local */
695:   PetscCall(VecGetLocalSize(vec, &N));
696:   if (N == xm * ym * zm * dof) {
697:     gxm = xm;
698:     gym = ym;
699:     gzm = zm;
700:     gxs = xs;
701:     gys = ys;
702:     gzs = zs;
703:   } else PetscCheck(N == gxm * gym * gzm * dof, PETSC_COMM_SELF, PETSC_ERR_ARG_INCOMP, "Vector local size %" PetscInt_FMT " is not compatible with DMDA local sizes %" PetscInt_FMT " %" PetscInt_FMT, N, xm * ym * zm * dof, gxm * gym * gzm * dof);

705:   if (dim == 1) {
706:     PetscCall(VecGetArray2dWrite(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
707:   } else if (dim == 2) {
708:     PetscCall(VecGetArray3dWrite(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
709:   } else if (dim == 3) {
710:     PetscCall(VecGetArray4dWrite(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
711:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
712:   PetscFunctionReturn(PETSC_SUCCESS);
713: }

715: /*@
716:    DMDAVecRestoreArrayDOFWrite - Restores a multiple dimension array obtained with `DMDAVecGetArrayDOFWrite()`

718:    Not Collective

720:    Input Parameters:
721: +  da - the distributed array
722: .  vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
723: -  array - the `array` pointer is zeroed

725:   Level: intermediate

727:   Fortran Note:
728:   Use `DMDAVecRestoreArrayWriteF90()`

730: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`, `DMDAVecRestoreArrayDOF()`,
731:           `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`
732: @*/
733: PetscErrorCode DMDAVecRestoreArrayDOFWrite(DM da, Vec vec, void *array)
734: {
735:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

737:   PetscFunctionBegin;
738:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
739:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
740:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

742:   /* Handle case where user passes in global vector as opposed to local */
743:   PetscCall(VecGetLocalSize(vec, &N));
744:   if (N == xm * ym * zm * dof) {
745:     gxm = xm;
746:     gym = ym;
747:     gzm = zm;
748:     gxs = xs;
749:     gys = ys;
750:     gzs = zs;
751:   }

753:   if (dim == 1) {
754:     PetscCall(VecRestoreArray2dWrite(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
755:   } else if (dim == 2) {
756:     PetscCall(VecRestoreArray3dWrite(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
757:   } else if (dim == 3) {
758:     PetscCall(VecRestoreArray4dWrite(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
759:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
760:   PetscFunctionReturn(PETSC_SUCCESS);
761: }