Actual source code: dagetarray.c
petsc-3.11.4 2019-09-28
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 Vec
10: Input Parameter:
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: Notes:
18: Call DMDAVecRestoreArray() once you have finished accessing the vector entries.
20: In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]!
22: If vec is a local vector (obtained with DMCreateLocalVector() etc) then the ghost point locations are accessible. If it is
23: a global vector then the ghost points are not accessible. Of course with the local vector you will have had to do the
25: appropriate DMGlobalToLocalBegin() and DMGlobalToLocalEnd() to have correct values in the ghost locations.
27: Fortran Notes:
28: From Fortran use DMDAVecGetArrayF90() and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
29: 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
30: 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
31: array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from
32: DMDAGetCorners() for a global array or DMDAGetGhostCorners() for a local array. Include petsc/finclude/petscdmda.h90 to access this routine.
34: Due to bugs in the compiler DMDAVecGetArrayF90() does not work with gfortran versions before 4.5
36: Level: intermediate
38: .keywords: distributed array, get, corners, nodes, local indices, coordinates
40: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecRestoreArrayDOF()
41: DMDAVecGetArrayDOF()
42: @*/
43: PetscErrorCode DMDAVecGetArray(DM da,Vec vec,void *array)
44: {
46: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
52: if (da->defaultSection) {
53: VecGetArray(vec,(PetscScalar**)array);
54: return(0);
55: }
56: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
57: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
58: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
60: /* Handle case where user passes in global vector as opposed to local */
61: 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 if (N != gxm*gym*gzm*dof) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Vector local size %D is not compatible with DMDA local sizes %D %D\n",N,xm*ym*zm*dof,gxm*gym*gzm*dof);
71: if (dim == 1) {
72: VecGetArray1d(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
73: } else if (dim == 2) {
74: VecGetArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
75: } else if (dim == 3) {
76: VecGetArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
77: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
78: return(0);
79: }
81: /*@
82: DMDAVecRestoreArray - Restores a multiple dimension array obtained with DMDAVecGetArray()
84: Logically collective on Vec
86: Input Parameter:
87: + da - the distributed array
88: . vec - the vector, either a vector the same size as one obtained with
89: DMCreateGlobalVector() or DMCreateLocalVector()
90: - array - the array, non-NULL pointer is zeroed
92: Level: intermediate
94: Fortran Notes:
95: From Fortran use DMDAVecRestoreArrayF90()
97: .keywords: distributed array, get, corners, nodes, local indices, coordinates
99: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray()
100: @*/
101: PetscErrorCode DMDAVecRestoreArray(DM da,Vec vec,void *array)
102: {
104: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
110: if (da->defaultSection) {
111: VecRestoreArray(vec,(PetscScalar**)array);
112: return(0);
113: }
114: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
115: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
116: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
118: /* Handle case where user passes in global vector as opposed to local */
119: VecGetLocalSize(vec,&N);
120: if (N == xm*ym*zm*dof) {
121: gxm = xm;
122: gym = ym;
123: gzm = zm;
124: gxs = xs;
125: gys = ys;
126: gzs = zs;
127: } else if (N != gxm*gym*gzm*dof) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Vector local size %D is not compatible with DMDA local sizes %D %D\n",N,xm*ym*zm*dof,gxm*gym*gzm*dof);
129: if (dim == 1) {
130: VecRestoreArray1d(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
131: } else if (dim == 2) {
132: VecRestoreArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
133: } else if (dim == 3) {
134: VecRestoreArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
135: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
136: return(0);
137: }
139: /*@C
140: DMDAVecGetArrayDOF - Returns a multiple dimension array that shares data with
141: the underlying vector and is indexed using the global dimensions.
143: Logically collective
145: Input Parameter:
146: + da - the distributed array
147: - vec - the vector, either a vector the same size as one obtained with
148: DMCreateGlobalVector() or DMCreateLocalVector()
150: Output Parameter:
151: . array - the array
153: Notes:
154: Call DMDAVecRestoreArrayDOF() once you have finished accessing the vector entries.
156: In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]!
158: In Fortran 90 you do not need a version of DMDAVecRestoreArrayDOF() just use DMDAVecRestoreArrayF90() and declare your array with one higher dimension,
159: see src/dm/examples/tutorials/ex11f90.F
161: Level: intermediate
163: .keywords: distributed array, get, corners, nodes, local indices, coordinates
165: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecGetArray(), DMDAVecRestoreArrayDOF()
166: @*/
167: PetscErrorCode DMDAVecGetArrayDOF(DM da,Vec vec,void *array)
168: {
170: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
173: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
174: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
175: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
177: /* Handle case where user passes in global vector as opposed to local */
178: VecGetLocalSize(vec,&N);
179: if (N == xm*ym*zm*dof) {
180: gxm = xm;
181: gym = ym;
182: gzm = zm;
183: gxs = xs;
184: gys = ys;
185: gzs = zs;
186: } else if (N != gxm*gym*gzm*dof) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Vector local size %D is not compatible with DMDA local sizes %D %D\n",N,xm*ym*zm*dof,gxm*gym*gzm*dof);
188: if (dim == 1) {
189: VecGetArray2d(vec,gxm,dof,gxs,0,(PetscScalar***)array);
190: } else if (dim == 2) {
191: VecGetArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
192: } else if (dim == 3) {
193: VecGetArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
194: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
195: return(0);
196: }
198: /*@
199: DMDAVecRestoreArrayDOF - Restores a multiple dimension array obtained with DMDAVecGetArrayDOF()
201: Logically collective
203: Input Parameter:
204: + da - the distributed array
205: . vec - the vector, either a vector the same size as one obtained with
206: DMCreateGlobalVector() or DMCreateLocalVector()
207: - array - the array
209: Level: intermediate
211: .keywords: distributed array, get, corners, nodes, local indices, coordinates
213: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(), DMDAVecRestoreArrayDOF()
214: @*/
215: PetscErrorCode DMDAVecRestoreArrayDOF(DM da,Vec vec,void *array)
216: {
218: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
221: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
222: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
223: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
225: /* Handle case where user passes in global vector as opposed to local */
226: VecGetLocalSize(vec,&N);
227: if (N == xm*ym*zm*dof) {
228: gxm = xm;
229: gym = ym;
230: gzm = zm;
231: gxs = xs;
232: gys = ys;
233: gzs = zs;
234: }
236: if (dim == 1) {
237: VecRestoreArray2d(vec,gxm,dof,gxs,0,(PetscScalar***)array);
238: } else if (dim == 2) {
239: VecRestoreArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
240: } else if (dim == 3) {
241: VecRestoreArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
242: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
243: return(0);
244: }
246: /*@C
247: DMDAVecGetArrayRead - Returns a multiple dimension array that shares data with
248: the underlying vector and is indexed using the global dimensions.
250: Not collective
252: Input Parameter:
253: + da - the distributed array
254: - vec - the vector, either a vector the same size as one obtained with DMCreateGlobalVector() or DMCreateLocalVector()
256: Output Parameter:
257: . array - the array
259: Notes:
260: Call DMDAVecRestoreArrayRead() once you have finished accessing the vector entries.
262: In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]!
264: If vec is a local vector (obtained with DMCreateLocalVector() etc) then the ghost point locations are accessible. If it is
265: a global vector then the ghost points are not accessible. Of course with the local vector you will have had to do the
267: appropriate DMGlobalToLocalBegin() and DMGlobalToLocalEnd() to have correct values in the ghost locations.
269: Fortran Notes:
270: From Fortran use DMDAVecGetArrayReadF90() and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
271: 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
272: 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
273: array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from
274: DMDAGetCorners() for a global array or DMDAGetGhostCorners() for a local array. Include petsc/finclude/petscdmda.h90 to access this routine.
276: Due to bugs in the compiler DMDAVecGetArrayReadF90() does not work with gfortran versions before 4.5
278: Level: intermediate
280: .keywords: distributed array, get, corners, nodes, local indices, coordinates
282: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecRestoreArrayDOF()
283: DMDAVecGetArrayDOF()
284: @*/
285: PetscErrorCode DMDAVecGetArrayRead(DM da,Vec vec,void *array)
286: {
288: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
294: if (da->defaultSection) {
295: VecGetArrayRead(vec,(const PetscScalar**)array);
296: return(0);
297: }
298: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
299: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
300: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
302: /* Handle case where user passes in global vector as opposed to local */
303: VecGetLocalSize(vec,&N);
304: if (N == xm*ym*zm*dof) {
305: gxm = xm;
306: gym = ym;
307: gzm = zm;
308: gxs = xs;
309: gys = ys;
310: gzs = zs;
311: } else if (N != gxm*gym*gzm*dof) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Vector local size %D is not compatible with DMDA local sizes %D %D\n",N,xm*ym*zm*dof,gxm*gym*gzm*dof);
313: if (dim == 1) {
314: VecGetArray1dRead(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
315: } else if (dim == 2) {
316: VecGetArray2dRead(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
317: } else if (dim == 3) {
318: VecGetArray3dRead(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
319: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
320: return(0);
321: }
323: /*@
324: DMDAVecRestoreArrayRead - Restores a multiple dimension array obtained with DMDAVecGetArrayRead()
326: Not collective
328: Input Parameter:
329: + da - the distributed array
330: . vec - the vector, either a vector the same size as one obtained with
331: DMCreateGlobalVector() or DMCreateLocalVector()
332: - array - the array, non-NULL pointer is zeroed
334: Level: intermediate
336: Fortran Notes:
337: From Fortran use DMDAVecRestoreArrayReadF90()
339: .keywords: distributed array, get, corners, nodes, local indices, coordinates
341: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray()
342: @*/
343: PetscErrorCode DMDAVecRestoreArrayRead(DM da,Vec vec,void *array)
344: {
346: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
352: if (da->defaultSection) {
353: VecRestoreArrayRead(vec,(const PetscScalar**)array);
354: return(0);
355: }
356: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
357: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
358: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
360: /* Handle case where user passes in global vector as opposed to local */
361: VecGetLocalSize(vec,&N);
362: if (N == xm*ym*zm*dof) {
363: gxm = xm;
364: gym = ym;
365: gzm = zm;
366: gxs = xs;
367: gys = ys;
368: gzs = zs;
369: } else if (N != gxm*gym*gzm*dof) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Vector local size %D is not compatible with DMDA local sizes %D %D\n",N,xm*ym*zm*dof,gxm*gym*gzm*dof);
371: if (dim == 1) {
372: VecRestoreArray1dRead(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
373: } else if (dim == 2) {
374: VecRestoreArray2dRead(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
375: } else if (dim == 3) {
376: VecRestoreArray3dRead(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
377: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
378: return(0);
379: }
381: /*@C
382: DMDAVecGetArrayDOFRead - Returns a multiple dimension array that shares data with
383: the underlying vector and is indexed using the global dimensions.
385: Not Collective
387: Input Parameter:
388: + da - the distributed array
389: - vec - the vector, either a vector the same size as one obtained with
390: DMCreateGlobalVector() or DMCreateLocalVector()
392: Output Parameter:
393: . array - the array
395: Notes:
396: Call DMDAVecRestoreArrayDOFRead() once you have finished accessing the vector entries.
398: In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]!
400: In Fortran 90 you do not need a version of DMDAVecRestoreArrayDOF() just use DMDAVecRestoreArrayReadF90() and declare your array with one higher dimension,
401: see src/dm/examples/tutorials/ex11f90.F
403: Level: intermediate
405: .keywords: distributed array, get, corners, nodes, local indices, coordinates
407: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecGetArray(), DMDAVecRestoreArrayDOF()
408: @*/
409: PetscErrorCode DMDAVecGetArrayDOFRead(DM da,Vec vec,void *array)
410: {
412: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
415: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
416: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
417: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
419: /* Handle case where user passes in global vector as opposed to local */
420: VecGetLocalSize(vec,&N);
421: if (N == xm*ym*zm*dof) {
422: gxm = xm;
423: gym = ym;
424: gzm = zm;
425: gxs = xs;
426: gys = ys;
427: gzs = zs;
428: } else if (N != gxm*gym*gzm*dof) SETERRQ3(PETSC_COMM_SELF,PETSC_ERR_ARG_INCOMP,"Vector local size %D is not compatible with DMDA local sizes %D %D\n",N,xm*ym*zm*dof,gxm*gym*gzm*dof);
430: if (dim == 1) {
431: VecGetArray2dRead(vec,gxm,dof,gxs,0,(PetscScalar***)array);
432: } else if (dim == 2) {
433: VecGetArray3dRead(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
434: } else if (dim == 3) {
435: VecGetArray4dRead(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
436: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
437: return(0);
438: }
440: /*@
441: DMDAVecRestoreArrayDOFRead - Restores a multiple dimension array obtained with DMDAVecGetArrayDOFRead()
443: Not Collective
445: Input Parameter:
446: + da - the distributed array
447: . vec - the vector, either a vector the same size as one obtained with
448: DMCreateGlobalVector() or DMCreateLocalVector()
449: - array - the array
451: Level: intermediate
453: .keywords: distributed array, get, corners, nodes, local indices, coordinates
455: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(), DMDAVecRestoreArrayDOF()
456: @*/
457: PetscErrorCode DMDAVecRestoreArrayDOFRead(DM da,Vec vec,void *array)
458: {
460: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
463: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
464: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
465: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
467: /* Handle case where user passes in global vector as opposed to local */
468: VecGetLocalSize(vec,&N);
469: if (N == xm*ym*zm*dof) {
470: gxm = xm;
471: gym = ym;
472: gzm = zm;
473: gxs = xs;
474: gys = ys;
475: gzs = zs;
476: }
478: if (dim == 1) {
479: VecRestoreArray2dRead(vec,gxm,dof,gxs,0,(PetscScalar***)array);
480: } else if (dim == 2) {
481: VecRestoreArray3dRead(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
482: } else if (dim == 3) {
483: VecRestoreArray4dRead(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
484: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
485: return(0);
486: }