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