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