Actual source code: dagetarray.c
petsc-3.6.1 2015-08-06
2: #include <petsc/private/dmdaimpl.h> /*I "petscdmda.h" I*/
6: /*@C
7: DMDAVecGetArray - Returns a multiple dimension array that shares data with
8: the underlying vector and is indexed using the global dimensions.
10: Logically collective on Vec
12: Input Parameter:
13: + da - the distributed array
14: - vec - the vector, either a vector the same size as one obtained with DMCreateGlobalVector() or DMCreateLocalVector()
16: Output Parameter:
17: . array - the array
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: From Fortran use DMDAVecGetArrayF90() and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
30: 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
31: 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
32: array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from
33: DMDAGetCorners() for a global array or DMDAGetGhostCorners() for a local array. Include petsc/finclude/petscdmda.h90 to access this routine.
35: Due to bugs in the compiler DMDAVecGetArrayF90() does not work with gfortran versions before 4.5
37: Level: intermediate
39: .keywords: distributed array, get, corners, nodes, local indices, coordinates
41: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecRestoreArrayDOF()
42: DMDAVecGetArrayDOF()
43: @*/
44: PetscErrorCode DMDAVecGetArray(DM da,Vec vec,void *array)
45: {
47: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
53: if (da->defaultSection) {
54: VecGetArray(vec,(PetscScalar**)array);
55: return(0);
56: }
57: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
58: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
59: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
61: /* Handle case where user passes in global vector as opposed to local */
62: 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 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);
72: if (dim == 1) {
73: VecGetArray1d(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
74: } else if (dim == 2) {
75: VecGetArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
76: } else if (dim == 3) {
77: VecGetArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
78: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
79: return(0);
80: }
84: /*@
85: DMDAVecRestoreArray - Restores a multiple dimension array obtained with DMDAVecGetArray()
87: Logically collective on Vec
89: Input Parameter:
90: + da - the distributed array
91: . vec - the vector, either a vector the same size as one obtained with
92: DMCreateGlobalVector() or DMCreateLocalVector()
93: - array - the array, non-NULL pointer is zeroed
95: Level: intermediate
97: Fortran Notes: From Fortran use DMDAVecRestoreArrayF90()
99: .keywords: distributed array, get, corners, nodes, local indices, coordinates
101: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray()
102: @*/
103: PetscErrorCode DMDAVecRestoreArray(DM da,Vec vec,void *array)
104: {
106: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
112: if (da->defaultSection) {
113: VecRestoreArray(vec,(PetscScalar**)array);
114: return(0);
115: }
116: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
117: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
118: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
120: /* Handle case where user passes in global vector as opposed to local */
121: VecGetLocalSize(vec,&N);
122: if (N == xm*ym*zm*dof) {
123: gxm = xm;
124: gym = ym;
125: gzm = zm;
126: gxs = xs;
127: gys = ys;
128: gzs = zs;
129: } 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);
131: if (dim == 1) {
132: VecRestoreArray1d(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
133: } else if (dim == 2) {
134: VecRestoreArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
135: } else if (dim == 3) {
136: VecRestoreArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
137: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
138: return(0);
139: }
143: /*@C
144: DMDAVecGetArrayDOF - Returns a multiple dimension array that shares data with
145: the underlying vector and is indexed using the global dimensions.
147: Logically collective
149: Input Parameter:
150: + da - the distributed array
151: - vec - the vector, either a vector the same size as one obtained with
152: DMCreateGlobalVector() or DMCreateLocalVector()
154: Output Parameter:
155: . array - the array
157: Notes:
158: Call DMDAVecRestoreArrayDOF() once you have finished accessing the vector entries.
160: In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]!
162: In Fortran 90 you do not need a version of DMDAVecRestoreArrayDOF() just use DMDAVecRestoreArrayF90() and declare your array with one higher dimension,
163: see src/dm/examples/tutorials/ex11f90.F
165: Level: intermediate
167: .keywords: distributed array, get, corners, nodes, local indices, coordinates
169: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecGetArray(), DMDAVecRestoreArrayDOF()
170: @*/
171: PetscErrorCode DMDAVecGetArrayDOF(DM da,Vec vec,void *array)
172: {
174: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
177: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
178: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
179: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
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;
190: } 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);
192: if (dim == 1) {
193: VecGetArray2d(vec,gxm,dof,gxs,0,(PetscScalar***)array);
194: } else if (dim == 2) {
195: VecGetArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
196: } else if (dim == 3) {
197: VecGetArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
198: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
199: return(0);
200: }
204: /*@
205: DMDAVecRestoreArrayDOF - Restores a multiple dimension array obtained with DMDAVecGetArrayDOF()
207: Logically collective
209: Input Parameter:
210: + da - the distributed array
211: . vec - the vector, either a vector the same size as one obtained with
212: DMCreateGlobalVector() or DMCreateLocalVector()
213: - array - the array
215: Level: intermediate
217: .keywords: distributed array, get, corners, nodes, local indices, coordinates
219: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(), DMDAVecRestoreArrayDOF()
220: @*/
221: PetscErrorCode DMDAVecRestoreArrayDOF(DM da,Vec vec,void *array)
222: {
224: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
227: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
228: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
229: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
231: /* Handle case where user passes in global vector as opposed to local */
232: VecGetLocalSize(vec,&N);
233: if (N == xm*ym*zm*dof) {
234: gxm = xm;
235: gym = ym;
236: gzm = zm;
237: gxs = xs;
238: gys = ys;
239: gzs = zs;
240: }
242: if (dim == 1) {
243: VecRestoreArray2d(vec,gxm,dof,gxs,0,(PetscScalar***)array);
244: } else if (dim == 2) {
245: VecRestoreArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
246: } else if (dim == 3) {
247: VecRestoreArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
248: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
249: return(0);
250: }
254: /*@C
255: DMDAVecGetArrayRead - Returns a multiple dimension array that shares data with
256: the underlying vector and is indexed using the global dimensions.
258: Not collective
260: Input Parameter:
261: + da - the distributed array
262: - vec - the vector, either a vector the same size as one obtained with DMCreateGlobalVector() or DMCreateLocalVector()
264: Output Parameter:
265: . array - the array
267: Notes:
268: Call DMDAVecRestoreArrayRead() once you have finished accessing the vector entries.
270: In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]!
272: If vec is a local vector (obtained with DMCreateLocalVector() etc) then the ghost point locations are accessible. If it is
273: a global vector then the ghost points are not accessible. Of course with the local vector you will have had to do the
275: appropriate DMGlobalToLocalBegin() and DMGlobalToLocalEnd() to have correct values in the ghost locations.
277: Fortran Notes: From Fortran use DMDAVecGetArrayReadF90() and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
278: 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
279: 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
280: array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from
281: DMDAGetCorners() for a global array or DMDAGetGhostCorners() for a local array. Include petsc/finclude/petscdmda.h90 to access this routine.
283: Due to bugs in the compiler DMDAVecGetArrayReadF90() does not work with gfortran versions before 4.5
285: Level: intermediate
287: .keywords: distributed array, get, corners, nodes, local indices, coordinates
289: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecRestoreArrayDOF()
290: DMDAVecGetArrayDOF()
291: @*/
292: PetscErrorCode DMDAVecGetArrayRead(DM da,Vec vec,void *array)
293: {
295: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
301: if (da->defaultSection) {
302: VecGetArrayRead(vec,(const PetscScalar**)array);
303: return(0);
304: }
305: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
306: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
307: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
309: /* Handle case where user passes in global vector as opposed to local */
310: VecGetLocalSize(vec,&N);
311: if (N == xm*ym*zm*dof) {
312: gxm = xm;
313: gym = ym;
314: gzm = zm;
315: gxs = xs;
316: gys = ys;
317: gzs = zs;
318: } 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);
320: if (dim == 1) {
321: VecGetArray1dRead(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
322: } else if (dim == 2) {
323: VecGetArray2dRead(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
324: } else if (dim == 3) {
325: VecGetArray3dRead(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
326: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
327: return(0);
328: }
332: /*@
333: DMDAVecRestoreArrayRead - Restores a multiple dimension array obtained with DMDAVecGetArrayRead()
335: Not collective
337: Input Parameter:
338: + da - the distributed array
339: . vec - the vector, either a vector the same size as one obtained with
340: DMCreateGlobalVector() or DMCreateLocalVector()
341: - array - the array, non-NULL pointer is zeroed
343: Level: intermediate
345: Fortran Notes: From Fortran use DMDAVecRestoreArrayReadF90()
347: .keywords: distributed array, get, corners, nodes, local indices, coordinates
349: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray()
350: @*/
351: PetscErrorCode DMDAVecRestoreArrayRead(DM da,Vec vec,void *array)
352: {
354: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
360: if (da->defaultSection) {
361: VecRestoreArrayRead(vec,(const PetscScalar**)array);
362: return(0);
363: }
364: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
365: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
366: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
368: /* Handle case where user passes in global vector as opposed to local */
369: VecGetLocalSize(vec,&N);
370: if (N == xm*ym*zm*dof) {
371: gxm = xm;
372: gym = ym;
373: gzm = zm;
374: gxs = xs;
375: gys = ys;
376: gzs = zs;
377: } 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);
379: if (dim == 1) {
380: VecRestoreArray1dRead(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
381: } else if (dim == 2) {
382: VecRestoreArray2dRead(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
383: } else if (dim == 3) {
384: VecRestoreArray3dRead(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
385: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
386: return(0);
387: }
391: /*@C
392: DMDAVecGetArrayDOFRead - Returns a multiple dimension array that shares data with
393: the underlying vector and is indexed using the global dimensions.
395: Not Collective
397: Input Parameter:
398: + da - the distributed array
399: - vec - the vector, either a vector the same size as one obtained with
400: DMCreateGlobalVector() or DMCreateLocalVector()
402: Output Parameter:
403: . array - the array
405: Notes:
406: Call DMDAVecRestoreArrayDOFRead() once you have finished accessing the vector entries.
408: In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]!
410: In Fortran 90 you do not need a version of DMDAVecRestoreArrayDOF() just use DMDAVecRestoreArrayReadF90() and declare your array with one higher dimension,
411: see src/dm/examples/tutorials/ex11f90.F
413: Level: intermediate
415: .keywords: distributed array, get, corners, nodes, local indices, coordinates
417: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecGetArray(), DMDAVecRestoreArrayDOF()
418: @*/
419: PetscErrorCode DMDAVecGetArrayDOFRead(DM da,Vec vec,void *array)
420: {
422: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
425: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
426: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
427: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
429: /* Handle case where user passes in global vector as opposed to local */
430: VecGetLocalSize(vec,&N);
431: if (N == xm*ym*zm*dof) {
432: gxm = xm;
433: gym = ym;
434: gzm = zm;
435: gxs = xs;
436: gys = ys;
437: gzs = zs;
438: } 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);
440: if (dim == 1) {
441: VecGetArray2dRead(vec,gxm,dof,gxs,0,(PetscScalar***)array);
442: } else if (dim == 2) {
443: VecGetArray3dRead(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
444: } else if (dim == 3) {
445: VecGetArray4dRead(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
446: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
447: return(0);
448: }
452: /*@
453: DMDAVecRestoreArrayDOFRead - Restores a multiple dimension array obtained with DMDAVecGetArrayDOFRead()
455: Not Collective
457: Input Parameter:
458: + da - the distributed array
459: . vec - the vector, either a vector the same size as one obtained with
460: DMCreateGlobalVector() or DMCreateLocalVector()
461: - array - the array
463: Level: intermediate
465: .keywords: distributed array, get, corners, nodes, local indices, coordinates
467: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(), DMDAVecRestoreArrayDOF()
468: @*/
469: PetscErrorCode DMDAVecRestoreArrayDOFRead(DM da,Vec vec,void *array)
470: {
472: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
475: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
476: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
477: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
479: /* Handle case where user passes in global vector as opposed to local */
480: VecGetLocalSize(vec,&N);
481: if (N == xm*ym*zm*dof) {
482: gxm = xm;
483: gym = ym;
484: gzm = zm;
485: gxs = xs;
486: gys = ys;
487: gzs = zs;
488: }
490: if (dim == 1) {
491: VecRestoreArray2dRead(vec,gxm,dof,gxs,0,(PetscScalar***)array);
492: } else if (dim == 2) {
493: VecRestoreArray3dRead(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
494: } else if (dim == 3) {
495: VecRestoreArray4dRead(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
496: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
497: return(0);
498: }