1 2 #include <petsc/private/dmdaimpl.h> /*I "petscdmda.h" I*/ 3 4 /*@C 5 DMDAVecGetArray - Returns a multiple dimension array that shares data with 6 the underlying vector and is indexed using the global dimensions. 7 8 Logically collective on da 9 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() 13 14 Output Parameter: 15 . array - the array 16 17 Notes: 18 Call DMDAVecRestoreArray() once you have finished accessing the vector entries. 19 20 In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]! 21 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 24 25 appropriate DMGlobalToLocalBegin() and DMGlobalToLocalEnd() to have correct values in the ghost locations. 26 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. 33 34 Due to bugs in the compiler DMDAVecGetArrayF90() does not work with gfortran versions before 4.5 35 36 Level: intermediate 37 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 { 43 PetscErrorCode ierr; 44 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; 45 46 PetscFunctionBegin; 47 PetscValidHeaderSpecificType(da, DM_CLASSID, 1,DMDA); 48 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 49 PetscValidPointer(array, 3); 50 if (da->defaultSection) { 51 ierr = VecGetArray(vec,(PetscScalar**)array);CHKERRQ(ierr); 52 PetscFunctionReturn(0); 53 } 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 PetscFunctionReturn(0); 77 } 78 79 /*@ 80 DMDAVecRestoreArray - Restores a multiple dimension array obtained with DMDAVecGetArray() 81 82 Logically collective on da 83 84 Input Parameter: 85 + da - the distributed array 86 . vec - the vector, either a vector the same size as one obtained with 87 DMCreateGlobalVector() or DMCreateLocalVector() 88 - array - the array, non-NULL pointer is zeroed 89 90 Level: intermediate 91 92 Fortran Notes: 93 From Fortran use DMDAVecRestoreArayF90() 94 95 .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(), 96 DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead() 97 @*/ 98 PetscErrorCode DMDAVecRestoreArray(DM da,Vec vec,void *array) 99 { 100 PetscErrorCode ierr; 101 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; 102 103 PetscFunctionBegin; 104 PetscValidHeaderSpecificType(da, DM_CLASSID, 1,DMDA); 105 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 106 PetscValidPointer(array, 3); 107 if (da->defaultSection) { 108 ierr = VecRestoreArray(vec,(PetscScalar**)array);CHKERRQ(ierr); 109 PetscFunctionReturn(0); 110 } 111 ierr = DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);CHKERRQ(ierr); 112 ierr = DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);CHKERRQ(ierr); 113 ierr = DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);CHKERRQ(ierr); 114 115 /* Handle case where user passes in global vector as opposed to local */ 116 ierr = VecGetLocalSize(vec,&N);CHKERRQ(ierr); 117 if (N == xm*ym*zm*dof) { 118 gxm = xm; 119 gym = ym; 120 gzm = zm; 121 gxs = xs; 122 gys = ys; 123 gzs = zs; 124 } 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); 125 126 if (dim == 1) { 127 ierr = VecRestoreArray1d(vec,gxm*dof,gxs*dof,(PetscScalar**)array);CHKERRQ(ierr); 128 } else if (dim == 2) { 129 ierr = VecRestoreArray2d(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);CHKERRQ(ierr); 130 } else if (dim == 3) { 131 ierr = VecRestoreArray3d(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);CHKERRQ(ierr); 132 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim); 133 PetscFunctionReturn(0); 134 } 135 136 /*@C 137 DMDAVecGetArrayWrite - Returns a multiple dimension array that shares data with 138 the underlying vector and is indexed using the global dimensions. 139 140 Logically collective on Vec 141 142 Input Parameter: 143 + da - the distributed array 144 - vec - the vector, either a vector the same size as one obtained with DMCreateGlobalVector() or DMCreateLocalVector() 145 146 Output Parameter: 147 . array - the array 148 149 Notes: 150 Call DMDAVecRestoreArray() once you have finished accessing the vector entries. 151 152 In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]! 153 154 If vec is a local vector (obtained with DMCreateLocalVector() etc) then the ghost point locations are accessible. If it is 155 a global vector then the ghost points are not accessible. Of course with the local vector you will have had to do the 156 157 appropriate DMGlobalToLocalBegin() and DMGlobalToLocalEnd() to have correct values in the ghost locations. 158 159 Fortran Notes: 160 From Fortran use DMDAVecGetArrayF90() and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate 161 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 162 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 163 array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from 164 DMDAGetCorners() for a global array or DMDAGetGhostCorners() for a local array. Include petsc/finclude/petscdmda.h90 to access this routine. 165 166 Due to bugs in the compiler DMDAVecGetArrayF90() does not work with gfortran versions before 4.5 167 168 Level: intermediate 169 170 Developer Notes: This has code duplication with DMDAVecGetArray() and DMDAVecGetArrayRead() 171 172 .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArrayWrite(), DMDAVecRestoreArrayDOF() 173 DMDAVecGetArrayDOF(), DMDAVecGetArray(), DMDAVecRestoreArray(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead() 174 @*/ 175 PetscErrorCode DMDAVecGetArrayWrite(DM da,Vec vec,void *array) 176 { 177 PetscErrorCode ierr; 178 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; 179 180 PetscFunctionBegin; 181 PetscValidHeaderSpecificType(da, DM_CLASSID, 1,DMDA); 182 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 183 PetscValidPointer(array, 3); 184 if (da->defaultSection) { 185 ierr = VecGetArrayWrite(vec,(PetscScalar**)array);CHKERRQ(ierr); 186 PetscFunctionReturn(0); 187 } 188 ierr = DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);CHKERRQ(ierr); 189 ierr = DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);CHKERRQ(ierr); 190 ierr = DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);CHKERRQ(ierr); 191 192 /* Handle case where user passes in global vector as opposed to local */ 193 ierr = VecGetLocalSize(vec,&N);CHKERRQ(ierr); 194 if (N == xm*ym*zm*dof) { 195 gxm = xm; 196 gym = ym; 197 gzm = zm; 198 gxs = xs; 199 gys = ys; 200 gzs = zs; 201 } 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); 202 203 if (dim == 1) { 204 ierr = VecGetArray1dWrite(vec,gxm*dof,gxs*dof,(PetscScalar**)array);CHKERRQ(ierr); 205 } else if (dim == 2) { 206 ierr = VecGetArray2dWrite(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);CHKERRQ(ierr); 207 } else if (dim == 3) { 208 ierr = VecGetArray3dWrite(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);CHKERRQ(ierr); 209 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim); 210 PetscFunctionReturn(0); 211 } 212 213 /*@ 214 DMDAVecRestoreArrayWrite - Restores a multiple dimension array obtained with DMDAVecGetArrayWrite() 215 216 Logically collective on Vec 217 218 Input Parameter: 219 + da - the distributed array 220 . vec - the vector, either a vector the same size as one obtained with 221 DMCreateGlobalVector() or DMCreateLocalVector() 222 - array - the array, non-NULL pointer is zeroed 223 224 Level: intermediate 225 226 Fortran Notes: 227 From Fortran use DMDAVecRestoreArayF90() 228 229 .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArrayWrite(), 230 DMDAVecGetArray(), DMDAVecRestoreArray(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead() 231 @*/ 232 PetscErrorCode DMDAVecRestoreArrayWrite(DM da,Vec vec,void *array) 233 { 234 PetscErrorCode ierr; 235 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; 236 237 PetscFunctionBegin; 238 PetscValidHeaderSpecificType(da, DM_CLASSID, 1,DMDA); 239 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 240 PetscValidPointer(array, 3); 241 if (da->defaultSection) { 242 ierr = VecRestoreArray(vec,(PetscScalar**)array);CHKERRQ(ierr); 243 PetscFunctionReturn(0); 244 } 245 ierr = DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);CHKERRQ(ierr); 246 ierr = DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);CHKERRQ(ierr); 247 ierr = DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);CHKERRQ(ierr); 248 249 /* Handle case where user passes in global vector as opposed to local */ 250 ierr = VecGetLocalSize(vec,&N);CHKERRQ(ierr); 251 if (N == xm*ym*zm*dof) { 252 gxm = xm; 253 gym = ym; 254 gzm = zm; 255 gxs = xs; 256 gys = ys; 257 gzs = zs; 258 } 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); 259 260 if (dim == 1) { 261 ierr = VecRestoreArray1dWrite(vec,gxm*dof,gxs*dof,(PetscScalar**)array);CHKERRQ(ierr); 262 } else if (dim == 2) { 263 ierr = VecRestoreArray2dWrite(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);CHKERRQ(ierr); 264 } else if (dim == 3) { 265 ierr = VecRestoreArray3dWrite(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);CHKERRQ(ierr); 266 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim); 267 PetscFunctionReturn(0); 268 } 269 270 /*@C 271 DMDAVecGetArrayDOF - Returns a multiple dimension array that shares data with 272 the underlying vector and is indexed using the global dimensions. 273 274 Logically collective 275 276 Input Parameter: 277 + da - the distributed array 278 - vec - the vector, either a vector the same size as one obtained with 279 DMCreateGlobalVector() or DMCreateLocalVector() 280 281 Output Parameter: 282 . array - the array 283 284 Notes: 285 Call DMDAVecRestoreArrayDOF() once you have finished accessing the vector entries. 286 287 In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]! 288 289 In Fortran 90 you do not need a version of DMDAVecRestoreArrayDOF() just use DMDAVecRestoreArrayF90() and declare your array with one higher dimension, 290 see src/dm/examples/tutorials/ex11f90.F 291 292 Level: intermediate 293 294 .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecGetArray(), DMDAVecRestoreArrayDOF(), 295 DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead() 296 @*/ 297 PetscErrorCode DMDAVecGetArrayDOF(DM da,Vec vec,void *array) 298 { 299 PetscErrorCode ierr; 300 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; 301 302 PetscFunctionBegin; 303 ierr = DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);CHKERRQ(ierr); 304 ierr = DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);CHKERRQ(ierr); 305 ierr = DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);CHKERRQ(ierr); 306 307 /* Handle case where user passes in global vector as opposed to local */ 308 ierr = VecGetLocalSize(vec,&N);CHKERRQ(ierr); 309 if (N == xm*ym*zm*dof) { 310 gxm = xm; 311 gym = ym; 312 gzm = zm; 313 gxs = xs; 314 gys = ys; 315 gzs = zs; 316 } 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); 317 318 if (dim == 1) { 319 ierr = VecGetArray2d(vec,gxm,dof,gxs,0,(PetscScalar***)array);CHKERRQ(ierr); 320 } else if (dim == 2) { 321 ierr = VecGetArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);CHKERRQ(ierr); 322 } else if (dim == 3) { 323 ierr = VecGetArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);CHKERRQ(ierr); 324 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim); 325 PetscFunctionReturn(0); 326 } 327 328 /*@ 329 DMDAVecRestoreArrayDOF - Restores a multiple dimension array obtained with DMDAVecGetArrayDOF() 330 331 Logically collective 332 333 Input Parameter: 334 + da - the distributed array 335 . vec - the vector, either a vector the same size as one obtained with 336 DMCreateGlobalVector() or DMCreateLocalVector() 337 - array - the array 338 339 Level: intermediate 340 341 .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(), DMDAVecRestoreArrayDOF(), 342 DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead() 343 @*/ 344 PetscErrorCode DMDAVecRestoreArrayDOF(DM da,Vec vec,void *array) 345 { 346 PetscErrorCode ierr; 347 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; 348 349 PetscFunctionBegin; 350 ierr = DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);CHKERRQ(ierr); 351 ierr = DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);CHKERRQ(ierr); 352 ierr = DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);CHKERRQ(ierr); 353 354 /* Handle case where user passes in global vector as opposed to local */ 355 ierr = VecGetLocalSize(vec,&N);CHKERRQ(ierr); 356 if (N == xm*ym*zm*dof) { 357 gxm = xm; 358 gym = ym; 359 gzm = zm; 360 gxs = xs; 361 gys = ys; 362 gzs = zs; 363 } 364 365 if (dim == 1) { 366 ierr = VecRestoreArray2d(vec,gxm,dof,gxs,0,(PetscScalar***)array);CHKERRQ(ierr); 367 } else if (dim == 2) { 368 ierr = VecRestoreArray3d(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);CHKERRQ(ierr); 369 } else if (dim == 3) { 370 ierr = VecRestoreArray4d(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);CHKERRQ(ierr); 371 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim); 372 PetscFunctionReturn(0); 373 } 374 375 /*@C 376 DMDAVecGetArrayRead - Returns a multiple dimension array that shares data with 377 the underlying vector and is indexed using the global dimensions. 378 379 Not collective 380 381 Input Parameter: 382 + da - the distributed array 383 - vec - the vector, either a vector the same size as one obtained with DMCreateGlobalVector() or DMCreateLocalVector() 384 385 Output Parameter: 386 . array - the array 387 388 Notes: 389 Call DMDAVecRestoreArrayRead() once you have finished accessing the vector entries. 390 391 In C, the indexing is "backwards" from what expects: array[k][j][i] NOT array[i][j][k]! 392 393 If vec is a local vector (obtained with DMCreateLocalVector() etc) then the ghost point locations are accessible. If it is 394 a global vector then the ghost points are not accessible. Of course with the local vector you will have had to do the 395 396 appropriate DMGlobalToLocalBegin() and DMGlobalToLocalEnd() to have correct values in the ghost locations. 397 398 Fortran Notes: 399 From Fortran use DMDAVecGetArrayReadF90() and pass for the array type PetscScalar,pointer :: array(:,...,:) of the appropriate 400 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 401 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 402 array(0:dof-1,xs:xs+xm-1,ys:ys+ym-1,zs:zs+zm-1) where the values are obtained from 403 DMDAGetCorners() for a global array or DMDAGetGhostCorners() for a local array. Include petsc/finclude/petscdmda.h90 to access this routine. 404 405 Due to bugs in the compiler DMDAVecGetArrayReadF90() does not work with gfortran versions before 4.5 406 407 Level: intermediate 408 409 .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArrayRead(), DMDAVecRestoreArrayDOF() 410 DMDAVecGetArrayDOF(), DMDAVecGetArray(), DMDAVecRestoreArray(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead() 411 @*/ 412 PetscErrorCode DMDAVecGetArrayRead(DM da,Vec vec,void *array) 413 { 414 PetscErrorCode ierr; 415 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; 416 417 PetscFunctionBegin; 418 PetscValidHeaderSpecificType(da, DM_CLASSID, 1,DMDA); 419 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 420 PetscValidPointer(array, 3); 421 if (da->defaultSection) { 422 ierr = VecGetArrayRead(vec,(const PetscScalar**)array);CHKERRQ(ierr); 423 PetscFunctionReturn(0); 424 } 425 ierr = DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);CHKERRQ(ierr); 426 ierr = DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);CHKERRQ(ierr); 427 ierr = DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);CHKERRQ(ierr); 428 429 /* Handle case where user passes in global vector as opposed to local */ 430 ierr = VecGetLocalSize(vec,&N);CHKERRQ(ierr); 431 if (N == xm*ym*zm*dof) { 432 gxm = xm; 433 gym = ym; 434 gzm = zm; 435 gxs = xs; 436 gys = ys; 437 gzs = zs; 438 } 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); 439 440 if (dim == 1) { 441 ierr = VecGetArray1dRead(vec,gxm*dof,gxs*dof,(PetscScalar**)array);CHKERRQ(ierr); 442 } else if (dim == 2) { 443 ierr = VecGetArray2dRead(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);CHKERRQ(ierr); 444 } else if (dim == 3) { 445 ierr = VecGetArray3dRead(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);CHKERRQ(ierr); 446 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim); 447 PetscFunctionReturn(0); 448 } 449 450 /*@ 451 DMDAVecRestoreArrayRead - Restores a multiple dimension array obtained with DMDAVecGetArrayRead() 452 453 Not collective 454 455 Input Parameter: 456 + da - the distributed array 457 . vec - the vector, either a vector the same size as one obtained with 458 DMCreateGlobalVector() or DMCreateLocalVector() 459 - array - the array, non-NULL pointer is zeroed 460 461 Level: intermediate 462 463 Fortran Notes: 464 From Fortran use DMDAVecRestoreArrayReadF90() 465 466 .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArrayRead(), 467 DMDAVecGetArray(), DMDAVecRestoreArray(), DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite() 468 @*/ 469 PetscErrorCode DMDAVecRestoreArrayRead(DM da,Vec vec,void *array) 470 { 471 PetscErrorCode ierr; 472 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; 473 474 PetscFunctionBegin; 475 PetscValidHeaderSpecificType(da, DM_CLASSID, 1,DMDA); 476 PetscValidHeaderSpecific(vec, VEC_CLASSID, 2); 477 PetscValidPointer(array, 3); 478 if (da->defaultSection) { 479 ierr = VecRestoreArrayRead(vec,(const PetscScalar**)array);CHKERRQ(ierr); 480 PetscFunctionReturn(0); 481 } 482 ierr = DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);CHKERRQ(ierr); 483 ierr = DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);CHKERRQ(ierr); 484 ierr = DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);CHKERRQ(ierr); 485 486 /* Handle case where user passes in global vector as opposed to local */ 487 ierr = VecGetLocalSize(vec,&N);CHKERRQ(ierr); 488 if (N == xm*ym*zm*dof) { 489 gxm = xm; 490 gym = ym; 491 gzm = zm; 492 gxs = xs; 493 gys = ys; 494 gzs = zs; 495 } 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); 496 497 if (dim == 1) { 498 ierr = VecRestoreArray1dRead(vec,gxm*dof,gxs*dof,(PetscScalar**)array);CHKERRQ(ierr); 499 } else if (dim == 2) { 500 ierr = VecRestoreArray2dRead(vec,gym,gxm*dof,gys,gxs*dof,(PetscScalar***)array);CHKERRQ(ierr); 501 } else if (dim == 3) { 502 ierr = VecRestoreArray3dRead(vec,gzm,gym,gxm*dof,gzs,gys,gxs*dof,(PetscScalar****)array);CHKERRQ(ierr); 503 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim); 504 PetscFunctionReturn(0); 505 } 506 507 /*@C 508 DMDAVecGetArrayDOFRead - Returns a multiple dimension array that shares data with 509 the underlying vector and is indexed using the global dimensions. 510 511 Not Collective 512 513 Input Parameter: 514 + da - the distributed array 515 - vec - the vector, either a vector the same size as one obtained with 516 DMCreateGlobalVector() or DMCreateLocalVector() 517 518 Output Parameter: 519 . array - the array 520 521 Notes: 522 Call DMDAVecRestoreArrayDOFRead() once you have finished accessing the vector entries. 523 524 In C, the indexing is "backwards" from what expects: array[k][j][i][DOF] NOT array[i][j][k][DOF]! 525 526 In Fortran 90 you do not need a version of DMDAVecRestoreArrayDOF() just use DMDAVecRestoreArrayReadF90() and declare your array with one higher dimension, 527 see src/dm/examples/tutorials/ex11f90.F 528 529 Level: intermediate 530 531 .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecRestoreArray(), DMDAVecGetArray(), DMDAVecRestoreArrayDOF(), 532 DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead() 533 @*/ 534 PetscErrorCode DMDAVecGetArrayDOFRead(DM da,Vec vec,void *array) 535 { 536 PetscErrorCode ierr; 537 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; 538 539 PetscFunctionBegin; 540 ierr = DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);CHKERRQ(ierr); 541 ierr = DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);CHKERRQ(ierr); 542 ierr = DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);CHKERRQ(ierr); 543 544 /* Handle case where user passes in global vector as opposed to local */ 545 ierr = VecGetLocalSize(vec,&N);CHKERRQ(ierr); 546 if (N == xm*ym*zm*dof) { 547 gxm = xm; 548 gym = ym; 549 gzm = zm; 550 gxs = xs; 551 gys = ys; 552 gzs = zs; 553 } 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); 554 555 if (dim == 1) { 556 ierr = VecGetArray2dRead(vec,gxm,dof,gxs,0,(PetscScalar***)array);CHKERRQ(ierr); 557 } else if (dim == 2) { 558 ierr = VecGetArray3dRead(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);CHKERRQ(ierr); 559 } else if (dim == 3) { 560 ierr = VecGetArray4dRead(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);CHKERRQ(ierr); 561 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim); 562 PetscFunctionReturn(0); 563 } 564 565 /*@ 566 DMDAVecRestoreArrayDOFRead - Restores a multiple dimension array obtained with DMDAVecGetArrayDOFRead() 567 568 Not Collective 569 570 Input Parameter: 571 + da - the distributed array 572 . vec - the vector, either a vector the same size as one obtained with 573 DMCreateGlobalVector() or DMCreateLocalVector() 574 - array - the array 575 576 Level: intermediate 577 578 .seealso: DMDAGetGhostCorners(), DMDAGetCorners(), VecGetArray(), VecRestoreArray(), DMDAVecGetArray(), DMDAVecGetArrayDOF(), DMDAVecRestoreArrayDOF(), 579 DMDAVecGetArrayWrite(), DMDAVecRestoreArrayWrite(), DMDAVecGetArrayRead(), DMDAVecRestoreArrayRead() 580 @*/ 581 PetscErrorCode DMDAVecRestoreArrayDOFRead(DM da,Vec vec,void *array) 582 { 583 PetscErrorCode ierr; 584 PetscInt xs,ys,zs,xm,ym,zm,gxs,gys,gzs,gxm,gym,gzm,N,dim,dof; 585 586 PetscFunctionBegin; 587 ierr = DMDAGetCorners(da,&xs,&ys,&zs,&xm,&ym,&zm);CHKERRQ(ierr); 588 ierr = DMDAGetGhostCorners(da,&gxs,&gys,&gzs,&gxm,&gym,&gzm);CHKERRQ(ierr); 589 ierr = DMDAGetInfo(da,&dim,0,0,0,0,0,0,&dof,0,0,0,0,0);CHKERRQ(ierr); 590 591 /* Handle case where user passes in global vector as opposed to local */ 592 ierr = VecGetLocalSize(vec,&N);CHKERRQ(ierr); 593 if (N == xm*ym*zm*dof) { 594 gxm = xm; 595 gym = ym; 596 gzm = zm; 597 gxs = xs; 598 gys = ys; 599 gzs = zs; 600 } 601 602 if (dim == 1) { 603 ierr = VecRestoreArray2dRead(vec,gxm,dof,gxs,0,(PetscScalar***)array);CHKERRQ(ierr); 604 } else if (dim == 2) { 605 ierr = VecRestoreArray3dRead(vec,gym,gxm,dof,gys,gxs,0,(PetscScalar****)array);CHKERRQ(ierr); 606 } else if (dim == 3) { 607 ierr = VecRestoreArray4dRead(vec,gzm,gym,gxm,dof,gzs,gys,gxs,0,(PetscScalar*****)array);CHKERRQ(ierr); 608 } else SETERRQ1(PETSC_COMM_SELF,PETSC_ERR_ARG_CORRUPT,"DMDA dimension not 1, 2, or 3, it is %D\n",dim); 609 PetscFunctionReturn(0); 610 } 611 612 613 614 615 616 617 618 619 620 621 622 623 624