Actual source code: dagetarray.c

  1: 
  2: #include <petscdmda.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:    Not Collective

 12:    Input Parameter:
 13: +  da - the distributed array
 14: -  vec - the vector, either a vector the same size as one obtained with 
 15:          DMCreateGlobalVector() or DMCreateLocalVector()
 16:    
 17:    Output Parameter:
 18: .  array - the array

 20:    Notes:
 21:     Call DMDAVecRestoreArray() once you have finished accessing the vector entries.

 23:     In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]!

 25:     If vec is a local vector (obtained with DMCreateLocalVector() etc) then they ghost point locations are accessable. If it is 
 26:     a global vector then the ghost points are not accessable. Of course with the local vector you will have had to do the 

 28:     appropriate DMLocalToGlobalBegin() and DMLocalToGlobalEnd() to have correct values in the ghost locations.

 30:   Fortran Notes: From Fortran use DMDAVecGetArrayF90() and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate 
 31:        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 
 32:        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
 33:        array(1:dof,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from 
 34:        DMDAGetCorners() for a global array or DMDAGetGhostCorners() for a local array. Include finclude/petscdmda.h90 to access this routine.

 36:   Due to bugs in the compiler DMDAVecGetArrayF90() does not work with gfortran versions before 2.5

 38:   Level: intermediate

 40: .keywords: distributed array, get, corners, nodes, local indices, coordinates

 42: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecRestoreArrayDOF()
 43:           DMDAVecGetArrayDOF()
 44: @*/
 45: PetscErrorCode  DMDAVecGetArray(DM da,Vec vec,void *array)
 46: {
 48:   PetscInt       xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;

 54:   DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
 55:   DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
 56:   DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);

 58:   /* Handle case where user passes in global vector as opposed to local */
 59:   VecGetLocalSize(vec,&N);
 60:   if (N == xm*ym*zm*dof) {
 61:     gxm = xm;
 62:     gym = ym;
 63:     gzm = zm;
 64:     gxs = xs;
 65:     gys = ys;
 66:     gzs = zs;
 67:   } else if (N != gxm*gym*gzm*dof) {
 68:     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);
 69:   }

 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 {
 78:     SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
 79:   }

 81:   return(0);
 82: }

 86: /*@
 87:    DMDAVecRestoreArray - Restores a multiple dimension array obtained with DMDAVecGetArray()

 89:    Not Collective

 91:    Input Parameter:
 92: +  da - the distributed array
 93: .  vec - the vector, either a vector the same size as one obtained with 
 94:          DMCreateGlobalVector() or DMCreateLocalVector()
 95: -  array - the array

 97:   Level: intermediate

 99:   Fortran Notes: From Fortran use DMDAVecRestoreArrayF90()

101: .keywords: distributed array, get, corners, nodes, local indices, coordinates

103: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray()
104: @*/
105: PetscErrorCode  DMDAVecRestoreArray(DM da,Vec vec,void *array)
106: {
108:   PetscInt       xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;

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) {
128:     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:   }

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 {
138:     SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
139:   }
140:   return(0);
141: }

145: /*@C
146:    DMDAVecGetArrayDOF - Returns a multiple dimension array that shares data with 
147:       the underlying vector and is indexed using the global dimensions.

149:    Not Collective

151:    Input Parameter:
152: +  da - the distributed array
153: -  vec - the vector, either a vector the same size as one obtained with 
154:          DMCreateGlobalVector() or DMCreateLocalVector()
155:    
156:    Output Parameter:
157: .  array - the array

159:    Notes:
160:     Call DMDAVecRestoreArrayDOF() once you have finished accessing the vector entries.

162:     In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]!

164:   Level: intermediate

166: .keywords: distributed array, get, corners, nodes, local indices, coordinates

168: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecGetArray(), DMDAVecRestoreArrayDOF()
169: @*/
170: PetscErrorCode  DMDAVecGetArrayDOF(DM da,Vec vec,void *array)
171: {
173:   PetscInt       xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;

176:   DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
177:   DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
178:   DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);

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;
189:   } else if (N != gxm*gym*gzm*dof) {
190:     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);
191:   }

193:   if (dim == 1) {
194:     VecGetArray2d(vec,gxm,dof,gxs,0,(PetscScalar ***)array);
195:   } else if (dim == 2) {
196:     VecGetArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
197:   } else if (dim == 3) {
198:     VecGetArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
199:   } else {
200:     SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
201:   }
202:   return(0);
203: }

207: /*@
208:    DMDAVecRestoreArrayDOF - Restores a multiple dimension array obtained with DMDAVecGetArrayDOF()

210:    Not Collective

212:    Input Parameter:
213: +  da - the distributed array
214: .  vec - the vector, either a vector the same size as one obtained with 
215:          DMCreateGlobalVector() or DMCreateLocalVector()
216: -  array - the array

218:   Level: intermediate

220: .keywords: distributed array, get, corners, nodes, local indices, coordinates

222: .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(), DMDAVecRestoreArrayDOF()
223: @*/
224: PetscErrorCode  DMDAVecRestoreArrayDOF(DM da,Vec vec,void *array)
225: {
227:   PetscInt       xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof;

230:   DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);
231:   DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);
232:   DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);

234:   /* Handle case where user passes in global vector as opposed to local */
235:   VecGetLocalSize(vec,&N);
236:   if (N == xm*ym*zm*dof) {
237:     gxm = xm;
238:     gym = ym;
239:     gzm = zm;
240:     gxs = xs;
241:     gys = ys;
242:     gzs = zs;
243:   }

245:   if (dim == 1) {
246:     VecRestoreArray2d(vec,gxm,dof,gxs,0,(PetscScalar***)array);
247:   } else if (dim == 2) {
248:     VecRestoreArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);
249:   } else if (dim == 3) {
250:     VecRestoreArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);
251:   } else {
252:     SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim);
253:   }
254:   return(0);
255: }