1 2 #include <petsc-private/daimpl.h> /*I "petscdmda.h" I*/ 3 4 #undef __FUNCT__ 5 #define __FUNCT__ "DMDAVecGetArray" 6 /*@C 7 DMDAVecGetArray - Returns a multiple dimension array that shares data with 8 the underlying vector and is indexed using the global dimensions. 9 10 Collective on Vec 11 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 19 20 Notes: 21 Call DMDAVecRestoreArray() once you have finished accessing the vector entries. 22 23 In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]! 24 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 27 28 appropriate DMLocalToGlobalBegin() and DMLocalToGlobalEnd() to have correct values in the ghost locations. 29 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(0:dof-1,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. 35 36 Due to bugs in the compiler DMDAVecGetArrayF90() does not work with gfortran versions before 4.5 37 38 Level: intermediate 39 40 .keywords: distributed array, get, corners, nodes, local indices, coordinates 41 42 .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecRestoreArrayDOF() 43 DMDAVecGetArrayDOF() 44 @*/ 45 PetscErrorCode DMDAVecGetArray(DM da,Vec vec,void *array) 46 { 47 PetscErrorCode ierr; 48 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; 49 50 PetscFunctionBegin; 51 PetscValidHeaderSpecific(da, DM_CLASSID, 1); 52 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 53 PetscValidPointer(array, 3); 54 if (da->defaultSection) { 55 ierr = VecGetArray(vec,(PetscScalar**)array);CHKERRQ(ierr); 56 PetscFunctionReturn(0); 57 } 58 ierr = DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);CHKERRQ(ierr); 59 ierr = DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);CHKERRQ(ierr); 60 ierr = DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);CHKERRQ(ierr); 61 62 /* Handle case where user passes in global vector as opposed to local */ 63 ierr = VecGetLocalSize(vec,&N);CHKERRQ(ierr); 64 if (N == xm*ym*zm*dof) { 65 gxm = xm; 66 gym = ym; 67 gzm = zm; 68 gxs = xs; 69 gys = ys; 70 gzs = zs; 71 } 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 73 if (dim == 1) { 74 ierr = VecGetArray1d(vec,gxm*dof,gxs*dof,(PetscScalar **)array);CHKERRQ(ierr); 75 } else if (dim == 2) { 76 ierr = VecGetArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);CHKERRQ(ierr); 77 } else if (dim == 3) { 78 ierr = VecGetArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);CHKERRQ(ierr); 79 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim); 80 81 PetscFunctionReturn(0); 82 } 83 84 #undef __FUNCT__ 85 #define __FUNCT__ "DMDAVecRestoreArray" 86 /*@ 87 DMDAVecRestoreArray - Restores a multiple dimension array obtained with DMDAVecGetArray() 88 89 Collective on Vec 90 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 96 97 Level: intermediate 98 99 Fortran Notes: From Fortran use DMDAVecRestoreArrayF90() 100 101 .keywords: distributed array, get, corners, nodes, local indices, coordinates 102 103 .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray() 104 @*/ 105 PetscErrorCode DMDAVecRestoreArray(DM da,Vec vec,void *array) 106 { 107 PetscErrorCode ierr; 108 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; 109 110 PetscFunctionBegin; 111 PetscValidHeaderSpecific(da, DM_CLASSID, 1); 112 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 113 PetscValidPointer(array, 3); 114 if (da->defaultSection) { 115 ierr = VecRestoreArray(vec,(PetscScalar**)array);CHKERRQ(ierr); 116 PetscFunctionReturn(0); 117 } 118 ierr = DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);CHKERRQ(ierr); 119 ierr = DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);CHKERRQ(ierr); 120 ierr = DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);CHKERRQ(ierr); 121 122 /* Handle case where user passes in global vector as opposed to local */ 123 ierr = VecGetLocalSize(vec,&N);CHKERRQ(ierr); 124 if (N == xm*ym*zm*dof) { 125 gxm = xm; 126 gym = ym; 127 gzm = zm; 128 gxs = xs; 129 gys = ys; 130 gzs = zs; 131 } 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); 132 133 if (dim == 1) { 134 ierr = VecRestoreArray1d(vec,gxm*dof,gxs*dof,(PetscScalar **)array);CHKERRQ(ierr); 135 } else if (dim == 2) { 136 ierr = VecRestoreArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);CHKERRQ(ierr); 137 } else if (dim == 3) { 138 ierr = VecRestoreArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);CHKERRQ(ierr); 139 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim); 140 PetscFunctionReturn(0); 141 } 142 143 #undef __FUNCT__ 144 #define __FUNCT__ "DMDAVecGetArrayDOF" 145 /*@C 146 DMDAVecGetArrayDOF - Returns a multiple dimension array that shares data with 147 the underlying vector and is indexed using the global dimensions. 148 149 Not Collective 150 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 158 159 Notes: 160 Call DMDAVecRestoreArrayDOF() once you have finished accessing the vector entries. 161 162 In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]! 163 164 In Fortran 90 you do not need a version of DMDAVecRestoreArrayDOF() just use DMDAVecRestoreArrayF90() and declare your array with one higher dimension, 165 see src/dm/examples/tutorials/ex11f90.F 166 167 Level: intermediate 168 169 .keywords: distributed array, get, corners, nodes, local indices, coordinates 170 171 .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecGetArray(), DMDAVecRestoreArrayDOF() 172 @*/ 173 PetscErrorCode DMDAVecGetArrayDOF(DM da,Vec vec,void *array) 174 { 175 PetscErrorCode ierr; 176 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; 177 178 PetscFunctionBegin; 179 ierr = DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);CHKERRQ(ierr); 180 ierr = DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);CHKERRQ(ierr); 181 ierr = DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);CHKERRQ(ierr); 182 183 /* Handle case where user passes in global vector as opposed to local */ 184 ierr = VecGetLocalSize(vec,&N);CHKERRQ(ierr); 185 if (N == xm*ym*zm*dof) { 186 gxm = xm; 187 gym = ym; 188 gzm = zm; 189 gxs = xs; 190 gys = ys; 191 gzs = zs; 192 } 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); 193 194 if (dim == 1) { 195 ierr = VecGetArray2d(vec,gxm,dof,gxs,0,(PetscScalar ***)array);CHKERRQ(ierr); 196 } else if (dim == 2) { 197 ierr = VecGetArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);CHKERRQ(ierr); 198 } else if (dim == 3) { 199 ierr = VecGetArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);CHKERRQ(ierr); 200 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim); 201 PetscFunctionReturn(0); 202 } 203 204 #undef __FUNCT__ 205 #define __FUNCT__ "DMDAVecRestoreArrayDOF" 206 /*@ 207 DMDAVecRestoreArrayDOF - Restores a multiple dimension array obtained with DMDAVecGetArrayDOF() 208 209 Not Collective 210 211 Input Parameter: 212 + da - the distributed array 213 . vec - the vector, either a vector the same size as one obtained with 214 DMCreateGlobalVector() or DMCreateLocalVector() 215 - array - the array 216 217 Level: intermediate 218 219 .keywords: distributed array, get, corners, nodes, local indices, coordinates 220 221 .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(), DMDAVecRestoreArrayDOF() 222 @*/ 223 PetscErrorCode DMDAVecRestoreArrayDOF(DM da,Vec vec,void *array) 224 { 225 PetscErrorCode ierr; 226 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; 227 228 PetscFunctionBegin; 229 ierr = DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);CHKERRQ(ierr); 230 ierr = DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);CHKERRQ(ierr); 231 ierr = DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);CHKERRQ(ierr); 232 233 /* Handle case where user passes in global vector as opposed to local */ 234 ierr = VecGetLocalSize(vec,&N);CHKERRQ(ierr); 235 if (N == xm*ym*zm*dof) { 236 gxm = xm; 237 gym = ym; 238 gzm = zm; 239 gxs = xs; 240 gys = ys; 241 gzs = zs; 242 } 243 244 if (dim == 1) { 245 ierr = VecRestoreArray2d(vec,gxm,dof,gxs,0,(PetscScalar***)array);CHKERRQ(ierr); 246 } else if (dim == 2) { 247 ierr = VecRestoreArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);CHKERRQ(ierr); 248 } else if (dim == 3) { 249 ierr = VecRestoreArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);CHKERRQ(ierr); 250 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim); 251 PetscFunctionReturn(0); 252 } 253 254 255 256 257 258 259 260 261 262 263 264 265 266