Actual source code: dagetarray.c
1: #include <petsc/private/dmdaimpl.h>
3: /*MC
4: DMDAVecGetArrayF90 - check Fortran Notes at `DMDAVecGetArray()`
6: Level: intermediate
7: M*/
9: /*@C
10: DMDAVecGetArray - Returns a multiple dimension array that shares data with
11: the underlying vector and is indexed using the global dimensions.
13: Logically Collective
15: Input Parameters:
16: + da - the distributed array
17: - vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
19: Output Parameter:
20: . array - the array
22: Level: intermediate
24: Notes:
25: Call `DMDAVecRestoreArray()` once you have finished accessing the vector entries.
27: In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]!
29: If `vec` is a local vector (obtained with DMCreateLocalVector() etc) then the ghost point locations are accessible. If it is
30: a global vector then the ghost points are not accessible. Of course, with a local vector you will have had to do the
31: appropriate `DMGlobalToLocalBegin()` and `DMGlobalToLocalEnd()` to have correct values in the ghost locations.
33: The accessible indices are `array[zs:zs+zm-1][ys:ys+ym-1][xs:xs+xm-1]` where the values are obtained from
34: `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.
36: Fortran Notes:
37: Use `DMDAVecGetArrayF90()` and pass for the array type `PetscScalar`,pointer :: array(:,...,:) of the appropriate
38: 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
39: dimension of the `DMDA`.
41: The order of the indices is `array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` (when dof is 1) otherwise
42: `array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` where the values are obtained from
43: `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.
45: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArray()`, `DMDAVecRestoreArrayDOF()`
46: `DMDAVecGetArrayDOF()`, `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`,
47: `DMStagVecGetArray()`
48: @*/
49: PetscErrorCode DMDAVecGetArray(DM da, Vec vec, void *array)
50: {
51: PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
53: PetscFunctionBegin;
56: PetscAssertPointer(array, 3);
57: PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
58: PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
59: PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));
61: /* Handle case where user passes in global vector as opposed to local */
62: PetscCall(VecGetLocalSize(vec, &N));
63: if (N == xm * ym * zm * dof) {
64: gxm = xm;
65: gym = ym;
66: gzm = zm;
67: gxs = xs;
68: gys = ys;
69: gzs = zs;
70: } 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);
72: if (dim == 1) {
73: PetscCall(VecGetArray1d(vec, gxm * dof, gxs * dof, (PetscScalar **)array));
74: } else if (dim == 2) {
75: PetscCall(VecGetArray2d(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array));
76: } else if (dim == 3) {
77: PetscCall(VecGetArray3d(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array));
78: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
79: PetscFunctionReturn(PETSC_SUCCESS);
80: }
82: /*MC
83: DMDAVecRestoreArrayF90 - check Fortran Notes at `DMDAVecRestoreArray()`
85: Level: intermediate
86: M*/
88: /*@
89: DMDAVecRestoreArray - Restores a multiple dimension array obtained with `DMDAVecGetArray()`
91: Logically Collective
93: Input Parameters:
94: + da - the distributed array
95: . vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
96: - array - the `array` pointer is zeroed
98: Level: intermediate
100: Fortran Notes:
101: 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;
114: PetscAssertPointer(array, 3);
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
152: Input Parameters:
153: + da - the distributed array
154: - vec - 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: The accessible indices are `array[zs:zs+zm-1][ys:ys+ym-1][xs:xs+xm-1]` where the values are obtained from
171: `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.
173: 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: The order of the indices is array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) (when dof is 1) otherwise
184: array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from
185: `DMDAGetCorners()` for a global array or `DMDAGetGhostCorners()` for a local array.
187: Developer Notes:
188: This has code duplication with `DMDAVecGetArray()` and `DMDAVecGetArrayRead()`
190: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecRestoreArrayDOF()`
191: `DMDAVecGetArrayDOF()`, `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`
192: @*/
193: PetscErrorCode DMDAVecGetArrayWrite(DM da, Vec vec, void *array)
194: {
195: PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
197: PetscFunctionBegin;
200: PetscAssertPointer(array, 3);
201: if (da->localSection) {
202: PetscCall(VecGetArrayWrite(vec, (PetscScalar **)array));
203: PetscFunctionReturn(PETSC_SUCCESS);
204: }
205: PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
206: PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
207: PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));
209: /* Handle case where user passes in global vector as opposed to local */
210: PetscCall(VecGetLocalSize(vec, &N));
211: if (N == xm * ym * zm * dof) {
212: gxm = xm;
213: gym = ym;
214: gzm = zm;
215: gxs = xs;
216: gys = ys;
217: gzs = zs;
218: } 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);
220: if (dim == 1) {
221: PetscCall(VecGetArray1dWrite(vec, gxm * dof, gxs * dof, (PetscScalar **)array));
222: } else if (dim == 2) {
223: PetscCall(VecGetArray2dWrite(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array));
224: } else if (dim == 3) {
225: PetscCall(VecGetArray3dWrite(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array));
226: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
227: PetscFunctionReturn(PETSC_SUCCESS);
228: }
230: /*MC
231: DMDAVecRestoreArrayWriteF90 - check Fortran Notes at `DMDAVecRestoreArrayWrite()`
233: Level: intermediate
234: M*/
236: /*@
237: DMDAVecRestoreArrayWrite - Restores a multiple dimension array obtained with `DMDAVecGetArrayWrite()`
239: Logically Collective
241: Input Parameters:
242: + da - the distributed array
243: . vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
244: - array - the `array` pointer is zeroed
246: Level: intermediate
248: Fortran Notes:
249: Use `DMDAVecRestoreArayWriteF90()`
251: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArrayWrite()`,
252: `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`
253: @*/
254: PetscErrorCode DMDAVecRestoreArrayWrite(DM da, Vec vec, void *array)
255: {
256: PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
258: PetscFunctionBegin;
261: PetscAssertPointer(array, 3);
262: if (da->localSection) {
263: PetscCall(VecRestoreArray(vec, (PetscScalar **)array));
264: PetscFunctionReturn(PETSC_SUCCESS);
265: }
266: PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
267: PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
268: PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));
270: /* Handle case where user passes in global vector as opposed to local */
271: PetscCall(VecGetLocalSize(vec, &N));
272: if (N == xm * ym * zm * dof) {
273: gxm = xm;
274: gym = ym;
275: gzm = zm;
276: gxs = xs;
277: gys = ys;
278: gzs = zs;
279: } 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);
281: if (dim == 1) {
282: PetscCall(VecRestoreArray1dWrite(vec, gxm * dof, gxs * dof, (PetscScalar **)array));
283: } else if (dim == 2) {
284: PetscCall(VecRestoreArray2dWrite(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array));
285: } else if (dim == 3) {
286: PetscCall(VecRestoreArray3dWrite(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array));
287: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
288: PetscFunctionReturn(PETSC_SUCCESS);
289: }
291: /*@C
292: DMDAVecGetArrayDOF - Returns a multiple dimension array that shares data with
293: the underlying vector and is indexed using the global dimensions.
295: Logically Collective
297: Input Parameters:
298: + da - the distributed array
299: - vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
301: Output Parameter:
302: . array - the `array` pointer is zeroed
304: Level: intermediate
306: Notes:
307: Call `DMDAVecRestoreArrayDOF()` once you have finished accessing the vector entries.
309: In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]
311: 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
312: `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.
314: Fortran Notes:
315: Use `DMDAVecGetArrayF90()` and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
316: 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
317: dimension of the `DMDA`.
319: The order of the indices is `array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` (when ndof is 1) otherwise
320: `array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` where the values are obtained from
321: `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.
323: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecRestoreArrayDOF()`,
324: `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`, `DMDAVecGetArrayDOFRead()`
325: @*/
326: PetscErrorCode DMDAVecGetArrayDOF(DM da, Vec vec, void *array)
327: {
328: PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
330: PetscFunctionBegin;
331: PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
332: PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
333: PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));
335: /* Handle case where user passes in global vector as opposed to local */
336: PetscCall(VecGetLocalSize(vec, &N));
337: if (N == xm * ym * zm * dof) {
338: gxm = xm;
339: gym = ym;
340: gzm = zm;
341: gxs = xs;
342: gys = ys;
343: gzs = zs;
344: } 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);
346: if (dim == 1) {
347: PetscCall(VecGetArray2d(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
348: } else if (dim == 2) {
349: PetscCall(VecGetArray3d(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
350: } else if (dim == 3) {
351: PetscCall(VecGetArray4d(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
352: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
353: PetscFunctionReturn(PETSC_SUCCESS);
354: }
356: /*@
357: DMDAVecRestoreArrayDOF - Restores a multiple dimension array obtained with `DMDAVecGetArrayDOF()`
359: Logically Collective
361: Input Parameters:
362: + da - the distributed array
363: . vec - vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
364: - array - the `array` point is zeroed
366: Level: intermediate
368: Fortran Notes:
369: Use `DMDAVecRestoreArayF90()`
371: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`,
372: `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`
373: @*/
374: PetscErrorCode DMDAVecRestoreArrayDOF(DM da, Vec vec, void *array)
375: {
376: PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
378: PetscFunctionBegin;
379: PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
380: PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
381: PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));
383: /* Handle case where user passes in global vector as opposed to local */
384: PetscCall(VecGetLocalSize(vec, &N));
385: if (N == xm * ym * zm * dof) {
386: gxm = xm;
387: gym = ym;
388: gzm = zm;
389: gxs = xs;
390: gys = ys;
391: gzs = zs;
392: }
394: if (dim == 1) {
395: PetscCall(VecRestoreArray2d(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
396: } else if (dim == 2) {
397: PetscCall(VecRestoreArray3d(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
398: } else if (dim == 3) {
399: PetscCall(VecRestoreArray4d(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
400: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
401: PetscFunctionReturn(PETSC_SUCCESS);
402: }
404: /*MC
405: DMDAVecGetArrayReadF90 - check Fortran Notes at `DMDAVecGetArrayRead()`
407: Level: intermediate
408: M*/
410: /*@C
411: DMDAVecGetArrayRead - Returns a multiple dimension array that shares data with
412: the underlying vector and is indexed using the global dimensions.
414: Not Collective
416: Input Parameters:
417: + da - the distributed array
418: - vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
420: Output Parameter:
421: . array - the array
423: Level: intermediate
425: Notes:
426: Call `DMDAVecRestoreArrayRead()` once you have finished accessing the vector entries.
428: In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]!
430: If `vec` is a local vector (obtained with `DMCreateLocalVector()` etc) then the ghost point locations are accessible. If it is
431: a global vector then the ghost points are not accessible. Of course with the local vector you will have had to do the
432: appropriate `DMGlobalToLocalBegin()` and `DMGlobalToLocalEnd()` to have correct values in the ghost locations.
434: The accessible indices are `array[zs:zs+zm-1][ys:ys+ym-1][xs:xs+xm-1]` where the values are obtained from
435: `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.
437: Fortran Notes:
438: Use `DMDAVecGetArrayReadF90()` and pass for the array type `PetscScalar`,pointer :: array(:,...,:) of the appropriate
439: 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
440: dimension of the `DMDA`.
442: The order of the indices is `array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` (when dof is 1) otherwise
443: `array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` where the values are obtained from
444: `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.
446: .seealso: `DM`, `DMDA`, `DMDAVecGetArrayReadF90()`, `DMDAGetGhostCorners()`,
447: `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArrayRead()`,
448: `DMDAVecRestoreArrayDOF()`, `DMDAVecGetArrayDOF()`, `DMDAVecGetArray()`,
449: `DMDAVecRestoreArray()`, `DMStagVecGetArrayRead()`
450: @*/
451: PetscErrorCode DMDAVecGetArrayRead(DM da, Vec vec, void *array)
452: {
453: PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
455: PetscFunctionBegin;
458: PetscAssertPointer(array, 3);
459: PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
460: PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
461: PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));
463: /* Handle case where user passes in global vector as opposed to local */
464: PetscCall(VecGetLocalSize(vec, &N));
465: if (N == xm * ym * zm * dof) {
466: gxm = xm;
467: gym = ym;
468: gzm = zm;
469: gxs = xs;
470: gys = ys;
471: gzs = zs;
472: } 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);
474: if (dim == 1) {
475: PetscCall(VecGetArray1dRead(vec, gxm * dof, gxs * dof, (PetscScalar **)array));
476: } else if (dim == 2) {
477: PetscCall(VecGetArray2dRead(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array));
478: } else if (dim == 3) {
479: PetscCall(VecGetArray3dRead(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array));
480: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
481: PetscFunctionReturn(PETSC_SUCCESS);
482: }
484: /*MC
485: DMDAVecRestoreArrayReadF90 - check Fortran Notes at `DMDAVecRestoreArrayRead()`
487: Level: intermediate
488: M*/
490: /*@
491: DMDAVecRestoreArrayRead - Restores a multiple dimension array obtained with `DMDAVecGetArrayRead()`
493: Not Collective
495: Input Parameters:
496: + da - the distributed array
497: . vec - vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
498: - array - the `array` pointer is zeroed
500: Level: intermediate
502: Fortran Notes:
503: Use `DMDAVecRestoreArrayReadF90()`
505: .seealso: `DM`, `DMDA`, `DMDAVecRestoreArrayReadF90()`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArrayRead()`,
506: `DMDAVecGetArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`,
507: `DMStagVecRestoreArrayRead()`
508: @*/
509: PetscErrorCode DMDAVecRestoreArrayRead(DM da, Vec vec, void *array)
510: {
511: PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
513: PetscFunctionBegin;
516: PetscAssertPointer(array, 3);
517: PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
518: PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
519: PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));
521: /* Handle case where user passes in global vector as opposed to local */
522: PetscCall(VecGetLocalSize(vec, &N));
523: if (N == xm * ym * zm * dof) {
524: gxm = xm;
525: gym = ym;
526: gzm = zm;
527: gxs = xs;
528: gys = ys;
529: gzs = zs;
530: } 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);
532: if (dim == 1) {
533: PetscCall(VecRestoreArray1dRead(vec, gxm * dof, gxs * dof, (PetscScalar **)array));
534: } else if (dim == 2) {
535: PetscCall(VecRestoreArray2dRead(vec, gym, gxm * dof, gys, gxs * dof, (PetscScalar ***)array));
536: } else if (dim == 3) {
537: PetscCall(VecRestoreArray3dRead(vec, gzm, gym, gxm * dof, gzs, gys, gxs * dof, (PetscScalar ****)array));
538: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
539: PetscFunctionReturn(PETSC_SUCCESS);
540: }
542: /*@C
543: DMDAVecGetArrayDOFRead - Returns a multiple dimension array that shares data with
544: the underlying vector and is indexed using the global dimensions.
546: Not Collective
548: Input Parameters:
549: + da - the distributed array
550: - vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
552: Output Parameter:
553: . array - the array
555: Level: intermediate
557: Notes:
558: Call `DMDAVecRestoreArrayDOFRead()` once you have finished accessing the vector entries.
560: In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]!
562: The accessible indices are `array[zs:zs+zm-1][ys:ys+ym-1][xs:xs+xm-1]` where the values are obtained from
563: `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.
565: Fortran Notes:
566: Use `DMDAVecGetArrayReadF90()` and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
567: 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
568: dimension of the `DMDA`.
570: The order of the indices is `array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` (when dof is 1) otherwise
571: `array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` where the values are obtained from
572: `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.
574: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`,
575: `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`
576: @*/
577: PetscErrorCode DMDAVecGetArrayDOFRead(DM da, Vec vec, void *array)
578: {
579: PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
581: PetscFunctionBegin;
582: PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
583: PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
584: PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));
586: /* Handle case where user passes in global vector as opposed to local */
587: PetscCall(VecGetLocalSize(vec, &N));
588: if (N == xm * ym * zm * dof) {
589: gxm = xm;
590: gym = ym;
591: gzm = zm;
592: gxs = xs;
593: gys = ys;
594: gzs = zs;
595: } 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);
597: if (dim == 1) {
598: PetscCall(VecGetArray2dRead(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
599: } else if (dim == 2) {
600: PetscCall(VecGetArray3dRead(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
601: } else if (dim == 3) {
602: PetscCall(VecGetArray4dRead(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
603: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
604: PetscFunctionReturn(PETSC_SUCCESS);
605: }
607: /*@
608: DMDAVecRestoreArrayDOFRead - Restores a multiple dimension array obtained with `DMDAVecGetArrayDOFRead()`
610: Not Collective
612: Input Parameters:
613: + da - the distributed array
614: . vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
615: - array - the `array` pointer is zeroed
617: Level: intermediate
619: Fortran Notes:
620: Use `DMDAVecRestoreArrayReadF90()`
622: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`, `DMDAVecRestoreArrayDOF()`,
623: `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`, `DMDAVecGetArrayRead()`, `DMDAVecRestoreArrayRead()`
624: @*/
625: PetscErrorCode DMDAVecRestoreArrayDOFRead(DM da, Vec vec, void *array)
626: {
627: PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
629: PetscFunctionBegin;
630: PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
631: PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
632: PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));
634: /* Handle case where user passes in global vector as opposed to local */
635: PetscCall(VecGetLocalSize(vec, &N));
636: if (N == xm * ym * zm * dof) {
637: gxm = xm;
638: gym = ym;
639: gzm = zm;
640: gxs = xs;
641: gys = ys;
642: gzs = zs;
643: }
645: if (dim == 1) {
646: PetscCall(VecRestoreArray2dRead(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
647: } else if (dim == 2) {
648: PetscCall(VecRestoreArray3dRead(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
649: } else if (dim == 3) {
650: PetscCall(VecRestoreArray4dRead(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
651: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
652: PetscFunctionReturn(PETSC_SUCCESS);
653: }
655: /*@C
656: DMDAVecGetArrayDOFWrite - Returns a multiple dimension array that shares data with
657: the underlying vector and is indexed using the global dimensions.
659: Not Collective
661: Input Parameters:
662: + da - the distributed array
663: - vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
665: Output Parameter:
666: . array - the array
668: Level: intermediate
670: Notes:
671: Call `DMDAVecRestoreArrayDOFWrite()` once you have finished accessing the vector entries.
673: In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]!
675: 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
676: `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.
678: Fortran Notes:
679: Use `DMDAVecGetArrayWriteF90()` and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
680: 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
681: dimension of the `DMDA`.
683: The order of the indices is `array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` (when dof is 1) otherwise
684: `array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1)` where the values are obtained from
685: `DMDAGetCorners()` for a global vector or `DMDAGetGhostCorners()` for a local vector.
687: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `DMDAVecRestoreArrayWriteF90()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`,
688: `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`
689: @*/
690: PetscErrorCode DMDAVecGetArrayDOFWrite(DM da, Vec vec, void *array)
691: {
692: PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
694: PetscFunctionBegin;
695: PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
696: PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
697: PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));
699: /* Handle case where user passes in global vector as opposed to local */
700: PetscCall(VecGetLocalSize(vec, &N));
701: if (N == xm * ym * zm * dof) {
702: gxm = xm;
703: gym = ym;
704: gzm = zm;
705: gxs = xs;
706: gys = ys;
707: gzs = zs;
708: } 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);
710: if (dim == 1) {
711: PetscCall(VecGetArray2dWrite(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
712: } else if (dim == 2) {
713: PetscCall(VecGetArray3dWrite(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
714: } else if (dim == 3) {
715: PetscCall(VecGetArray4dWrite(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
716: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
717: PetscFunctionReturn(PETSC_SUCCESS);
718: }
720: /*@
721: DMDAVecRestoreArrayDOFWrite - Restores a multiple dimension array obtained with `DMDAVecGetArrayDOFWrite()`
723: Not Collective
725: Input Parameters:
726: + da - the distributed array
727: . vec - a vector the same size as one obtained with `DMCreateGlobalVector()` or `DMCreateLocalVector()`
728: - array - the `array` pointer is zeroed
730: Level: intermediate
732: Fortran Notes:
733: Use `DMDAVecRestoreArrayWriteF90()`
735: .seealso: `DM`, `DMDA`, `DMDAGetGhostCorners()`, `DMDAGetCorners()`, `VecGetArray()`, `VecRestoreArray()`, `DMDAVecGetArray()`, `DMDAVecGetArrayDOF()`, `DMDAVecRestoreArrayDOF()`,
736: `DMDAVecGetArrayWrite()`, `DMDAVecRestoreArrayWrite()`
737: @*/
738: PetscErrorCode DMDAVecRestoreArrayDOFWrite(DM da, Vec vec, void *array)
739: {
740: PetscInt xs, ys, zs, xm, ym, zm, gxs, gys, gzs, gxm, gym, gzm, N, dim, dof;
742: PetscFunctionBegin;
743: PetscCall(DMDAGetCorners(da, &xs, &ys, &zs, &xm, &ym, &zm));
744: PetscCall(DMDAGetGhostCorners(da, &gxs, &gys, &gzs, &gxm, &gym, &gzm));
745: PetscCall(DMDAGetInfo(da, &dim, NULL, NULL, NULL, NULL, NULL, NULL, &dof, NULL, NULL, NULL, NULL, NULL));
747: /* Handle case where user passes in global vector as opposed to local */
748: PetscCall(VecGetLocalSize(vec, &N));
749: if (N == xm * ym * zm * dof) {
750: gxm = xm;
751: gym = ym;
752: gzm = zm;
753: gxs = xs;
754: gys = ys;
755: gzs = zs;
756: }
758: if (dim == 1) {
759: PetscCall(VecRestoreArray2dWrite(vec, gxm, dof, gxs, 0, (PetscScalar ***)array));
760: } else if (dim == 2) {
761: PetscCall(VecRestoreArray3dWrite(vec, gym, gxm, dof, gys, gxs, 0, (PetscScalar ****)array));
762: } else if (dim == 3) {
763: PetscCall(VecRestoreArray4dWrite(vec, gzm, gym, gxm, dof, gzs, gys, gxs, 0, (PetscScalar *****)array));
764: } else SETERRQ(PETSC_COMM_SELF, PETSC_ERR_ARG_CORRUPT, "DMDA dimension not 1, 2, or 3, it is %" PetscInt_FMT, dim);
765: PetscFunctionReturn(PETSC_SUCCESS);
766: }