Actual source code: dagetarray.c


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

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

  8:    Logically collective on da

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

 14:    Output Parameter:
 15: .  array - the array

 17:   Level: intermediate

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

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

 24:     If vec is a local vector (obtained with DMCreateLocalVector() etc) then the ghost point locations are accessible. If it is
 25:     a global vector then the ghost points are not accessible. Of course with the local vector you will have had to do the

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

 29:   Fortran Notes:
 30:     From Fortran use `DMDAVecGetArrayF90()` and pass for the array type `PetscScalar`,pointer :: array(:,...,:) of the appropriate
 31:        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
 32:        dimension of the `DMDA`. The order of the indices is array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) (when dof is 1) otherwise
 33:        array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from
 34:        DMDAGetCorners() for a global array or DMDAGetGhostCorners() for a local array. Include petsc/finclude/petscdmda.h90 to access this routine.

 36:   Due to bugs in the compiler `DMDAVecGetArrayF90()` does not work with gfortran versions before 4.5

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

 49:   DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm);
 50:   DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
 51:   DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);

 53:   /* Handle case where user passes in global vector as opposed to local */
 54:   VecGetLocalSize(vec, &N);
 55:   if (N == xm * ym * zm * dof) {
 56:     gxm = xm;
 57:     gym = ym;
 58:     gzm = zm;
 59:     gxs = xs;
 60:     gys = ys;
 61:     gzs = zs;

 64:   if (dim == 1) {
 65:     VecGetArray1d(vec, gxm * dof, gxs * dof, (PetscScalar **)array);
 66:   } else if (dim == 2) {
 67:     VecGetArray2d(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array);
 68:   } else if (dim == 3) {
 69:     VecGetArray3d(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array);
 70:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
 71:   return 0;
 72: }

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

 77:    Logically collective on da

 79:    Input Parameters:
 80: +  da - the distributed array
 81: .  vec - the vector, either a vector the same size as one obtained with
 82:          DMCreateGlobalVector() or DMCreateLocalVector()
 83: -  array - the array, non-NULL pointer is zeroed

 85:   Level: intermediate

 87:   Fortran Note:
 88:     From Fortran use `DMDAVecRestoreArayF90()`

 90: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArray()`,
 91:           `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`,
 92:           `DMDStagVecRestoreArray()`
 93: @*/
 94: PetscErrorCode DMDAVecRestoreArray(DM da, Vec vec, void *array)
 95: {
 96:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

101:   DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm);
102:   DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
103:   DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);

105:   /* Handle case where user passes in global vector as opposed to local */
106:   VecGetLocalSize(vec, &N);
107:   if (N == xm * ym * zm * dof) {
108:     gxm = xm;
109:     gym = ym;
110:     gzm = zm;
111:     gxs = xs;
112:     gys = ys;
113:     gzs = zs;

116:   if (dim == 1) {
117:     VecRestoreArray1d(vec, gxm * dof, gxs * dof, (PetscScalar **)array);
118:   } else if (dim == 2) {
119:     VecRestoreArray2d(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array);
120:   } else if (dim == 3) {
121:     VecRestoreArray3d(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array);
122:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
123:   return 0;
124: }

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

130:    Logically collective on Vec

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

136:    Output Parameter:
137: .  array - the array

139:   Level: intermediate

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

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

146:     If vec is a local vector (obtained with `DMCreateLocalVector()` etc) then the ghost point locations are accessible. If it is
147:     a global vector then the ghost points are not accessible. Of course with the local vector you will have had to do the

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

151:   Fortran Notes:
152:     From Fortran use `DMDAVecGetArrayF90()` and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
153:        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
154:        dimension of the `DMDA`. The order of the indices is array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) (when dof is 1) otherwise
155:        array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from
156:        DMDAGetCorners() for a global array or `DMDAGetGhostCorners()` for a local array. Include petsc/finclude/petscdmda.h90 to access this routine.

158:   Due to bugs in the compiler `DMDAVecGetArrayF90()` does not work with gfortran versions before 4.5

160:   Developer Note:
161:   This has code duplication with `DMDAVecGetArray()` and `DMDAVecGetArrayRead()`

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

173:   if (da->localSection) {
174:     VecGetArrayWrite(vec, (PetscScalar **)array);
175:     return 0;
176:   }
177:   DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm);
178:   DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
179:   DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);

181:   /* Handle case where user passes in global vector as opposed to local */
182:   VecGetLocalSize(vec, &N);
183:   if (N == xm * ym * zm * dof) {
184:     gxm = xm;
185:     gym = ym;
186:     gzm = zm;
187:     gxs = xs;
188:     gys = ys;
189:     gzs = zs;

192:   if (dim == 1) {
193:     VecGetArray1dWrite(vec, gxm * dof, gxs * dof, (PetscScalar **)array);
194:   } else if (dim == 2) {
195:     VecGetArray2dWrite(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array);
196:   } else if (dim == 3) {
197:     VecGetArray3dWrite(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array);
198:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
199:   return 0;
200: }

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

205:    Logically collective on vec

207:    Input Parameters:
208: +  da - the distributed array
209: .  vec - the vector, either a vector the same size as one obtained with
210:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`
211: -  array - the array, non-NULL pointer is zeroed

213:   Level: intermediate

215:   Fortran Note:
216:     From Fortran use `DMDAVecRestoreArayF90()`

218: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArrayWrite()`,
219:           `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`
220: @*/
221: PetscErrorCode DMDAVecRestoreArrayWrite(DM da, Vec vec, void *array)
222: {
223:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

228:   if (da->localSection) {
229:     VecRestoreArray(vec, (PetscScalar **)array);
230:     return 0;
231:   }
232:   DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm);
233:   DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
234:   DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);

236:   /* Handle case where user passes in global vector as opposed to local */
237:   VecGetLocalSize(vec, &N);
238:   if (N == xm * ym * zm * dof) {
239:     gxm = xm;
240:     gym = ym;
241:     gzm = zm;
242:     gxs = xs;
243:     gys = ys;
244:     gzs = zs;

247:   if (dim == 1) {
248:     VecRestoreArray1dWrite(vec, gxm * dof, gxs * dof, (PetscScalar **)array);
249:   } else if (dim == 2) {
250:     VecRestoreArray2dWrite(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array);
251:   } else if (dim == 3) {
252:     VecRestoreArray3dWrite(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array);
253:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
254:   return 0;
255: }

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

261:    Logically collective

263:    Input Parameters:
264: +  da - the distributed array
265: -  vec - the vector, either a vector the same size as one obtained with
266:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`

268:    Output Parameter:
269: .  array - the array

271:   Level: intermediate

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

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

278:     In Fortran 90 you do not need a version of `DMDAVecRestoreArrayDOF()` just use `DMDAVecRestoreArrayF90()` and declare your array with one higher dimension,
279:     see src/dm/tutorials/ex11f90.F

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

288:   DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm);
289:   DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
290:   DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);

292:   /* Handle case where user passes in global vector as opposed to local */
293:   VecGetLocalSize(vec, &N);
294:   if (N == xm * ym * zm * dof) {
295:     gxm = xm;
296:     gym = ym;
297:     gzm = zm;
298:     gxs = xs;
299:     gys = ys;
300:     gzs = zs;

303:   if (dim == 1) {
304:     VecGetArray2d(vec, gxm, dof, gxs, 0, (PetscScalar ***)array);
305:   } else if (dim == 2) {
306:     VecGetArray3d(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array);
307:   } else if (dim == 3) {
308:     VecGetArray4d(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array);
309:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
310:   return 0;
311: }

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

316:    Logically collective

318:    Input Parameters:
319: +  da - the distributed array
320: .  vec - the vector, either a vector the same size as one obtained with
321:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`
322: -  array - the array

324:   Level: intermediate

326: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`, `DMDAVecRestoreArrayDOF()`,
327:           `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`
328: @*/
329: PetscErrorCode DMDAVecRestoreArrayDOF(DM da, Vec vec, void *array)
330: {
331:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

333:   DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm);
334:   DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
335:   DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);

337:   /* Handle case where user passes in global vector as opposed to local */
338:   VecGetLocalSize(vec, &N);
339:   if (N == xm * ym * zm * dof) {
340:     gxm = xm;
341:     gym = ym;
342:     gzm = zm;
343:     gxs = xs;
344:     gys = ys;
345:     gzs = zs;
346:   }

348:   if (dim == 1) {
349:     VecRestoreArray2d(vec, gxm, dof, gxs, 0, (PetscScalar ***)array);
350:   } else if (dim == 2) {
351:     VecRestoreArray3d(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array);
352:   } else if (dim == 3) {
353:     VecRestoreArray4d(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array);
354:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
355:   return 0;
356: }

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

362:    Not collective

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

368:    Output Parameter:
369: .  array - the array

371:   Level: intermediate

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

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

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

382:   Fortran Notes:
383:     From Fortran use `DMDAVecGetArrayReadF90()` and pass for the array type `PetscScalar`,pointer :: array(:,...,:) of the appropriate
384:        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
385:        dimension of the `DMDA`. The order of the indices is array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) (when dof is 1) otherwise
386:        array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from
387:        DMDAGetCorners() for a global array or `DMDAGetGhostCorners()` for a local array. Include petsc/finclude/petscdmda.h90 to access this routine.

389:   Due to bugs in the compiler `DMDAVecGetArrayReadF90()` does not work with gfortran versions before 4.5

391: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArrayRead()`, `DMDAVecRestoreArrayDOF()`
392:           `DMDAVecGetArrayDOF()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`,
393:           `DMStagVecGetArrayRead()`
394: @*/
395: PetscErrorCode DMDAVecGetArrayRead(DM da, Vec vec, void *array)
396: {
397:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

402:   DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm);
403:   DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
404:   DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);

406:   /* Handle case where user passes in global vector as opposed to local */
407:   VecGetLocalSize(vec, &N);
408:   if (N == xm * ym * zm * dof) {
409:     gxm = xm;
410:     gym = ym;
411:     gzm = zm;
412:     gxs = xs;
413:     gys = ys;
414:     gzs = zs;

417:   if (dim == 1) {
418:     VecGetArray1dRead(vec, gxm * dof, gxs * dof, (PetscScalar **)array);
419:   } else if (dim == 2) {
420:     VecGetArray2dRead(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array);
421:   } else if (dim == 3) {
422:     VecGetArray3dRead(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array);
423:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
424:   return 0;
425: }

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

430:    Not collective

432:    Input Parameters:
433: +  da - the distributed array
434: .  vec - the vector, either a vector the same size as one obtained with
435:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`
436: -  array - the array, non-NULL pointer is zeroed

438:   Level: intermediate

440:   Fortran Note:
441:     From Fortran use `DMDAVecRestoreArrayReadF90()`

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

454:   DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm);
455:   DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
456:   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:   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;

469:   if (dim == 1) {
470:     VecRestoreArray1dRead(vec, gxm * dof, gxs * dof, (PetscScalar **)array);
471:   } else if (dim == 2) {
472:     VecRestoreArray2dRead(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array);
473:   } else if (dim == 3) {
474:     VecRestoreArray3dRead(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:   return 0;
477: }

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

483:    Not Collective

485:    Input Parameters:
486: +  da - the distributed array
487: -  vec - the vector, either a vector the same size as one obtained with
488:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`

490:    Output Parameter:
491: .  array - the array

493:   Level: intermediate

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

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

500:    Fortran Note:
501:     In Fortran you do not need a version of `DMDAVecRestoreArrayDOF()` just use  `DMDAVecRestoreArrayReadF90()` and declare your array with one higher dimension,
502:     see src/dm/tutorials/ex11f90.F

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

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

515:   /* Handle case where user passes in global vector as opposed to local */
516:   VecGetLocalSize(vec, &N);
517:   if (N == xm * ym * zm * dof) {
518:     gxm = xm;
519:     gym = ym;
520:     gzm = zm;
521:     gxs = xs;
522:     gys = ys;
523:     gzs = zs;

526:   if (dim == 1) {
527:     VecGetArray2dRead(vec, gxm, dof, gxs, 0, (PetscScalar ***)array);
528:   } else if (dim == 2) {
529:     VecGetArray3dRead(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array);
530:   } else if (dim == 3) {
531:     VecGetArray4dRead(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array);
532:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
533:   return 0;
534: }

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

539:    Not Collective

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

547:   Level: intermediate

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

556:   DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm);
557:   DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
558:   DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);

560:   /* Handle case where user passes in global vector as opposed to local */
561:   VecGetLocalSize(vec, &N);
562:   if (N == xm * ym * zm * dof) {
563:     gxm = xm;
564:     gym = ym;
565:     gzm = zm;
566:     gxs = xs;
567:     gys = ys;
568:     gzs = zs;
569:   }

571:   if (dim == 1) {
572:     VecRestoreArray2dRead(vec, gxm, dof, gxs, 0, (PetscScalar ***)array);
573:   } else if (dim == 2) {
574:     VecRestoreArray3dRead(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array);
575:   } else if (dim == 3) {
576:     VecRestoreArray4dRead(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array);
577:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
578:   return 0;
579: }

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

585:    Not Collective

587:    Input Parameters:
588: +  da - the distributed array
589: -  vec - the vector, either a vector the same size as one obtained with
590:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`

592:    Output Parameter:
593: .  array - the array

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

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

600:    Fortran Noe:
601:     In Fortran you do not need a version of `DMDAVecRestoreArrayDOF()` just use  `DMDAVecRestoreArrayWriteF90()` and declare your array with one higher dimension,
602:     see src/dm/tutorials/ex11f90.F

604:   Level: intermediate

606: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`,
607:           `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`
608: @*/
609: PetscErrorCode DMDAVecGetArrayDOFWrite(DM da, Vec vec, void *array)
610: {
611:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

613:   DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm);
614:   DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
615:   DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);

617:   /* Handle case where user passes in global vector as opposed to local */
618:   VecGetLocalSize(vec, &N);
619:   if (N == xm * ym * zm * dof) {
620:     gxm = xm;
621:     gym = ym;
622:     gzm = zm;
623:     gxs = xs;
624:     gys = ys;
625:     gzs = zs;

628:   if (dim == 1) {
629:     VecGetArray2dWrite(vec, gxm, dof, gxs, 0, (PetscScalar ***)array);
630:   } else if (dim == 2) {
631:     VecGetArray3dWrite(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array);
632:   } else if (dim == 3) {
633:     VecGetArray4dWrite(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array);
634:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
635:   return 0;
636: }

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

641:    Not Collective

643:    Input Parameters:
644: +  da - the distributed array
645: .  vec - the vector, either a vector the same size as one obtained with
646:          `DMCreateGlobalVector()` or `DMCreateLocalVector()`
647: -  array - the array

649:   Level: intermediate

651: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`, `DMDAVecRestoreArrayDOF()`,
652:           `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`
653: @*/
654: PetscErrorCode DMDAVecRestoreArrayDOFWrite(DM da, Vec vec, void *array)
655: {
656:   PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;

658:   DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm);
659:   DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm);
660:   DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL);

662:   /* Handle case where user passes in global vector as opposed to local */
663:   VecGetLocalSize(vec, &N);
664:   if (N == xm * ym * zm * dof) {
665:     gxm = xm;
666:     gym = ym;
667:     gzm = zm;
668:     gxs = xs;
669:     gys = ys;
670:     gzs = zs;
671:   }

673:   if (dim == 1) {
674:     VecRestoreArray2dWrite(vec, gxm, dof, gxs, 0, (PetscScalar ***)array);
675:   } else if (dim == 2) {
676:     VecRestoreArray3dWrite(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array);
677:   } else if (dim == 3) {
678:     VecRestoreArray4dWrite(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array);
679:   } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
680:   return 0;
681: }