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