1# -------------------------------------------------------------------- 2 3class DMDAStencilType(object): 4 """Stencil types.""" 5 STAR = DMDA_STENCIL_STAR 6 BOX = DMDA_STENCIL_BOX 7 8 9class DMDAInterpolationType(object): 10 """Interpolation types.""" 11 Q0 = DMDA_INTERPOLATION_Q0 12 Q1 = DMDA_INTERPOLATION_Q1 13 14 15class DMDAElementType(object): 16 """Element types.""" 17 P1 = DMDA_ELEMENT_P1 18 Q1 = DMDA_ELEMENT_Q1 19 20# -------------------------------------------------------------------- 21 22 23cdef class DMDA(DM): 24 """A DM object that is used to manage data for a structured grid.""" 25 26 StencilType = DMDAStencilType 27 InterpolationType = DMDAInterpolationType 28 ElementType = DMDAElementType 29 30 # 31 32 def create( 33 self, 34 dim: int | None = None, 35 dof: int | None = None, 36 sizes: DimsSpec | None = None, 37 proc_sizes: DimsSpec | None = None, 38 boundary_type: tuple[DM.BoundaryType | int | str | bool, ...] | None = None, 39 stencil_type: StencilType | None = None, 40 stencil_width: int | None = None, 41 bint setup: bool = True, 42 ownership_ranges: tuple[Sequence[int], ...] | None = None, 43 comm: Comm | None = None) -> Self: 44 """Create a ``DMDA`` object. 45 46 Collective. 47 48 This routine performs the following steps of the C API: 49 - ``petsc.DMDACreate`` 50 - ``petsc.DMSetDimension`` 51 - ``petsc.DMDASetDof`` 52 - ``petsc.DMDASetSizes`` 53 - ``petsc.DMDASetNumProcs`` 54 - ``petsc.DMDASetOwnershipRanges`` 55 - ``petsc.DMDASetBoundaryType`` 56 - ``petsc.DMDASetStencilType`` 57 - ``petsc.DMDASetStencilWidth`` 58 - ``petsc.DMSetUp`` (optionally) 59 60 Parameters 61 ---------- 62 dim 63 The number of dimensions. 64 dof 65 The number of degrees of freedom. 66 sizes 67 The number of elements in each dimension. 68 proc_sizes 69 The number of processes in x, y, z dimensions. 70 boundary_type 71 The boundary types. 72 stencil_type 73 The ghost/halo stencil type. 74 stencil_width 75 The width of the ghost/halo region. 76 setup 77 Whether to call the setup routine after creating the object. 78 ownership_ranges 79 Local x, y, z element counts, of length equal to ``proc_sizes``, 80 summing to ``sizes``. 81 comm 82 MPI communicator, defaults to `Sys.getDefaultComm`. 83 84 See Also 85 -------- 86 petsc.DMDACreate, petsc.DMSetDimension, petsc.DMDASetDof 87 petsc.DMDASetSizes, petsc.DMDASetNumProcs 88 petsc.DMDASetOwnershipRanges, petsc.DMDASetBoundaryType 89 petsc.DMDASetStencilType, petsc.DMDASetStencilWidth, petsc.DMSetUp 90 91 """ 92 # 93 cdef object arg = None 94 try: arg = tuple(dim) 95 except TypeError: pass 96 else: dim, sizes = None, arg 97 # 98 cdef PetscInt ndim = PETSC_DECIDE 99 cdef PetscInt ndof = PETSC_DECIDE 100 cdef PetscInt M = 1, m = PETSC_DECIDE, *lx = NULL 101 cdef PetscInt N = 1, n = PETSC_DECIDE, *ly = NULL 102 cdef PetscInt P = 1, p = PETSC_DECIDE, *lz = NULL 103 cdef PetscDMBoundaryType btx = DM_BOUNDARY_NONE 104 cdef PetscDMBoundaryType bty = DM_BOUNDARY_NONE 105 cdef PetscDMBoundaryType btz = DM_BOUNDARY_NONE 106 cdef PetscDMDAStencilType stype = DMDA_STENCIL_BOX 107 cdef PetscInt swidth = PETSC_DECIDE 108 # grid and proc sizes 109 cdef object gsizes = sizes 110 cdef object psizes = proc_sizes 111 cdef PetscInt gdim = PETSC_DECIDE 112 cdef PetscInt pdim = PETSC_DECIDE 113 if sizes is not None: 114 gdim = asDims(gsizes, &M, &N, &P) 115 if psizes is not None: 116 pdim = asDims(psizes, &m, &n, &p) 117 if gdim>=0 and pdim>=0: 118 assert gdim == pdim 119 # dim and dof 120 if dim is not None: ndim = asInt(dim) 121 if dof is not None: ndof = asInt(dof) 122 if ndim==PETSC_DECIDE: ndim = gdim 123 if ndof==PETSC_DECIDE: ndof = 1 124 # vertex distribution 125 if ownership_ranges is not None: 126 ownership_ranges = asOwnershipRanges(ownership_ranges, 127 ndim, &m, &n, &p, 128 &lx, &ly, &lz) 129 # periodicity, stencil type & width 130 if boundary_type is not None: 131 asBoundary(boundary_type, &btx, &bty, &btz) 132 if stencil_type is not None: 133 stype = asStencil(stencil_type) 134 if stencil_width is not None: 135 swidth = asInt(stencil_width) 136 if setup and swidth == PETSC_DECIDE: swidth = 0 137 # create the DMDA object 138 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 139 cdef PetscDM newda = NULL 140 CHKERR(DMDACreateND(ccomm, ndim, ndof, 141 M, N, P, m, n, p, lx, ly, lz, 142 btx, bty, btz, stype, swidth, 143 &newda)) 144 if setup and ndim > 0: CHKERR(DMSetUp(newda)) 145 CHKERR(PetscCLEAR(self.obj)); self.dm = newda 146 return self 147 148 def duplicate( 149 self, 150 dof: int | None = None, 151 boundary_type: tuple[DM.BoundaryType | int | str | bool, ...] | None = None, 152 stencil_type: StencilType | None = None, 153 stencil_width: int | None = None) -> DMDA: 154 """Duplicate a DMDA. 155 156 Collective. 157 158 This routine retrieves the information from the DMDA and recreates it. 159 Parameters ``dof``, ``boundary_type``, ``stencil_type``, 160 ``stencil_width`` will be overwritten, if provided. 161 162 Parameters 163 ---------- 164 dof 165 The number of degrees of freedom. 166 boundary_type 167 Boundary types. 168 stencil_type 169 The ghost/halo stencil type. 170 stencil_width 171 The width of the ghost/halo region. 172 173 See Also 174 -------- 175 create, petsc.DMDAGetInfo, petsc.DMSetUp 176 177 """ 178 cdef PetscInt ndim = 0, ndof = 0 179 cdef PetscInt M = 1, N = 1, P = 1 180 cdef PetscInt m = 1, n = 1, p = 1 181 cdef PetscDMBoundaryType btx = DM_BOUNDARY_NONE 182 cdef PetscDMBoundaryType bty = DM_BOUNDARY_NONE 183 cdef PetscDMBoundaryType btz = DM_BOUNDARY_NONE 184 cdef PetscDMDAStencilType stype = DMDA_STENCIL_BOX 185 cdef PetscInt swidth = PETSC_DECIDE 186 CHKERR(DMDAGetInfo(self.dm, 187 &ndim, 188 &M, &N, &P, 189 &m, &n, &p, 190 &ndof, &swidth, 191 &btx, &bty, &btz, 192 &stype)) 193 cdef const PetscInt *lx = NULL, *ly = NULL, *lz = NULL 194 CHKERR(DMDAGetOwnershipRanges(self.dm, &lx, &ly, &lz)) 195 cdef MPI_Comm comm = MPI_COMM_NULL 196 CHKERR(PetscObjectGetComm(<PetscObject>self.dm, &comm)) 197 # 198 if dof is not None: 199 ndof = asInt(dof) 200 if boundary_type is not None: 201 asBoundary(boundary_type, &btx, &bty, &btz) 202 if stencil_type is not None: 203 stype = asStencil(stencil_type) 204 if stencil_width is not None: 205 swidth = asInt(stencil_width) 206 # 207 cdef DMDA da = DMDA() 208 CHKERR(DMDACreateND(comm, ndim, ndof, 209 M, N, P, m, n, p, lx, ly, lz, 210 btx, bty, btz, stype, swidth, 211 &da.dm)) 212 CHKERR(DMSetUp(da.dm)) 213 return da 214 215 # 216 217 def setDim(self, dim: int) -> None: 218 """Set the topological dimension. 219 220 Collective. 221 222 Parameters 223 ---------- 224 dim 225 Topological dimension. 226 227 See Also 228 -------- 229 getDim, petsc.DMSetDimension 230 231 """ 232 return self.setDimension(dim) 233 234 def getDim(self) -> int: 235 """Return the topological dimension. 236 237 Not collective. 238 239 See Also 240 -------- 241 setDim, petsc.DMGetDimension 242 243 """ 244 return self.getDimension() 245 246 def setDof(self, dof: int) -> None: 247 """Set the number of degrees of freedom per vertex. 248 249 Not collective. 250 251 Parameters 252 ---------- 253 dof 254 The number of degrees of freedom. 255 256 See Also 257 -------- 258 getDof, petsc.DMDASetDof 259 260 """ 261 cdef PetscInt ndof = asInt(dof) 262 CHKERR(DMDASetDof(self.dm, ndof)) 263 264 def getDof(self) -> int: 265 """Return the number of degrees of freedom per node. 266 267 Not collective. 268 269 See Also 270 -------- 271 setDof, petsc.DMDAGetInfo 272 273 """ 274 cdef PetscInt dof = 0 275 CHKERR(DMDAGetInfo(self.dm, 276 NULL, 277 NULL, NULL, NULL, 278 NULL, NULL, NULL, 279 &dof, NULL, 280 NULL, NULL, NULL, 281 NULL)) 282 return toInt(dof) 283 284 def setSizes(self, sizes: DimsSpec) -> None: 285 """Set the number of grid points in each dimension. 286 287 Logically collective. 288 289 Parameters 290 ---------- 291 sizes 292 The global (x,), (x, y), or (x, y, z) size. 293 294 See Also 295 -------- 296 getSizes, petsc.DMDASetSizes 297 298 """ 299 cdef tuple gsizes = tuple(sizes) 300 cdef PetscInt gdim = PETSC_DECIDE 301 cdef PetscInt M = 1 302 cdef PetscInt N = 1 303 cdef PetscInt P = 1 304 gdim = asDims(gsizes, &M, &N, &P) 305 cdef PetscInt dim = PETSC_DECIDE 306 CHKERR(DMDAGetDim(self.dm, &dim)) 307 if dim == PETSC_DECIDE: 308 CHKERR(DMSetDimension(self.dm, gdim)) 309 CHKERR(DMDASetSizes(self.dm, M, N, P)) 310 311 def getSizes(self) -> tuple[int, ...]: 312 """Return the number of grid points in each dimension. 313 314 Not collective. 315 316 See Also 317 -------- 318 setSizes, petsc.DMDAGetInfo 319 320 """ 321 cdef PetscInt dim = 0 322 cdef PetscInt M = PETSC_DECIDE 323 cdef PetscInt N = PETSC_DECIDE 324 cdef PetscInt P = PETSC_DECIDE 325 CHKERR(DMDAGetInfo(self.dm, 326 &dim, 327 &M, &N, &P, 328 NULL, NULL, NULL, 329 NULL, NULL, 330 NULL, NULL, NULL, 331 NULL)) 332 return toDims(dim, M, N, P) 333 334 def setProcSizes(self, proc_sizes: DimsSpec) -> None: 335 """Set the number of processes in each dimension. 336 337 Logically collective. 338 339 Parameters 340 ---------- 341 proc_sizes 342 The number of processes in (x,), (x, y), or (x, y, z) dimensions. 343 344 See Also 345 -------- 346 getProcSizes, petsc.DMDASetNumProcs 347 348 """ 349 cdef tuple psizes = tuple(proc_sizes) 350 cdef PetscInt pdim = PETSC_DECIDE 351 cdef PetscInt m = PETSC_DECIDE 352 cdef PetscInt n = PETSC_DECIDE 353 cdef PetscInt p = PETSC_DECIDE 354 pdim = asDims(psizes, &m, &n, &p) 355 cdef PetscInt dim = PETSC_DECIDE 356 CHKERR(DMDAGetDim(self.dm, &dim)) 357 if dim == PETSC_DECIDE: 358 CHKERR(DMSetDimension(self.dm, pdim)) 359 CHKERR(DMDASetNumProcs(self.dm, m, n, p)) 360 361 def getProcSizes(self) -> tuple[int, ...]: 362 """Return the number of processes in each dimension. 363 364 Not collective. 365 366 See Also 367 -------- 368 setProcSizes, petsc.DMDAGetInfo 369 370 """ 371 cdef PetscInt dim = 0 372 cdef PetscInt m = PETSC_DECIDE 373 cdef PetscInt n = PETSC_DECIDE 374 cdef PetscInt p = PETSC_DECIDE 375 CHKERR(DMDAGetInfo(self.dm, 376 &dim, 377 NULL, NULL, NULL, 378 &m, &n, &p, 379 NULL, NULL, 380 NULL, NULL, NULL, 381 NULL)) 382 return toDims(dim, m, n, p) 383 384 def setBoundaryType( 385 self, 386 boundary_type: tuple[DM.BoundaryType | int | str | bool, ...]) -> None: 387 """Set the type of ghost nodes on domain boundaries. 388 389 Not collective. 390 391 Parameters 392 ---------- 393 boundary_type 394 The boundary type in (x), (x, y), or (x, y, z) dimensions. 395 396 See Also 397 -------- 398 getBoundaryType, petsc.DMDASetBoundaryType 399 400 """ 401 cdef PetscDMBoundaryType btx = DM_BOUNDARY_NONE 402 cdef PetscDMBoundaryType bty = DM_BOUNDARY_NONE 403 cdef PetscDMBoundaryType btz = DM_BOUNDARY_NONE 404 asBoundary(boundary_type, &btx, &bty, &btz) 405 CHKERR(DMDASetBoundaryType(self.dm, btx, bty, btz)) 406 407 def getBoundaryType(self) -> tuple[DM.BoundaryType, ...]: 408 """Return the type of ghost nodes at boundary in each dimension. 409 410 Not collective. 411 412 See Also 413 -------- 414 setBoundaryType 415 416 """ 417 cdef PetscInt dim = 0 418 cdef PetscDMBoundaryType btx = DM_BOUNDARY_NONE 419 cdef PetscDMBoundaryType bty = DM_BOUNDARY_NONE 420 cdef PetscDMBoundaryType btz = DM_BOUNDARY_NONE 421 CHKERR(DMDAGetInfo(self.dm, 422 &dim, 423 NULL, NULL, NULL, 424 NULL, NULL, NULL, 425 NULL, NULL, 426 &btx, &bty, &btz, 427 NULL)) 428 return toDims(dim, btx, bty, btz) 429 430 def setStencilType(self, stencil_type: StencilType) -> None: 431 """Set the stencil type. 432 433 Logically collective. 434 435 Parameters 436 ---------- 437 stype 438 The stencil type. 439 440 See Also 441 -------- 442 getStencilType, setStencil, petsc.DMDASetStencilType 443 444 """ 445 cdef PetscDMDAStencilType stype = asStencil(stencil_type) 446 CHKERR(DMDASetStencilType(self.dm, stype)) 447 448 def getStencilType(self) -> StencilType: 449 """Return the stencil type. 450 451 Not collective. 452 453 See Also 454 -------- 455 setStencilType, petsc.DMDAGetInfo 456 457 """ 458 cdef PetscDMDAStencilType stype = DMDA_STENCIL_BOX 459 CHKERR(DMDAGetInfo(self.dm, 460 NULL, 461 NULL, NULL, NULL, 462 NULL, NULL, NULL, 463 NULL, NULL, 464 NULL, NULL, NULL, 465 &stype)) 466 return stype 467 468 def setStencilWidth(self, stencil_width: int) -> None: 469 """Set the stencil width. 470 471 Logically collective. 472 473 Parameters 474 ---------- 475 stencil_width 476 The stencil width. 477 478 See Also 479 -------- 480 getStencilWidth, setStencil, petsc.DMDASetStencilWidth 481 482 """ 483 cdef PetscInt swidth = asInt(stencil_width) 484 CHKERR(DMDASetStencilWidth(self.dm, swidth)) 485 486 def getStencilWidth(self) -> int: 487 """Return the stencil width. 488 489 Not collective. 490 491 See Also 492 -------- 493 setStencilWidth 494 495 """ 496 cdef PetscInt swidth = 0 497 CHKERR(DMDAGetInfo(self.dm, 498 NULL, 499 NULL, NULL, NULL, 500 NULL, NULL, NULL, 501 NULL, &swidth, 502 NULL, NULL, NULL, 503 NULL)) 504 return toInt(swidth) 505 506 def setStencil( 507 self, 508 stencil_type: StencilType, 509 stencil_width: int) -> None: 510 """Set the stencil type and width. 511 512 Not collective. 513 514 Parameters 515 ---------- 516 stencil_type 517 The stencil type. 518 stencil_width 519 The stencil width. 520 521 See Also 522 -------- 523 setStencilWidth, setStencilType, petsc.DMDASetStencilType 524 petsc.DMDASetStencilWidth 525 526 """ 527 cdef PetscDMDAStencilType stype = asStencil(stencil_type) 528 cdef PetscInt swidth = asInt(stencil_width) 529 CHKERR(DMDASetStencilType(self.dm, stype)) 530 CHKERR(DMDASetStencilWidth(self.dm, swidth)) 531 532 def getStencil(self) -> tuple[StencilType, int]: 533 """Return the stencil type and width. 534 535 Not collective. 536 537 See Also 538 -------- 539 getStencilType, getStencilWidth 540 541 """ 542 cdef PetscDMDAStencilType stype = DMDA_STENCIL_BOX 543 cdef PetscInt swidth = 0 544 CHKERR(DMDAGetInfo(self.dm, 545 NULL, 546 NULL, NULL, NULL, 547 NULL, NULL, NULL, 548 NULL, &swidth, 549 NULL, NULL, NULL, 550 &stype)) 551 return (toStencil(stype), toInt(swidth)) 552 553 # 554 555 def getRanges(self) -> tuple[tuple[int, int], ...]: 556 """Return the ranges of the owned local region in each dimension. 557 558 Not collective. 559 560 Excluding ghost nodes. 561 562 See Also 563 -------- 564 getGhostRanges, getOwnershipRanges, getCorners, getGhostCorners 565 petsc.DMDAGetCorners 566 567 """ 568 cdef PetscInt dim=0, x=0, y=0, z=0, m=0, n=0, p=0 569 CHKERR(DMDAGetDim(self.dm, &dim)) 570 CHKERR(DMDAGetCorners(self.dm, 571 &x, &y, &z, 572 &m, &n, &p)) 573 return ((toInt(x), toInt(x+m)), 574 (toInt(y), toInt(y+n)), 575 (toInt(z), toInt(z+p)))[:<Py_ssize_t>dim] 576 577 def getGhostRanges(self) -> tuple[tuple[int, int], ...]: 578 """Return the ranges of the local region in each dimension, including ghost nodes. 579 580 Not collective. 581 582 See Also 583 -------- 584 getRanges, getOwnershipRanges, getCorners, getGhostCorners 585 petsc.DMDAGetGhostCorners 586 587 """ 588 cdef PetscInt dim=0, x=0, y=0, z=0, m=0, n=0, p=0 589 CHKERR(DMDAGetDim(self.dm, &dim)) 590 CHKERR(DMDAGetGhostCorners(self.dm, 591 &x, &y, &z, 592 &m, &n, &p)) 593 return ((toInt(x), toInt(x+m)), 594 (toInt(y), toInt(y+n)), 595 (toInt(z), toInt(z+p)))[:<Py_ssize_t>dim] 596 597 def getOwnershipRanges(self) -> tuple[ArrayInt, ...]: 598 """Return the ranges of indices in each dimension owned by each process. 599 600 Not collective. 601 602 These numbers are not multiplied by the number of DOFs per node. 603 604 See Also 605 -------- 606 getRanges, getGhostRanges, getCorners, getGhostCorners 607 petsc.DMDAGetOwnershipRanges 608 609 """ 610 cdef PetscInt dim=0, m=0, n=0, p=0 611 cdef const PetscInt *lx = NULL, *ly = NULL, *lz = NULL 612 CHKERR(DMDAGetInfo(self.dm, 613 &dim, 614 NULL, NULL, NULL, 615 &m, &n, &p, 616 NULL, NULL, 617 NULL, NULL, NULL, 618 NULL)) 619 CHKERR(DMDAGetOwnershipRanges(self.dm, &lx, &ly, &lz)) 620 return toOwnershipRanges(dim, m, n, p, lx, ly, lz) 621 622 def getCorners(self) -> tuple[tuple[int, ...], tuple[int, ...]]: 623 """Return the lower left corner and the sizes of the owned local region. 624 625 Not collective. 626 627 Returns the global (x,y,z) indices of the lower left corner (first 628 tuple) and size of the local region (second tuple). 629 630 Excluding ghost points. 631 632 The corner information is independent of the number of degrees of 633 freedom per node. Thus the returned values can be thought of as 634 coordinates on a logical grid, where each grid point has (potentially) 635 several degrees of freedom. 636 637 See Also 638 -------- 639 getRanges, getGhostRanges, getOwnershipRanges, getGhostCorners 640 petsc.DMDAGetCorners 641 642 """ 643 cdef PetscInt dim=0, x=0, y=0, z=0, m=0, n=0, p=0 644 CHKERR(DMDAGetDim(self.dm, &dim)) 645 CHKERR(DMDAGetCorners(self.dm, 646 &x, &y, &z, 647 &m, &n, &p)) 648 return ((toInt(x), toInt(y), toInt(z))[:<Py_ssize_t>dim], 649 (toInt(m), toInt(n), toInt(p))[:<Py_ssize_t>dim]) 650 651 def getGhostCorners(self) -> tuple[tuple[int, ...], tuple[int, ...]]: 652 """Return the lower left corner and the size of the ghosted local region. 653 654 Not collective. 655 656 Returns the global (x,y,z) indices of the lower left corner (first 657 tuple) and size of the local region (second tuple). 658 659 See Also 660 -------- 661 getRanges, getGhostRanges, getOwnershipRanges, getCorners 662 petsc.DMDAGetGhostCorners 663 664 """ 665 cdef PetscInt dim=0, x=0, y=0, z=0, m=0, n=0, p=0 666 CHKERR(DMDAGetDim(self.dm, &dim)) 667 CHKERR(DMDAGetGhostCorners(self.dm, 668 &x, &y, &z, 669 &m, &n, &p)) 670 return ((toInt(x), toInt(y), toInt(z))[:<Py_ssize_t>dim], 671 (toInt(m), toInt(n), toInt(p))[:<Py_ssize_t>dim]) 672 673 # 674 675 def setFieldName(self, field: int, name: str) -> None: 676 """Set the name of individual field components. 677 678 Logically collective. 679 680 Parameters 681 ---------- 682 field 683 The field number for the DMDA (``0``, ``1``, ..., ``dof-1``), 684 where ``dof`` indicates the number of degrees of freedom per node 685 within the `DMDA`. 686 name 687 The name of the field (component). 688 689 See Also 690 -------- 691 getFieldName, petsc.DMDASetFieldName 692 693 """ 694 cdef PetscInt ival = asInt(field) 695 cdef const char *cval = NULL 696 name = str2bytes(name, &cval) 697 CHKERR(DMDASetFieldName(self.dm, ival, cval)) 698 699 def getFieldName(self, field: int) -> str: 700 """Return the name of an individual field component. 701 702 Not collective. 703 704 Parameters 705 ---------- 706 field 707 The field number for the DMDA (``0``, ``1``, ..., ``dof-1``), 708 where ``dof`` indicates the number of degrees of freedom per node 709 within the `DMDA`. 710 711 See Also 712 -------- 713 setFieldName, petsc.DMDAGetFieldName 714 715 """ 716 cdef PetscInt ival = asInt(field) 717 cdef const char *cval = NULL 718 CHKERR(DMDAGetFieldName(self.dm, ival, &cval)) 719 return bytes2str(cval) 720 721 # 722 723 def getVecArray(self, Vec vec, readonly: bool = False) -> Any: 724 """Get access to the vector as laid out on a N-d grid. 725 726 Logically collective. 727 728 Parameters 729 ---------- 730 vec 731 The vector to which access is being requested. 732 readonly 733 Request read-only access. 734 735 See Also 736 -------- 737 Vec.getArray, petsc.DMDAVecGetArray, petsc.DMDAVecGetArrayDOF 738 739 """ 740 return _DMDA_Vec_array(self, vec, readonly) 741 742 # 743 744 def setUniformCoordinates( 745 self, 746 xmin: float = 0, 747 xmax: float = 1, 748 ymin: float = 0, 749 ymax: float = 1, 750 zmin: float = 0, 751 zmax: float = 1) -> None: 752 """Set the DMDA coordinates to be a uniform grid. 753 754 Collective. 755 756 Parameters 757 ---------- 758 xmin 759 The minimum in the ``x`` dimension. 760 xmax 761 The maximum in the ``x`` dimension. 762 ymin 763 The minimum in the ``y`` dimension (value ignored for 1 dimensional 764 problems). 765 ymax 766 The maximum in the ``y`` dimension (value ignored for 1 dimensional 767 problems). 768 zmin 769 The minimum in the ``z`` dimension (value ignored for 1 or 2 770 dimensional problems). 771 zmax 772 The maximum in the ``z`` dimension (value ignored for 1 or 2 773 dimensional problems). 774 775 See Also 776 -------- 777 petsc.DMDASetUniformCoordinates 778 779 """ 780 cdef PetscReal _xmin = asReal(xmin), _xmax = asReal(xmax) 781 cdef PetscReal _ymin = asReal(ymin), _ymax = asReal(ymax) 782 cdef PetscReal _zmin = asReal(zmin), _zmax = asReal(zmax) 783 CHKERR(DMDASetUniformCoordinates(self.dm, 784 _xmin, _xmax, 785 _ymin, _ymax, 786 _zmin, _zmax)) 787 788 def setCoordinateName(self, index: int, name: str) -> None: 789 """Set the name of the coordinate dimension. 790 791 Logically collective. 792 793 Parameters 794 ---------- 795 index 796 The coordinate number for the DMDA (``0``, ``1``, ..., ``dim-1``). 797 name 798 The name of the coordinate. 799 800 See Also 801 -------- 802 getCoordinateName, petsc.DMDASetCoordinateName 803 804 """ 805 cdef PetscInt ival = asInt(index) 806 cdef const char *cval = NULL 807 name = str2bytes(name, &cval) 808 CHKERR(DMDASetCoordinateName(self.dm, ival, cval)) 809 810 def getCoordinateName(self, index: int) -> str: 811 """Return the name of a coordinate dimension. 812 813 Not collective. 814 815 Parameters 816 ---------- 817 index 818 The coordinate number for the DMDA (``0``, ``1``, ..., ``dim-1``). 819 820 See Also 821 -------- 822 setCoordinateName, petsc.DMDAGetCoordinateName 823 824 """ 825 cdef PetscInt ival = asInt(index) 826 cdef const char *cval = NULL 827 CHKERR(DMDAGetCoordinateName(self.dm, ival, &cval)) 828 return bytes2str(cval) 829 830 # 831 832 def createNaturalVec(self) -> Vec: 833 """Create a vector that will hold values in the natural numbering. 834 835 Collective. 836 837 The number of local entries in the vector on each process is the same 838 as in a vector created with `DM.createGlobalVec`. 839 840 See Also 841 -------- 842 petsc.DMDACreateNaturalVector 843 844 """ 845 cdef Vec vn = Vec() 846 CHKERR(DMDACreateNaturalVector(self.dm, &vn.vec)) 847 return vn 848 849 def globalToNatural( 850 self, 851 Vec vg, 852 Vec vn, 853 addv: InsertMode | None = None) -> None: 854 """Map values to the "natural" grid ordering. 855 856 Neighborwise collective. 857 858 You must call `createNaturalVec` before using this routine. 859 860 Parameters 861 ---------- 862 vg 863 The global vector in a grid ordering. 864 vn 865 The global vector in a "natural" ordering. 866 addv 867 The insertion mode. 868 869 See Also 870 -------- 871 naturalToGlobal, petsc.DMDAGlobalToNaturalBegin 872 petsc.DMDAGlobalToNaturalEnd 873 874 """ 875 cdef PetscInsertMode im = insertmode(addv) 876 CHKERR(DMDAGlobalToNaturalBegin(self.dm, vg.vec, im, vn.vec)) 877 CHKERR(DMDAGlobalToNaturalEnd(self.dm, vg.vec, im, vn.vec)) 878 879 def naturalToGlobal( 880 self, 881 Vec vn, 882 Vec vg, 883 addv: InsertMode | None = None) -> None: 884 """Map values the to grid ordering. 885 886 Neighborwise collective. 887 888 Parameters 889 ---------- 890 vn 891 The global vector in a natural ordering. 892 vg 893 the global vector in a grid ordering. 894 addv 895 The insertion mode. 896 897 See Also 898 -------- 899 globalToNatural, petsc.DMDANaturalToGlobalBegin 900 petsc.DMDANaturalToGlobalEnd 901 902 """ 903 cdef PetscInsertMode im = insertmode(addv) 904 CHKERR(DMDANaturalToGlobalBegin(self.dm, vn.vec, im, vg.vec)) 905 CHKERR(DMDANaturalToGlobalEnd(self.dm, vn.vec, im, vg.vec)) 906 907 # 908 909 def getAO(self) -> AO: 910 """Return the application ordering context for a distributed array. 911 912 Collective. 913 914 The returned `AO` maps to the natural grid ordering that would be 915 used for the `DMDA` if only 1 processor were employed (ordering most 916 rapidly in the x-dimension, then y, then z). Multiple degrees of 917 freedom are numbered for each node (rather than 1 component for the 918 whole grid, then the next component, etc.). 919 920 See Also 921 -------- 922 petsc.DMDAGetAO 923 924 """ 925 cdef AO ao = AO() 926 CHKERR(DMDAGetAO(self.dm, &ao.ao)) 927 CHKERR(PetscINCREF(ao.obj)) 928 return ao 929 930 def getScatter(self) -> tuple[Scatter, Scatter]: 931 """Return the global-to-local, and local-to-local scatter contexts. 932 933 Collective. 934 935 See Also 936 -------- 937 petsc.DMDAGetScatter 938 939 """ 940 cdef Scatter l2g = Scatter() 941 cdef Scatter g2l = Scatter() 942 CHKERR(DMDAGetScatter(self.dm, &l2g.sct, &g2l.sct)) 943 CHKERR(PetscINCREF(l2g.obj)) 944 CHKERR(PetscINCREF(g2l.obj)) 945 return (l2g, g2l) 946 947 # 948 949 def setRefinementFactor( 950 self, 951 refine_x: int = 2, 952 refine_y: int = 2, 953 refine_z: int = 2) -> None: 954 """Set the ratios for the DMDA grid refinement. 955 956 Logically collective. 957 958 Parameters 959 ---------- 960 refine_x 961 Ratio of fine grid to coarse in x dimension. 962 refine_y 963 Ratio of fine grid to coarse in y dimension. 964 refine_z 965 Ratio of fine grid to coarse in z dimension. 966 967 See Also 968 -------- 969 getRefinementFactor, petsc.DMDASetRefinementFactor 970 971 """ 972 cdef PetscInt refine[3] 973 refine[0] = asInt(refine_x) 974 refine[1] = asInt(refine_y) 975 refine[2] = asInt(refine_z) 976 CHKERR(DMDASetRefinementFactor(self.dm, 977 refine[0], 978 refine[1], 979 refine[2])) 980 981 def getRefinementFactor(self) -> tuple[int, ...]: 982 """Return the ratios that the DMDA grid is refined in each dimension. 983 984 Not collective. 985 986 See Also 987 -------- 988 setRefinementFactor, petsc.DMDAGetRefinementFactor 989 990 """ 991 cdef PetscInt dim = 0, refine[3] 992 CHKERR(DMDAGetDim(self.dm, &dim)) 993 CHKERR(DMDAGetRefinementFactor(self.dm, 994 &refine[0], 995 &refine[1], 996 &refine[2])) 997 return tuple([toInt(refine[i]) for 0 <= i < dim]) 998 999 def setInterpolationType(self, interp_type: InterpolationType) -> None: 1000 """Set the type of interpolation. 1001 1002 Logically collective. 1003 1004 You should call this on the coarser of the two DMDAs you pass to 1005 `DM.createInterpolation`. 1006 1007 Parameters 1008 ---------- 1009 interp_type 1010 The interpolation type. 1011 1012 See Also 1013 -------- 1014 getInterpolationType, petsc.DMDASetInterpolationType 1015 1016 """ 1017 cdef PetscDMDAInterpolationType ival = dainterpolationtype(interp_type) 1018 CHKERR(DMDASetInterpolationType(self.dm, ival)) 1019 1020 def getInterpolationType(self) -> InterpolationType: 1021 """Return the type of interpolation. 1022 1023 Not collective. 1024 1025 See Also 1026 -------- 1027 setInterpolationType, petsc.DMDAGetInterpolationType 1028 1029 """ 1030 cdef PetscDMDAInterpolationType ival = DMDA_INTERPOLATION_Q0 1031 CHKERR(DMDAGetInterpolationType(self.dm, &ival)) 1032 return <long>ival 1033 1034 # 1035 1036 def setElementType(self, elem_type: ElementType | str) -> None: 1037 """Set the element type to be returned by `getElements`. 1038 1039 Not collective. 1040 1041 See Also 1042 -------- 1043 getElementType, petsc.DMDASetElementType 1044 1045 """ 1046 cdef PetscDMDAElementType ival = daelementtype(elem_type) 1047 CHKERR(DMDASetElementType(self.dm, ival)) 1048 1049 # FIXME: Return type 1050 def getElementType(self) -> ElementType: 1051 """Return the element type to be returned by `getElements`. 1052 1053 Not collective. 1054 1055 See Also 1056 -------- 1057 setElementType, petsc.DMDAGetElementType 1058 1059 """ 1060 cdef PetscDMDAElementType ival = DMDA_ELEMENT_Q1 1061 CHKERR(DMDAGetElementType(self.dm, &ival)) 1062 return <long>ival 1063 1064 def getElements(self, elem_type: ElementType | None = None) -> ArrayInt: 1065 """Return an array containing the indices of all the local elements. 1066 1067 Not collective. 1068 1069 The elements are in local coordinates. 1070 1071 Each process uniquely owns a subset of the elements. That is, no 1072 element is owned by two or more processes. 1073 1074 Parameters 1075 ---------- 1076 elem_type 1077 The element type. 1078 1079 See Also 1080 -------- 1081 petsc.DMDAGetElements 1082 1083 """ 1084 cdef PetscInt dim=0 1085 cdef PetscDMDAElementType etype 1086 cdef PetscInt nel=0, nen=0 1087 cdef const PetscInt *elems=NULL 1088 cdef object elements 1089 CHKERR(DMDAGetDim(self.dm, &dim)) 1090 if elem_type is not None: 1091 etype = daelementtype(elem_type) 1092 CHKERR(DMDASetElementType(self.dm, etype)) 1093 try: 1094 CHKERR(DMDAGetElements(self.dm, &nel, &nen, &elems)) 1095 elements = array_i(nel*nen, elems) 1096 elements.shape = (toInt(nel), toInt(nen)) 1097 finally: 1098 CHKERR(DMDARestoreElements(self.dm, &nel, &nen, &elems)) 1099 return elements 1100 1101 # 1102 1103 property dim: 1104 """The grid dimension.""" 1105 def __get__(self) -> int: 1106 return self.getDim() 1107 1108 property dof: 1109 """The number of DOFs associated with each stratum of the grid.""" 1110 def __get__(self) -> int: 1111 return self.getDof() 1112 1113 property sizes: 1114 """The global dimension.""" 1115 def __get__(self) -> tuple[int, ...]: 1116 return self.getSizes() 1117 1118 property proc_sizes: 1119 """The number of processes in each dimension in the global decomposition.""" 1120 def __get__(self) -> tuple[int, ...]: 1121 return self.getProcSizes() 1122 1123 property boundary_type: 1124 """Boundary types in each dimension.""" 1125 def __get__(self) -> tuple[DM.BoundaryType, ...]: 1126 return self.getBoundaryType() 1127 1128 property stencil: 1129 """Stencil type and width.""" 1130 def __get__(self) -> tuple[StencilType, int]: 1131 return self.getStencil() 1132 1133 property stencil_type: 1134 """Stencil type.""" 1135 def __get__(self) -> str: 1136 return self.getStencilType() 1137 1138 property stencil_width: 1139 """Elementwise stencil width.""" 1140 def __get__(self) -> int: 1141 return self.getStencilWidth() 1142 1143 property ranges: 1144 """Ranges of the local region in each dimension.""" 1145 def __get__(self) -> tuple[tuple[int, int], ...]: 1146 return self.getRanges() 1147 1148 property ghost_ranges: 1149 """Ranges of local region, including ghost nodes.""" 1150 def __get__(self) -> tuple[tuple[int, int], ...]: 1151 return self.getGhostRanges() 1152 1153 property corners: 1154 """The lower left corner and size of local region in each dimension.""" 1155 def __get__(self) -> tuple[tuple[int, ...], tuple[int, ...]]: 1156 return self.getCorners() 1157 1158 property ghost_corners: 1159 """The lower left corner and size of local region in each dimension.""" 1160 def __get__(self) -> tuple[tuple[int, ...], tuple[int, ...]]: 1161 return self.getGhostCorners() 1162 1163 # backward compatibility 1164 createNaturalVector = createNaturalVec 1165 1166 1167# backward compatibility alias 1168DA = DMDA 1169 1170# -------------------------------------------------------------------- 1171 1172del DMDAStencilType 1173del DMDAInterpolationType 1174del DMDAElementType 1175 1176# -------------------------------------------------------------------- 1177