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