Actual source code: dagetarray.c
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 Parameters:
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: DMStagVecGetArray()
41: @*/
42: PetscErrorCode DMDAVecGetArray(DM da,Vec vec,void *array)
43: {
44: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
49: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
50: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
51: DMDAGetInfo(da,&dim,NULL,NULL,NULL,NULL,NULL,NULL,&dof,NULL,NULL,NULL,NULL,NULL);
53: /* Handle case where user passes in global vector as opposed to local */
54: VecGetLocalSize(vec,&N);
55: if (N == xm*ym*zm*dof) {
56: gxm = xm;
57: gym = ym;
58: gzm = zm;
59: gxs = xs;
60: gys = ys;
61: gzs = zs;
64: if (dim == 1) {
65: VecGetArray1d(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
66: } else if (dim == 2) {
67: VecGetArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
68: } else if (dim == 3) {
69: VecGetArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
70: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D",dim);
71: return 0;
72: }
74: /*@
75: DMDAVecRestoreArray - Restores a multiple dimension array obtained with DMDAVecGetArray()
77: Logically collective on da
79: Input Parameters:
80: + da - the distributed array
81: . vec - the vector, either a vector the same size as one obtained with
82: DMCreateGlobalVector() or DMCreateLocalVector()
83: - array - the array, non-NULL pointer is zeroed
85: Level: intermediate
87: Fortran Notes:
88: From Fortran use DMDAVecRestoreArayF90()
90: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(),
91: DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead(),
92: DMDStagVecRestoreArray()
93: @*/
94: PetscErrorCode DMDAVecRestoreArray(DM da,Vec vec,void *array)
95: {
96: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
101: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
102: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
103: DMDAGetInfo(da,&dim,NULL,NULL,NULL,NULL,NULL,NULL,&dof,NULL,NULL,NULL,NULL,NULL);
105: /* Handle case where user passes in global vector as opposed to local */
106: VecGetLocalSize(vec,&N);
107: if (N == xm*ym*zm*dof) {
108: gxm = xm;
109: gym = ym;
110: gzm = zm;
111: gxs = xs;
112: gys = ys;
113: gzs = zs;
116: if (dim == 1) {
117: VecRestoreArray1d(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
118: } else if (dim == 2) {
119: VecRestoreArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
120: } else if (dim == 3) {
121: VecRestoreArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
122: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D",dim);
123: return 0;
124: }
126: /*@C
127: DMDAVecGetArrayWrite - Returns a multiple dimension array that shares data with
128: the underlying vector and is indexed using the global dimensions.
130: Logically collective on Vec
132: Input Parameters:
133: + da - the distributed array
134: - vec - the vector, either a vector the same size as one obtained with DMCreateGlobalVector() or DMCreateLocalVector()
136: Output Parameter:
137: . array - the array
139: Notes:
140: Call DMDAVecRestoreArray() once you have finished accessing the vector entries.
142: In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]!
144: If vec is a local vector (obtained with DMCreateLocalVector() etc) then the ghost point locations are accessible. If it is
145: a global vector then the ghost points are not accessible. Of course with the local vector you will have had to do the
147: appropriate DMGlobalToLocalBegin() and DMGlobalToLocalEnd() to have correct values in the ghost locations.
149: Fortran Notes:
150: From Fortran use DMDAVecGetArrayF90() and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
151: 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
152: 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
153: array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from
154: DMDAGetCorners() for a global array or DMDAGetGhostCorners() for a local array. Include petsc/finclude/petscdmda.h90 to access this routine.
156: Due to bugs in the compiler DMDAVecGetArrayF90() does not work with gfortran versions before 4.5
158: Level: intermediate
160: Developer Notes: This has code duplication with DMDAVecGetArray() and DMDAVecGetArrayRead()
162: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArrayWrite(), DMDAVecRestoreArrayDOF()
163: DMDAVecGetArrayDOF(), DMDAVecGetArray(), DMDAVecRestoreArray(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead()
164: @*/
165: PetscErrorCode DMDAVecGetArrayWrite(DM da,Vec vec,void *array)
166: {
167: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
172: if (da->localSection) {
173: VecGetArrayWrite(vec,(PetscScalar**)array);
174: return 0;
175: }
176: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
177: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
178: DMDAGetInfo(da,&dim,NULL,NULL,NULL,NULL,NULL,NULL,&dof,NULL,NULL,NULL,NULL,NULL);
180: /* Handle case where user passes in global vector as opposed to local */
181: VecGetLocalSize(vec,&N);
182: if (N == xm*ym*zm*dof) {
183: gxm = xm;
184: gym = ym;
185: gzm = zm;
186: gxs = xs;
187: gys = ys;
188: gzs = zs;
191: if (dim == 1) {
192: VecGetArray1dWrite(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
193: } else if (dim == 2) {
194: VecGetArray2dWrite(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
195: } else if (dim == 3) {
196: VecGetArray3dWrite(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
197: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D",dim);
198: return 0;
199: }
201: /*@
202: DMDAVecRestoreArrayWrite - Restores a multiple dimension array obtained with DMDAVecGetArrayWrite()
204: Logically collective on Vec
206: Input Parameters:
207: + da - the distributed array
208: . vec - the vector, either a vector the same size as one obtained with
209: DMCreateGlobalVector() or DMCreateLocalVector()
210: - array - the array, non-NULL pointer is zeroed
212: Level: intermediate
214: Fortran Notes:
215: From Fortran use DMDAVecRestoreArayF90()
217: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArrayWrite(),
218: DMDAVecGetArray(), DMDAVecRestoreArray(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead()
219: @*/
220: PetscErrorCode DMDAVecRestoreArrayWrite(DM da,Vec vec,void *array)
221: {
222: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
227: if (da->localSection) {
228: VecRestoreArray(vec,(PetscScalar**)array);
229: return 0;
230: }
231: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
232: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
233: DMDAGetInfo(da,&dim,NULL,NULL,NULL,NULL,NULL,NULL,&dof,NULL,NULL,NULL,NULL,NULL);
235: /* Handle case where user passes in global vector as opposed to local */
236: VecGetLocalSize(vec,&N);
237: if (N == xm*ym*zm*dof) {
238: gxm = xm;
239: gym = ym;
240: gzm = zm;
241: gxs = xs;
242: gys = ys;
243: gzs = zs;
246: if (dim == 1) {
247: VecRestoreArray1dWrite(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
248: } else if (dim == 2) {
249: VecRestoreArray2dWrite(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
250: } else if (dim == 3) {
251: VecRestoreArray3dWrite(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
252: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D",dim);
253: return 0;
254: }
256: /*@C
257: DMDAVecGetArrayDOF - Returns a multiple dimension array that shares data with
258: the underlying vector and is indexed using the global dimensions.
260: Logically collective
262: Input Parameters:
263: + da - the distributed array
264: - vec - the vector, either a vector the same size as one obtained with
265: DMCreateGlobalVector() or DMCreateLocalVector()
267: Output Parameter:
268: . array - the array
270: Notes:
271: Call DMDAVecRestoreArrayDOF() once you have finished accessing the vector entries.
273: In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]!
275: In Fortran 90 you do not need a version of DMDAVecRestoreArrayDOF() just use DMDAVecRestoreArrayF90() and declare your array with one higher dimension,
276: see src/dm/tutorials/ex11f90.F
278: Level: intermediate
280: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecGetArray(), DMDAVecRestoreArrayDOF(),
281: DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead(), DMStagVecGetArrayDOF(), DMDAVecGetArrayDOFRead()
282: @*/
283: PetscErrorCode DMDAVecGetArrayDOF(DM da,Vec vec,void *array)
284: {
285: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
287: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
288: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
289: DMDAGetInfo(da,&dim,NULL,NULL,NULL,NULL,NULL,NULL,&dof,NULL,NULL,NULL,NULL,NULL);
291: /* Handle case where user passes in global vector as opposed to local */
292: VecGetLocalSize(vec,&N);
293: if (N == xm*ym*zm*dof) {
294: gxm = xm;
295: gym = ym;
296: gzm = zm;
297: gxs = xs;
298: gys = ys;
299: gzs = zs;
302: if (dim == 1) {
303: VecGetArray2d(vec,gxm,dof,gxs,0,(PetscScalar***)array);
304: } else if (dim == 2) {
305: VecGetArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
306: } else if (dim == 3) {
307: VecGetArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
308: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D",dim);
309: return 0;
310: }
312: /*@
313: DMDAVecRestoreArrayDOF - Restores a multiple dimension array obtained with DMDAVecGetArrayDOF()
315: Logically collective
317: Input Parameters:
318: + da - the distributed array
319: . vec - the vector, either a vector the same size as one obtained with
320: DMCreateGlobalVector() or DMCreateLocalVector()
321: - array - the array
323: Level: intermediate
325: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(), DMDAVecRestoreArrayDOF(),
326: DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead(), DMStagVecRestoreArrayDOF(), DMDAVecRestoreArrayDOFRead()
327: @*/
328: PetscErrorCode DMDAVecRestoreArrayDOF(DM da,Vec vec,void *array)
329: {
330: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
332: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
333: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
334: DMDAGetInfo(da,&dim,NULL,NULL,NULL,NULL,NULL,NULL,&dof,NULL,NULL,NULL,NULL,NULL);
336: /* Handle case where user passes in global vector as opposed to local */
337: VecGetLocalSize(vec,&N);
338: if (N == xm*ym*zm*dof) {
339: gxm = xm;
340: gym = ym;
341: gzm = zm;
342: gxs = xs;
343: gys = ys;
344: gzs = zs;
345: }
347: if (dim == 1) {
348: VecRestoreArray2d(vec,gxm,dof,gxs,0,(PetscScalar***)array);
349: } else if (dim == 2) {
350: VecRestoreArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
351: } else if (dim == 3) {
352: VecRestoreArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
353: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D",dim);
354: return 0;
355: }
357: /*@C
358: DMDAVecGetArrayRead - Returns a multiple dimension array that shares data with
359: the underlying vector and is indexed using the global dimensions.
361: Not collective
363: Input Parameters:
364: + da - the distributed array
365: - vec - the vector, either a vector the same size as one obtained with DMCreateGlobalVector() or DMCreateLocalVector()
367: Output Parameter:
368: . array - the array
370: Notes:
371: Call DMDAVecRestoreArrayRead() once you have finished accessing the vector entries.
373: In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]!
375: If vec is a local vector (obtained with DMCreateLocalVector() etc) then the ghost point locations are accessible. If it is
376: a global vector then the ghost points are not accessible. Of course with the local vector you will have had to do the
378: appropriate DMGlobalToLocalBegin() and DMGlobalToLocalEnd() to have correct values in the ghost locations.
380: Fortran Notes:
381: From Fortran use DMDAVecGetArrayReadF90() and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate
382: 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
383: 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
384: array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from
385: DMDAGetCorners() for a global array or DMDAGetGhostCorners() for a local array. Include petsc/finclude/petscdmda.h90 to access this routine.
387: Due to bugs in the compiler DMDAVecGetArrayReadF90() does not work with gfortran versions before 4.5
389: Level: intermediate
391: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArrayRead(), DMDAVecRestoreArrayDOF()
392: DMDAVecGetArrayDOF(), DMDAVecGetArray(), DMDAVecRestoreArray(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead(),
393: DMStagVecGetArrayRead()
394: @*/
395: PetscErrorCode DMDAVecGetArrayRead(DM da,Vec vec,void *array)
396: {
397: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
402: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
403: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
404: DMDAGetInfo(da,&dim,NULL,NULL,NULL,NULL,NULL,NULL,&dof,NULL,NULL,NULL,NULL,NULL);
406: /* Handle case where user passes in global vector as opposed to local */
407: VecGetLocalSize(vec,&N);
408: if (N == xm*ym*zm*dof) {
409: gxm = xm;
410: gym = ym;
411: gzm = zm;
412: gxs = xs;
413: gys = ys;
414: gzs = zs;
417: if (dim == 1) {
418: VecGetArray1dRead(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
419: } else if (dim == 2) {
420: VecGetArray2dRead(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
421: } else if (dim == 3) {
422: VecGetArray3dRead(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
423: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D",dim);
424: return 0;
425: }
427: /*@
428: DMDAVecRestoreArrayRead - Restores a multiple dimension array obtained with DMDAVecGetArrayRead()
430: Not collective
432: Input Parameters:
433: + da - the distributed array
434: . vec - the vector, either a vector the same size as one obtained with
435: DMCreateGlobalVector() or DMCreateLocalVector()
436: - array - the array, non-NULL pointer is zeroed
438: Level: intermediate
440: Fortran Notes:
441: From Fortran use DMDAVecRestoreArrayReadF90()
443: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArrayRead(),
444: DMDAVecGetArray(), DMDAVecRestoreArray(), DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(),
445: DMStagVecRestoreArrayRead()
446: @*/
447: PetscErrorCode DMDAVecRestoreArrayRead(DM da,Vec vec,void *array)
448: {
449: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
454: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
455: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
456: DMDAGetInfo(da,&dim,NULL,NULL,NULL,NULL,NULL,NULL,&dof,NULL,NULL,NULL,NULL,NULL);
458: /* Handle case where user passes in global vector as opposed to local */
459: VecGetLocalSize(vec,&N);
460: if (N == xm*ym*zm*dof) {
461: gxm = xm;
462: gym = ym;
463: gzm = zm;
464: gxs = xs;
465: gys = ys;
466: gzs = zs;
469: if (dim == 1) {
470: VecRestoreArray1dRead(vec,gxm*dof,gxs*dof,(PetscScalar**)array);
471: } else if (dim == 2) {
472: VecRestoreArray2dRead(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);
473: } else if (dim == 3) {
474: VecRestoreArray3dRead(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);
475: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D",dim);
476: return 0;
477: }
479: /*@C
480: DMDAVecGetArrayDOFRead - Returns a multiple dimension array that shares data with
481: the underlying vector and is indexed using the global dimensions.
483: Not Collective
485: Input Parameters:
486: + da - the distributed array
487: - vec - the vector, either a vector the same size as one obtained with
488: DMCreateGlobalVector() or DMCreateLocalVector()
490: Output Parameter:
491: . array - the array
493: Notes:
494: Call DMDAVecRestoreArrayDOFRead() once you have finished accessing the vector entries.
496: In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]!
498: In Fortran 90 you do not need a version of DMDAVecRestoreArrayDOF() just use DMDAVecRestoreArrayReadF90() and declare your array with one higher dimension,
499: see src/dm/tutorials/ex11f90.F
501: Level: intermediate
503: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(),
504: DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead(), DMStagVecGetArrayDOFRead()
505: @*/
506: PetscErrorCode DMDAVecGetArrayDOFRead(DM da,Vec vec,void *array)
507: {
508: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
510: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
511: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
512: DMDAGetInfo(da,&dim,NULL,NULL,NULL,NULL,NULL,NULL,&dof,NULL,NULL,NULL,NULL,NULL);
514: /* Handle case where user passes in global vector as opposed to local */
515: VecGetLocalSize(vec,&N);
516: if (N == xm*ym*zm*dof) {
517: gxm = xm;
518: gym = ym;
519: gzm = zm;
520: gxs = xs;
521: gys = ys;
522: gzs = zs;
525: if (dim == 1) {
526: VecGetArray2dRead(vec,gxm,dof,gxs,0,(PetscScalar***)array);
527: } else if (dim == 2) {
528: VecGetArray3dRead(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
529: } else if (dim == 3) {
530: VecGetArray4dRead(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
531: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D",dim);
532: return 0;
533: }
535: /*@
536: DMDAVecRestoreArrayDOFRead - Restores a multiple dimension array obtained with DMDAVecGetArrayDOFRead()
538: Not Collective
540: Input Parameters:
541: + da - the distributed array
542: . vec - the vector, either a vector the same size as one obtained with
543: DMCreateGlobalVector() or DMCreateLocalVector()
544: - array - the array
546: Level: intermediate
548: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(), DMDAVecRestoreArrayDOF(),
549: DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead(), DMStagVecRestoreArrayDOFRead()
550: @*/
551: PetscErrorCode DMDAVecRestoreArrayDOFRead(DM da,Vec vec,void *array)
552: {
553: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
555: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
556: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
557: DMDAGetInfo(da,&dim,NULL,NULL,NULL,NULL,NULL,NULL,&dof,NULL,NULL,NULL,NULL,NULL);
559: /* Handle case where user passes in global vector as opposed to local */
560: VecGetLocalSize(vec,&N);
561: if (N == xm*ym*zm*dof) {
562: gxm = xm;
563: gym = ym;
564: gzm = zm;
565: gxs = xs;
566: gys = ys;
567: gzs = zs;
568: }
570: if (dim == 1) {
571: VecRestoreArray2dRead(vec,gxm,dof,gxs,0,(PetscScalar***)array);
572: } else if (dim == 2) {
573: VecRestoreArray3dRead(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
574: } else if (dim == 3) {
575: VecRestoreArray4dRead(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
576: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D",dim);
577: return 0;
578: }
580: /*@C
581: DMDAVecGetArrayDOFWrite - Returns a multiple dimension array that shares data with
582: the underlying vector and is indexed using the global dimensions.
584: Not Collective
586: Input Parameters:
587: + da - the distributed array
588: - vec - the vector, either a vector the same size as one obtained with
589: DMCreateGlobalVector() or DMCreateLocalVector()
591: Output Parameter:
592: . array - the array
594: Notes:
595: Call DMDAVecRestoreArrayDOFWrite() once you have finished accessing the vector entries.
597: In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]!
599: In Fortran 90 you do not need a version of DMDAVecRestoreArrayDOF() just use DMDAVecRestoreArrayWriteF90() and declare your array with one higher dimension,
600: see src/dm/tutorials/ex11f90.F
602: Level: intermediate
604: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(),
605: DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMStagVecGetArrayDOFWrite(),
606: DMStagVecRestoreArrayDOFRead(), DMStagVecRestoreArrayDOFRead()
607: @*/
608: PetscErrorCode DMDAVecGetArrayDOFWrite(DM da,Vec vec,void *array)
609: {
610: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
612: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
613: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
614: DMDAGetInfo(da,&dim,NULL,NULL,NULL,NULL,NULL,NULL,&dof,NULL,NULL,NULL,NULL,NULL);
616: /* Handle case where user passes in global vector as opposed to local */
617: VecGetLocalSize(vec,&N);
618: if (N == xm*ym*zm*dof) {
619: gxm = xm;
620: gym = ym;
621: gzm = zm;
622: gxs = xs;
623: gys = ys;
624: gzs = zs;
627: if (dim == 1) {
628: VecGetArray2dWrite(vec,gxm,dof,gxs,0,(PetscScalar***)array);
629: } else if (dim == 2) {
630: VecGetArray3dWrite(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
631: } else if (dim == 3) {
632: VecGetArray4dWrite(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
633: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D",dim);
634: return 0;
635: }
637: /*@
638: DMDAVecRestoreArrayDOFWrite - Restores a multiple dimension array obtained with DMDAVecGetArrayDOFWrite()
640: Not Collective
642: Input Parameters:
643: + da - the distributed array
644: . vec - the vector, either a vector the same size as one obtained with
645: DMCreateGlobalVector() or DMCreateLocalVector()
646: - array - the array
648: Level: intermediate
650: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(), DMDAVecRestoreArrayDOF(),
651: DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMStagVecRestoreArrayDOFWrite(),
652: DMStagVecRestoreArrayDOFRead(), DMStagVecRestoreArrayDOFRead()
653: @*/
654: PetscErrorCode DMDAVecRestoreArrayDOFWrite(DM da,Vec vec,void *array)
655: {
656: PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;
658: DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
659: DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
660: DMDAGetInfo(da,&dim,NULL,NULL,NULL,NULL,NULL,NULL,&dof,NULL,NULL,NULL,NULL,NULL);
662: /* Handle case where user passes in global vector as opposed to local */
663: VecGetLocalSize(vec,&N);
664: if (N == xm*ym*zm*dof) {
665: gxm = xm;
666: gym = ym;
667: gzm = zm;
668: gxs = xs;
669: gys = ys;
670: gzs = zs;
671: }
673: if (dim == 1) {
674: VecRestoreArray2dWrite(vec,gxm,dof,gxs,0,(PetscScalar***)array);
675: } else if (dim == 2) {
676: VecRestoreArray3dWrite(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
677: } else if (dim == 3) {
678: VecRestoreArray4dWrite(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
679: } else SETERRQ(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D",dim);
680: return 0;
681: }