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