Actual source code: dagetarray.c
petsc-3.12.0 2019-09-29
2: #include <petsc/private/dmdaimpl.h>
4: /*@C
5: DMDAVecGetArray - Returns a multiple dimension array that shares data with
6: the underlying vector and is indexed using the global dimensions.
8: Logically collective on da
10: Input 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: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecRestoreArrayDOF()
39: DMDAVecGetArrayDOF(), DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead()
40: @*/
41: PetscErrorCode DMDAVecGetArray(DM da,Vec vec,void *array)
42: {
44: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
50: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
51: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
52: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
54: /* Handle case where user passes in global vector as opposed to local */
55: VecGetLocalSize(vec,&N);
56: if (N == xm*ym*zm*dof) {
57: gxm = xm;
58: gym = ym;
59: gzm = zm;
60: gxs = xs;
61: gys = ys;
62: gzs = zs;
63: } 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);
65: if (dim == 1) {
66: VecGetArray1d(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
67: } else if (dim == 2) {
68: VecGetArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
69: } else if (dim == 3) {
70: VecGetArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
71: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
72: return(0);
73: }
75: /*@
76: DMDAVecRestoreArray - Restores a multiple dimension array obtained with DMDAVecGetArray()
78: Logically collective on da
80: Input Parameter:
81: + da - the distributed array
82: . vec - the vector, either a vector the same size as one obtained with
83: DMCreateGlobalVector() or DMCreateLocalVector()
84: - array - the array, non-NULL pointer is zeroed
86: Level: intermediate
88: Fortran Notes:
89: From Fortran use DMDAVecRestoreArayF90()
91: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(),
92: DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead()
93: @*/
94: PetscErrorCode DMDAVecRestoreArray(DM da,Vec vec,void *array)
95: {
97: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
103: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
104: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
105: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
107: /* Handle case where user passes in global vector as opposed to local */
108: VecGetLocalSize(vec,&N);
109: if (N == xm*ym*zm*dof) {
110: gxm = xm;
111: gym = ym;
112: gzm = zm;
113: gxs = xs;
114: gys = ys;
115: gzs = zs;
116: } 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);
118: if (dim == 1) {
119: VecRestoreArray1d(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
120: } else if (dim == 2) {
121: VecRestoreArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
122: } else if (dim == 3) {
123: VecRestoreArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
124: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
125: return(0);
126: }
128: /*@C
129: DMDAVecGetArrayWrite - Returns a multiple dimension array that shares data with
130: the underlying vector and is indexed using the global dimensions.
132: Logically collective on Vec
134: Input Parameter:
135: + da - the distributed array
136: - vec - the vector, either a vector the same size as one obtained with DMCreateGlobalVector() or DMCreateLocalVector()
138: Output Parameter:
139: . array - the array
141: Notes:
142: Call DMDAVecRestoreArray() once you have finished accessing the vector entries.
144: In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]!
146: If vec is a local vector (obtained with DMCreateLocalVector() etc) then the ghost point locations are accessible. If it is
147: a global vector then the ghost points are not accessible. Of course with the local vector you will have had to do the
149: appropriate DMGlobalToLocalBegin() and DMGlobalToLocalEnd() to have correct values in the ghost locations.
151: Fortran Notes:
152: From Fortran use DMDAVecGetArrayF90() and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
153: dimension. For a DMDA created with a dof of 1 use the dimension of the DMDA, for a DMDA created with a dof greater than 1 use one more than the
154: dimension of the DMDA. The order of the indices is array(xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) (when dof is 1) otherwise
155: array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from
156: DMDAGetCorners() for a global array or DMDAGetGhostCorners() for a local array. Include petsc/finclude/petscdmda.h90 to access this routine.
158: Due to bugs in the compiler DMDAVecGetArrayF90() does not work with gfortran versions before 4.5
160: Level: intermediate
162: Developer Notes: This has code duplication with DMDAVecGetArray() and DMDAVecGetArrayRead()
164: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArrayWrite(), DMDAVecRestoreArrayDOF()
165: DMDAVecGetArrayDOF(), DMDAVecGetArray(), DMDAVecRestoreArray(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead()
166: @*/
167: PetscErrorCode DMDAVecGetArrayWrite(DM da,Vec vec,void *array)
168: {
170: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
176: if (da->localSection) {
177: VecGetArrayWrite(vec,(PetscScalar**)array);
178: return(0);
179: }
180: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
181: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
182: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
184: /* Handle case where user passes in global vector as opposed to local */
185: VecGetLocalSize(vec,&N);
186: if (N == xm*ym*zm*dof) {
187: gxm = xm;
188: gym = ym;
189: gzm = zm;
190: gxs = xs;
191: gys = ys;
192: gzs = zs;
193: } 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);
195: if (dim == 1) {
196: VecGetArray1dWrite(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
197: } else if (dim == 2) {
198: VecGetArray2dWrite(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
199: } else if (dim == 3) {
200: VecGetArray3dWrite(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
201: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
202: return(0);
203: }
205: /*@
206: DMDAVecRestoreArrayWrite - Restores a multiple dimension array obtained with DMDAVecGetArrayWrite()
208: Logically collective on Vec
210: Input Parameter:
211: + da - the distributed array
212: . vec - the vector, either a vector the same size as one obtained with
213: DMCreateGlobalVector() or DMCreateLocalVector()
214: - array - the array, non-NULL pointer is zeroed
216: Level: intermediate
218: Fortran Notes:
219: From Fortran use DMDAVecRestoreArayF90()
221: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArrayWrite(),
222: DMDAVecGetArray(), DMDAVecRestoreArray(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead()
223: @*/
224: PetscErrorCode DMDAVecRestoreArrayWrite(DM da,Vec vec,void *array)
225: {
227: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
233: if (da->localSection) {
234: VecRestoreArray(vec,(PetscScalar**)array);
235: return(0);
236: }
237: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
238: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
239: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
241: /* Handle case where user passes in global vector as opposed to local */
242: VecGetLocalSize(vec,&N);
243: if (N == xm*ym*zm*dof) {
244: gxm = xm;
245: gym = ym;
246: gzm = zm;
247: gxs = xs;
248: gys = ys;
249: gzs = zs;
250: } 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);
252: if (dim == 1) {
253: VecRestoreArray1dWrite(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
254: } else if (dim == 2) {
255: VecRestoreArray2dWrite(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
256: } else if (dim == 3) {
257: VecRestoreArray3dWrite(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
258: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
259: return(0);
260: }
262: /*@C
263: DMDAVecGetArrayDOF - Returns a multiple dimension array that shares data with
264: the underlying vector and is indexed using the global dimensions.
266: Logically collective
268: Input Parameter:
269: + da - the distributed array
270: - vec - the vector, either a vector the same size as one obtained with
271: DMCreateGlobalVector() or DMCreateLocalVector()
273: Output Parameter:
274: . array - the array
276: Notes:
277: Call DMDAVecRestoreArrayDOF() once you have finished accessing the vector entries.
279: In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]!
281: In Fortran 90 you do not need a version of DMDAVecRestoreArrayDOF() just use DMDAVecRestoreArrayF90() and declare your array with one higher dimension,
282: see src/dm/examples/tutorials/ex11f90.F
284: Level: intermediate
286: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecGetArray(), DMDAVecRestoreArrayDOF(),
287: DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead()
288: @*/
289: PetscErrorCode DMDAVecGetArrayDOF(DM da,Vec vec,void *array)
290: {
292: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
295: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
296: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
297: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
299: /* Handle case where user passes in global vector as opposed to local */
300: VecGetLocalSize(vec,&N);
301: if (N == xm*ym*zm*dof) {
302: gxm = xm;
303: gym = ym;
304: gzm = zm;
305: gxs = xs;
306: gys = ys;
307: gzs = zs;
308: } 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);
310: if (dim == 1) {
311: VecGetArray2d(vec,gxm,dof,gxs,0,(PetscScalar***)array);
312: } else if (dim == 2) {
313: VecGetArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
314: } else if (dim == 3) {
315: VecGetArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
316: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
317: return(0);
318: }
320: /*@
321: DMDAVecRestoreArrayDOF - Restores a multiple dimension array obtained with DMDAVecGetArrayDOF()
323: Logically collective
325: Input Parameter:
326: + da - the distributed array
327: . vec - the vector, either a vector the same size as one obtained with
328: DMCreateGlobalVector() or DMCreateLocalVector()
329: - array - the array
331: Level: intermediate
333: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(), DMDAVecRestoreArrayDOF(),
334: DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead()
335: @*/
336: PetscErrorCode DMDAVecRestoreArrayDOF(DM da,Vec vec,void *array)
337: {
339: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
342: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
343: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
344: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
346: /* Handle case where user passes in global vector as opposed to local */
347: VecGetLocalSize(vec,&N);
348: if (N == xm*ym*zm*dof) {
349: gxm = xm;
350: gym = ym;
351: gzm = zm;
352: gxs = xs;
353: gys = ys;
354: gzs = zs;
355: }
357: if (dim == 1) {
358: VecRestoreArray2d(vec,gxm,dof,gxs,0,(PetscScalar***)array);
359: } else if (dim == 2) {
360: VecRestoreArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
361: } else if (dim == 3) {
362: VecRestoreArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
363: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
364: return(0);
365: }
367: /*@C
368: DMDAVecGetArrayRead - Returns a multiple dimension array that shares data with
369: the underlying vector and is indexed using the global dimensions.
371: Not collective
373: Input Parameter:
374: + da - the distributed array
375: - vec - the vector, either a vector the same size as one obtained with DMCreateGlobalVector() or DMCreateLocalVector()
377: Output Parameter:
378: . array - the array
380: Notes:
381: Call DMDAVecRestoreArrayRead() once you have finished accessing the vector entries.
383: In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]!
385: If vec is a local vector (obtained with DMCreateLocalVector() etc) then the ghost point locations are accessible. If it is
386: a global vector then the ghost points are not accessible. Of course with the local vector you will have had to do the
388: appropriate DMGlobalToLocalBegin() and DMGlobalToLocalEnd() to have correct values in the ghost locations.
390: Fortran Notes:
391: From Fortran use DMDAVecGetArrayReadF90() and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
392: 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
393: 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
394: array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from
395: DMDAGetCorners() for a global array or DMDAGetGhostCorners() for a local array. Include petsc/finclude/petscdmda.h90 to access this routine.
397: Due to bugs in the compiler DMDAVecGetArrayReadF90() does not work with gfortran versions before 4.5
399: Level: intermediate
401: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArrayRead(), DMDAVecRestoreArrayDOF()
402: DMDAVecGetArrayDOF(), DMDAVecGetArray(), DMDAVecRestoreArray(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead()
403: @*/
404: PetscErrorCode DMDAVecGetArrayRead(DM da,Vec vec,void *array)
405: {
407: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
413: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
414: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
415: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
417: /* Handle case where user passes in global vector as opposed to local */
418: VecGetLocalSize(vec,&N);
419: if (N == xm*ym*zm*dof) {
420: gxm = xm;
421: gym = ym;
422: gzm = zm;
423: gxs = xs;
424: gys = ys;
425: gzs = zs;
426: } 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);
428: if (dim == 1) {
429: VecGetArray1dRead(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
430: } else if (dim == 2) {
431: VecGetArray2dRead(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
432: } else if (dim == 3) {
433: VecGetArray3dRead(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
434: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
435: return(0);
436: }
438: /*@
439: DMDAVecRestoreArrayRead - Restores a multiple dimension array obtained with DMDAVecGetArrayRead()
441: Not collective
443: Input Parameter:
444: + da - the distributed array
445: . vec - the vector, either a vector the same size as one obtained with
446: DMCreateGlobalVector() or DMCreateLocalVector()
447: - array - the array, non-NULL pointer is zeroed
449: Level: intermediate
451: Fortran Notes:
452: From Fortran use DMDAVecRestoreArrayReadF90()
454: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArrayRead(),
455: DMDAVecGetArray(), DMDAVecRestoreArray(), DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite()
456: @*/
457: PetscErrorCode DMDAVecRestoreArrayRead(DM da,Vec vec,void *array)
458: {
460: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
466: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
467: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
468: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
470: /* Handle case where user passes in global vector as opposed to local */
471: VecGetLocalSize(vec,&N);
472: if (N == xm*ym*zm*dof) {
473: gxm = xm;
474: gym = ym;
475: gzm = zm;
476: gxs = xs;
477: gys = ys;
478: gzs = zs;
479: } 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);
481: if (dim == 1) {
482: VecRestoreArray1dRead(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
483: } else if (dim == 2) {
484: VecRestoreArray2dRead(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
485: } else if (dim == 3) {
486: VecRestoreArray3dRead(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
487: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
488: return(0);
489: }
491: /*@C
492: DMDAVecGetArrayDOFRead - Returns a multiple dimension array that shares data with
493: the underlying vector and is indexed using the global dimensions.
495: Not Collective
497: Input Parameter:
498: + da - the distributed array
499: - vec - the vector, either a vector the same size as one obtained with
500: DMCreateGlobalVector() or DMCreateLocalVector()
502: Output Parameter:
503: . array - the array
505: Notes:
506: Call DMDAVecRestoreArrayDOFRead() once you have finished accessing the vector entries.
508: In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]!
510: In Fortran 90 you do not need a version of DMDAVecRestoreArrayDOF() just use DMDAVecRestoreArrayReadF90() and declare your array with one higher dimension,
511: see src/dm/examples/tutorials/ex11f90.F
513: Level: intermediate
515: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecGetArray(), DMDAVecRestoreArrayDOF(),
516: DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead()
517: @*/
518: PetscErrorCode DMDAVecGetArrayDOFRead(DM da,Vec vec,void *array)
519: {
521: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
524: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
525: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
526: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
528: /* Handle case where user passes in global vector as opposed to local */
529: VecGetLocalSize(vec,&N);
530: if (N == xm*ym*zm*dof) {
531: gxm = xm;
532: gym = ym;
533: gzm = zm;
534: gxs = xs;
535: gys = ys;
536: gzs = zs;
537: } 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);
539: if (dim == 1) {
540: VecGetArray2dRead(vec,gxm,dof,gxs,0,(PetscScalar***)array);
541: } else if (dim == 2) {
542: VecGetArray3dRead(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
543: } else if (dim == 3) {
544: VecGetArray4dRead(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
545: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
546: return(0);
547: }
549: /*@
550: DMDAVecRestoreArrayDOFRead - Restores a multiple dimension array obtained with DMDAVecGetArrayDOFRead()
552: Not Collective
554: Input Parameter:
555: + da - the distributed array
556: . vec - the vector, either a vector the same size as one obtained with
557: DMCreateGlobalVector() or DMCreateLocalVector()
558: - array - the array
560: Level: intermediate
562: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(), DMDAVecRestoreArrayDOF(),
563: DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead()
564: @*/
565: PetscErrorCode DMDAVecRestoreArrayDOFRead(DM da,Vec vec,void *array)
566: {
568: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
571: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
572: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
573: DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);
575: /* Handle case where user passes in global vector as opposed to local */
576: VecGetLocalSize(vec,&N);
577: if (N == xm*ym*zm*dof) {
578: gxm = xm;
579: gym = ym;
580: gzm = zm;
581: gxs = xs;
582: gys = ys;
583: gzs = zs;
584: }
586: if (dim == 1) {
587: VecRestoreArray2dRead(vec,gxm,dof,gxs,0,(PetscScalar***)array);
588: } else if (dim == 2) {
589: VecRestoreArray3dRead(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
590: } else if (dim == 3) {
591: VecRestoreArray4dRead(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
592: } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
593: return(0);
594: }