Actual source code: dagetarray.c
2: #include <petsc/private/dmdaimpl.h>
4: /*MC
5: DMDAVecGetArrayF90 - check Fortran Notes at `DMDAVecGetArray()`
7: Level: intermediate
8: M*/
10: /*@C
11: DMDAVecGetArray - Returns a multiple dimension array that shares data with
12: the underlying vector and is indexed using the global dimensions.
14: Logically Collective
16: Input Parameters:
17: + da - the distributed array
18: - vec - 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
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
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
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: }