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 on da

 16:    Input Parameters:
 17: +  da - the distributed array
 18: -  vec - the vector, either 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 the local vector you will have had to do the

 33:     appropriate `DMGlobalToLocalBegin()` and `DMGlobalToLocalEnd()` to have correct values in the ghost locations.

 35:   Fortran Notes:
 36:   Use `DMDAVecGetArrayF90()` and pass for the array type `PetscScalar`,pointer :: array(:,...,:) of the appropriate
 37:   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
 38:   dimension of the `DMDA`.

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

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

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

 60:   /* Handle case where user passes in global vector as opposed to local */
 61:   PetscCall(VecGetLocalSize(vec, &N));
 62:   if (N == xm * ym * zm * dof) {
 63:     gxm = xm;
 64:     gym = ym;
 65:     gzm = zm;
 66:     gxs = xs;
 67:     gys = ys;
 68:     gzs = zs;
 69:   } 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);

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

 81: /*MC
 82:   DMDAVecRestoreArrayF90 - check Fortran Notes at `DMDAVecRestoreArray()`

 84:   Level: intermediate
 85: M*/

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

 90:    Logically collective on da

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

 98:   Level: intermediate

100:   Fortran Note:
101:   From Fortran use `DMDAVecRestoreArayF90()`

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

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

119:   /* Handle case where user passes in global vector as opposed to local */
120:   PetscCall(VecGetLocalSize(vec, &N));
121:   if (N == xm * ym * zm * dof) {
122:     gxm = xm;
123:     gym = ym;
124:     gzm = zm;
125:     gxs = xs;
126:     gys = ys;
127:     gzs = zs;
128:   } 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);

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

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

143:   Level: intermediate
144: M*/

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

150:    Logically collective on Vec

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

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

159:   Level: intermediate

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

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

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

170:    Fortran Notes:
171:    From Fortran use `DMDAVecGetArrayWriteF90()` and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
172:    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
173:    dimension of the `DMDA`.

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

179:   Developer Note:
180:   This has code duplication with `DMDAVecGetArray()` and `DMDAVecGetArrayRead()`

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

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

201:   /* Handle case where user passes in global vector as opposed to local */
202:   PetscCall(VecGetLocalSize(vec, &N));
203:   if (N == xm * ym * zm * dof) {
204:     gxm = xm;
205:     gym = ym;
206:     gzm = zm;
207:     gxs = xs;
208:     gys = ys;
209:     gzs = zs;
210:   } 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);

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

222: /*MC
223:   DMDAVecRestoreArrayWriteF90 - check Fortran Notes at `DMDAVecRestoreArrayWrite()`

225:   Level: intermediate
226: M*/

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

231:    Logically collective on vec

233:    Input Parameters:
234: +  da - the distributed array
235: .  vec - the vector, either a vector the same size as one obtained with
236:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`
237: -  array - the array, non-NULL pointer is zeroed

239:   Level: intermediate

241:   Fortran Note:
242:   Use `DMDAVecRestoreArayWriteF90()`

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

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

263:   /* Handle case where user passes in global vector as opposed to local */
264:   PetscCall(VecGetLocalSize(vec, &N));
265:   if (N == xm * ym * zm * dof) {
266:     gxm = xm;
267:     gym = ym;
268:     gzm = zm;
269:     gxs = xs;
270:     gys = ys;
271:     gzs = zs;
272:   } 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);

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

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

288:    Logically collective

290:    Input Parameters:
291: +  da - the distributed array
292: -  vec - the vector, either a vector the same size as one obtained with
293:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`

295:    Output Parameter:
296: .  array - the array

298:   Level: intermediate

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

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

305:    Fortran Notes:
306:     Use `DMDAVecGetArrayF90()` and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
307:    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
308:    dimension of the `DMDA`.

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

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

321:   PetscFunctionBegin;
322:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
323:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
324:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

326:   /* Handle case where user passes in global vector as opposed to local */
327:   PetscCall(VecGetLocalSize(vec, &N));
328:   if (N == xm * ym * zm * dof) {
329:     gxm = xm;
330:     gym = ym;
331:     gzm = zm;
332:     gxs = xs;
333:     gys = ys;
334:     gzs = zs;
335:   } 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);

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

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

350:    Logically collective

352:    Input Parameters:
353: +  da - the distributed array
354: .  vec - the vector, either a vector the same size as one obtained with
355:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`
356: -  array - the array

358:   Level: intermediate

360:   Fortran Note:
361:   Use `DMDAVecRestoreArayF90()`

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

370:   PetscFunctionBegin;
371:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
372:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
373:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

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

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

396: /*MC
397:   DMDAVecGetArrayReadF90 - check Fortran Notes at `DMDAVecGetArrayRead()`

399:   Level: intermediate
400: M*/

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

406:    Not collective

408:    Input Parameters:
409: +  da - the distributed array
410: -  vec - the vector, either a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`

412:    Output Parameter:
413: .  array - the array

415:   Level: intermediate

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

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

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

426:   Fortran Notes:
427:   Use `DMDAVecGetArrayReadF90()` and pass for the array type `PetscScalar`,pointer :: array(:,...,:) of the appropriate
428:   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
429:   dimension of the `DMDA`.

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

435: .seealso: `DM`, `DMDA`, `DMDAVecGetArrayReadF90()`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArrayRead()`, `DMDAVecRestoreArrayDOF()`
436:           `DMDAVecGetArrayDOF()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`,
437:           `DMStagVecGetArrayRead()`
438: @*/
439: PetscErrorCode DMDAVecGetArrayRead(DM da, Vec vec, void *array)
440: {
441:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

443:   PetscFunctionBegin;
447:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
448:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
449:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

451:   /* Handle case where user passes in global vector as opposed to local */
452:   PetscCall(VecGetLocalSize(vec, &N));
453:   if (N == xm * ym * zm * dof) {
454:     gxm = xm;
455:     gym = ym;
456:     gzm = zm;
457:     gxs = xs;
458:     gys = ys;
459:     gzs = zs;
460:   } 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);

462:   if (dim == 1) {
463:     PetscCall(VecGetArray1dRead(vec, gxm * dof, gxs * dof, (PetscScalar **)array));
464:   } else if (dim == 2) {
465:     PetscCall(VecGetArray2dRead(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array));
466:   } else if (dim == 3) {
467:     PetscCall(VecGetArray3dRead(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array));
468:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
469:   PetscFunctionReturn(PETSC_SUCCESS);
470: }

472: /*MC
473:   DMDAVecRestoreArrayReadF90 - check Fortran Notes at `DMDAVecRestoreArrayRead()`

475:   Level: intermediate
476: M*/

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

481:    Not collective

483:    Input Parameters:
484: +  da - the distributed array
485: .  vec - the vector, either a vector the same size as one obtained with
486:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`
487: -  array - the array, non-NULL pointer is zeroed

489:   Level: intermediate

491:   Fortran Note:
492:   Use `DMDAVecRestoreArrayReadF90()`

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

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

510:   /* Handle case where user passes in global vector as opposed to local */
511:   PetscCall(VecGetLocalSize(vec, &N));
512:   if (N == xm * ym * zm * dof) {
513:     gxm = xm;
514:     gym = ym;
515:     gzm = zm;
516:     gxs = xs;
517:     gys = ys;
518:     gzs = zs;
519:   } 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);

521:   if (dim == 1) {
522:     PetscCall(VecRestoreArray1dRead(vec, gxm * dof, gxs * dof, (PetscScalar **)array));
523:   } else if (dim == 2) {
524:     PetscCall(VecRestoreArray2dRead(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array));
525:   } else if (dim == 3) {
526:     PetscCall(VecRestoreArray3dRead(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array));
527:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
528:   PetscFunctionReturn(PETSC_SUCCESS);
529: }

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

535:    Not Collective

537:    Input Parameters:
538: +  da - the distributed array
539: -  vec - the vector, either a vector the same size as one obtained with
540:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`

542:    Output Parameter:
543: .  array - the array

545:   Level: intermediate

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

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

552:    Fortran Note:
553:    Use  `DMDAVecGetArrayReadF90()` and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
554:    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
555:    dimension of the `DMDA`.

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

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

568:   PetscFunctionBegin;
569:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
570:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
571:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

573:   /* Handle case where user passes in global vector as opposed to local */
574:   PetscCall(VecGetLocalSize(vec, &N));
575:   if (N == xm * ym * zm * dof) {
576:     gxm = xm;
577:     gym = ym;
578:     gzm = zm;
579:     gxs = xs;
580:     gys = ys;
581:     gzs = zs;
582:   } 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);

584:   if (dim == 1) {
585:     PetscCall(VecGetArray2dRead(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
586:   } else if (dim == 2) {
587:     PetscCall(VecGetArray3dRead(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
588:   } else if (dim == 3) {
589:     PetscCall(VecGetArray4dRead(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
590:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
591:   PetscFunctionReturn(PETSC_SUCCESS);
592: }

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

597:    Not Collective

599:    Input Parameters:
600: +  da - the distributed array
601: .  vec - the vector, either a vector the same size as one obtained with
602:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`
603: -  array - the array

605:   Level: intermediate

607:   Fortran Note:
608:   Use `DMDAVecRestoreArrayReadF90()`

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

617:   PetscFunctionBegin;
618:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
619:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
620:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

622:   /* Handle case where user passes in global vector as opposed to local */
623:   PetscCall(VecGetLocalSize(vec, &N));
624:   if (N == xm * ym * zm * dof) {
625:     gxm = xm;
626:     gym = ym;
627:     gzm = zm;
628:     gxs = xs;
629:     gys = ys;
630:     gzs = zs;
631:   }

633:   if (dim == 1) {
634:     PetscCall(VecRestoreArray2dRead(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
635:   } else if (dim == 2) {
636:     PetscCall(VecRestoreArray3dRead(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
637:   } else if (dim == 3) {
638:     PetscCall(VecRestoreArray4dRead(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
639:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
640:   PetscFunctionReturn(PETSC_SUCCESS);
641: }

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

647:    Not Collective

649:    Input Parameters:
650: +  da - the distributed array
651: -  vec - the vector, either a vector the same size as one obtained with
652:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`

654:    Output Parameter:
655: .  array - the array

657:    Level: intermediate

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

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

664:    Fortran Note:
665:    Use  `DMDAVecGetArrayWriteF90()` and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
666:    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
667:    dimension of the `DMDA`.

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

673: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `DMDAVecRestoreArrayWriteF90()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`,
674:           `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`
675: @*/
676: PetscErrorCode DMDAVecGetArrayDOFWrite(DM da, Vec vec, void *array)
677: {
678:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

680:   PetscFunctionBegin;
681:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
682:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
683:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

685:   /* Handle case where user passes in global vector as opposed to local */
686:   PetscCall(VecGetLocalSize(vec, &N));
687:   if (N == xm * ym * zm * dof) {
688:     gxm = xm;
689:     gym = ym;
690:     gzm = zm;
691:     gxs = xs;
692:     gys = ys;
693:     gzs = zs;
694:   } 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);

696:   if (dim == 1) {
697:     PetscCall(VecGetArray2dWrite(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
698:   } else if (dim == 2) {
699:     PetscCall(VecGetArray3dWrite(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
700:   } else if (dim == 3) {
701:     PetscCall(VecGetArray4dWrite(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
702:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
703:   PetscFunctionReturn(PETSC_SUCCESS);
704: }

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

709:    Not Collective

711:    Input Parameters:
712: +  da - the distributed array
713: .  vec - the vector, either a vector the same size as one obtained with
714:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`
715: -  array - the array

717:   Level: intermediate

719:   Fortran Note:
720:   Use `DMDAVecRestoreArrayWriteF90()`

722: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`, `DMDAVecRestoreArrayDOF()`,
723:           `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`
724: @*/
725: PetscErrorCode DMDAVecRestoreArrayDOFWrite(DM da, Vec vec, void *array)
726: {
727:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

729:   PetscFunctionBegin;
730:   PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
731:   PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
732:   PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));

734:   /* Handle case where user passes in global vector as opposed to local */
735:   PetscCall(VecGetLocalSize(vec, &N));
736:   if (N == xm * ym * zm * dof) {
737:     gxm = xm;
738:     gym = ym;
739:     gzm = zm;
740:     gxs = xs;
741:     gys = ys;
742:     gzs = zs;
743:   }

745:   if (dim == 1) {
746:     PetscCall(VecRestoreArray2dWrite(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
747:   } else if (dim == 2) {
748:     PetscCall(VecRestoreArray3dWrite(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
749:   } else if (dim == 3) {
750:     PetscCall(VecRestoreArray4dWrite(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
751:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
752:   PetscFunctionReturn(PETSC_SUCCESS);
753: }