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