1# -------------------------------------------------------------------- 2 3class DMSwarmType(object): 4 """Swarm types.""" 5 BASIC = DMSWARM_BASIC 6 PIC = DMSWARM_PIC 7 8 9class DMSwarmMigrateType(object): 10 """Swarm migration types.""" 11 MIGRATE_BASIC = DMSWARM_MIGRATE_BASIC 12 MIGRATE_DMCELLNSCATTER = DMSWARM_MIGRATE_DMCELLNSCATTER 13 MIGRATE_DMCELLEXACT = DMSWARM_MIGRATE_DMCELLEXACT 14 MIGRATE_USER = DMSWARM_MIGRATE_USER 15 16 17class DMSwarmCollectType(object): 18 """Swarm collection types.""" 19 COLLECT_BASIC = DMSWARM_COLLECT_BASIC 20 COLLECT_DMDABOUNDINGBOX = DMSWARM_COLLECT_DMDABOUNDINGBOX 21 COLLECT_GENERAL = DMSWARM_COLLECT_GENERAL 22 COLLECT_USER = DMSWARM_COLLECT_USER 23 24 25class DMSwarmPICLayoutType(object): 26 """Swarm PIC layout types.""" 27 LAYOUT_REGULAR = DMSWARMPIC_LAYOUT_REGULAR 28 LAYOUT_GAUSS = DMSWARMPIC_LAYOUT_GAUSS 29 LAYOUT_SUBDIVISION = DMSWARMPIC_LAYOUT_SUBDIVISION 30 31 32cdef class DMSwarm(DM): 33 """A `DM` object used to represent arrays of data (fields) of arbitrary type.""" 34 Type = DMSwarmType 35 MigrateType = DMSwarmMigrateType 36 CollectType = DMSwarmCollectType 37 PICLayoutType = DMSwarmPICLayoutType 38 39 def create(self, comm: Comm | None = None) -> Self: 40 """Create an empty DM object and set its type to `DM.Type.SWARM`. 41 42 Collective. 43 44 DMs are the abstract objects in PETSc that mediate between meshes and 45 discretizations and the algebraic solvers, time integrators, and 46 optimization algorithms. 47 48 Parameters 49 ---------- 50 comm 51 MPI communicator, defaults to `Sys.getDefaultComm`. 52 53 See Also 54 -------- 55 petsc.DMCreate, petsc.DMSetType 56 57 """ 58 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 59 cdef PetscDM newdm = NULL 60 CHKERR(DMCreate(ccomm, &newdm)) 61 CHKERR(PetscCLEAR(self.obj)); self.dm = newdm 62 CHKERR(DMSetType(self.dm, DMSWARM)) 63 return self 64 65 def createGlobalVectorFromField(self, fieldname: str) -> Vec: 66 """Create a global `Vec` object associated with a given field. 67 68 Collective. 69 70 The vector must be returned to the `DMSwarm` using a matching call to 71 `destroyGlobalVectorFromField`. 72 73 Parameters 74 ---------- 75 fieldname 76 The textual name given to a registered field. 77 78 See Also 79 -------- 80 destroyGlobalVectorFromField, petsc.DMSwarmCreateGlobalVectorFromField 81 82 """ 83 cdef const char *cfieldname = NULL 84 cdef Vec vg = Vec() 85 fieldname = str2bytes(fieldname, &cfieldname) 86 CHKERR(DMSwarmCreateGlobalVectorFromField(self.dm, cfieldname, &vg.vec)) 87 return vg 88 89 def destroyGlobalVectorFromField(self, fieldname: str) -> None: 90 """Destroy the global `Vec` object associated with a given field. 91 92 Collective. 93 94 Parameters 95 ---------- 96 fieldname 97 The textual name given to a registered field. 98 99 See Also 100 -------- 101 createGlobalVectorFromField, petsc.DMSwarmDestroyGlobalVectorFromField 102 103 """ 104 cdef const char *cfieldname = NULL 105 cdef PetscVec vec = NULL 106 fieldname = str2bytes(fieldname, &cfieldname) 107 CHKERR(DMSwarmDestroyGlobalVectorFromField(self.dm, cfieldname, &vec)) 108 109 def createGlobalVectorFromFields(self, fieldnames: Sequence[str]) -> Vec: 110 """Create a global `Vec` object associated with a given set of fields. 111 112 Collective. 113 114 The vector must be returned to the `DMSwarm` using a matching call to 115 `destroyGlobalVectorFromFields`. 116 117 Parameters 118 ---------- 119 fieldnames 120 The textual name given to each registered field. 121 122 See Also 123 -------- 124 destroyGlobalVectorFromFields, petsc.DMSwarmCreateGlobalVectorFromFields 125 126 """ 127 cdef const char *cval = NULL 128 cdef PetscInt i = 0, nf = <PetscInt> len(fieldnames) 129 cdef const char **cfieldnames = NULL 130 cdef object unused = oarray_p(empty_p(nf), NULL, <void**>&cfieldnames) 131 fieldnames = list(fieldnames) 132 for i from 0 <= i < nf: 133 str2bytes(fieldnames[i], &cval) 134 cfieldnames[i] = cval 135 cdef Vec vg = Vec() 136 CHKERR(DMSwarmCreateGlobalVectorFromFields(self.dm, nf, cfieldnames, &vg.vec)) 137 return vg 138 139 def destroyGlobalVectorFromFields(self, fieldnames: Sequence[str]) -> None: 140 """Destroy the global `Vec` object associated with a given set of fields. 141 142 Collective. 143 144 Parameters 145 ---------- 146 fieldnames 147 The textual name given to each registered field. 148 149 See Also 150 -------- 151 createGlobalVectorFromFields, petsc.DMSwarmDestroyGlobalVectorFromFields 152 153 """ 154 cdef const char *cval = NULL 155 cdef PetscInt i = 0, nf = <PetscInt> len(fieldnames) 156 cdef const char **cfieldnames = NULL 157 cdef object unused = oarray_p(empty_p(nf), NULL, <void**>&cfieldnames) 158 fieldnames = list(fieldnames) 159 for i from 0 <= i < nf: 160 str2bytes(fieldnames[i], &cval) 161 cfieldnames[i] = cval 162 cdef PetscVec vec = NULL 163 CHKERR(DMSwarmDestroyGlobalVectorFromFields(self.dm, nf, cfieldnames, &vec)) 164 165 def createLocalVectorFromField(self, fieldname: str) -> Vec: 166 """Create a local `Vec` object associated with a given field. 167 168 Collective. 169 170 The vector must be returned to the `DMSwarm` using a matching call 171 to `destroyLocalVectorFromField`. 172 173 Parameters 174 ---------- 175 fieldname 176 The textual name given to a registered field. 177 178 See Also 179 -------- 180 destroyLocalVectorFromField, petsc.DMSwarmCreateLocalVectorFromField 181 182 """ 183 cdef const char *cfieldname = NULL 184 cdef Vec vl = Vec() 185 fieldname = str2bytes(fieldname, &cfieldname) 186 CHKERR(DMSwarmCreateLocalVectorFromField(self.dm, cfieldname, &vl.vec)) 187 return vl 188 189 def destroyLocalVectorFromField(self, fieldname: str) -> None: 190 """Destroy the local `Vec` object associated with a given field. 191 192 Collective. 193 194 Parameters 195 ---------- 196 fieldname 197 The textual name given to a registered field. 198 199 See Also 200 -------- 201 createLocalVectorFromField, petsc.DMSwarmDestroyLocalVectorFromField 202 203 """ 204 cdef const char *cfieldname = NULL 205 cdef PetscVec vec = NULL 206 fieldname = str2bytes(fieldname, &cfieldname) 207 CHKERR(DMSwarmDestroyLocalVectorFromField(self.dm, cfieldname, &vec)) 208 209 def createLocalVectorFromFields(self, fieldnames: Sequence[str]) -> Vec: 210 """Create a local `Vec` object associated with a given set of fields. 211 212 Collective. 213 214 The vector must be returned to the `DMSwarm` using a matching call to 215 `destroyLocalVectorFromFields`. 216 217 Parameters 218 ---------- 219 fieldnames 220 The textual name given to each registered field. 221 222 See Also 223 -------- 224 destroyLocalVectorFromFields, petsc.DMSwarmCreateLocalVectorFromFields 225 226 """ 227 cdef const char *cval = NULL 228 cdef PetscInt i = 0, nf = <PetscInt> len(fieldnames) 229 cdef const char **cfieldnames = NULL 230 cdef object unused = oarray_p(empty_p(nf), NULL, <void**>&cfieldnames) 231 fieldnames = list(fieldnames) 232 for i from 0 <= i < nf: 233 str2bytes(fieldnames[i], &cval) 234 cfieldnames[i] = cval 235 cdef Vec vl = Vec() 236 CHKERR(DMSwarmCreateLocalVectorFromFields(self.dm, nf, cfieldnames, &vl.vec)) 237 return vl 238 239 def destroyLocalVectorFromFields(self, fieldnames: Sequence[str]) -> None: 240 """Destroy the local `Vec` object associated with a given set of fields. 241 242 Collective. 243 244 Parameters 245 ---------- 246 fieldnames 247 The textual name given to each registered field. 248 249 See Also 250 -------- 251 createLocalVectorFromFields, petsc.DMSwarmDestroyLocalVectorFromFields 252 253 """ 254 cdef const char *cval = NULL 255 cdef PetscInt i = 0, nf = <PetscInt> len(fieldnames) 256 cdef const char **cfieldnames = NULL 257 cdef object unused = oarray_p(empty_p(nf), NULL, <void**>&cfieldnames) 258 fieldnames = list(fieldnames) 259 for i from 0 <= i < nf: 260 str2bytes(fieldnames[i], &cval) 261 cfieldnames[i] = cval 262 cdef PetscVec vec = NULL 263 CHKERR(DMSwarmDestroyLocalVectorFromFields(self.dm, nf, cfieldnames, &vec)) 264 265 def initializeFieldRegister(self) -> None: 266 """Initiate the registration of fields to a `DMSwarm`. 267 268 Collective. 269 270 After all fields have been registered, you must call `finalizeFieldRegister`. 271 272 See Also 273 -------- 274 finalizeFieldRegister, petsc.DMSwarmInitializeFieldRegister 275 276 """ 277 CHKERR(DMSwarmInitializeFieldRegister(self.dm)) 278 279 def finalizeFieldRegister(self) -> None: 280 """Finalize the registration of fields to a `DMSwarm`. 281 282 Collective. 283 284 See Also 285 -------- 286 initializeFieldRegister, petsc.DMSwarmFinalizeFieldRegister 287 288 """ 289 CHKERR(DMSwarmFinalizeFieldRegister(self.dm)) 290 291 def setLocalSizes(self, nlocal: int, buffer: int) -> Self: 292 """Set the length of all registered fields on the `DMSwarm`. 293 294 Not collective. 295 296 Parameters 297 ---------- 298 nlocal 299 The length of each registered field. 300 buffer 301 The length of the buffer used for efficient dynamic resizing. 302 303 See Also 304 -------- 305 petsc.DMSwarmSetLocalSizes 306 307 """ 308 cdef PetscInt cnlocal = asInt(nlocal) 309 cdef PetscInt cbuffer = asInt(buffer) 310 CHKERR(DMSwarmSetLocalSizes(self.dm, cnlocal, cbuffer)) 311 return self 312 313 def registerField(self, fieldname: str, blocksize: int, dtype: type | dtype = ScalarType) -> None: 314 """Register a field to a `DMSwarm` with a native PETSc data type. 315 316 Collective. 317 318 Parameters 319 ---------- 320 fieldname 321 The textual name to identify this field. 322 blocksize 323 The number of each data type. 324 dtype 325 A valid PETSc data type. 326 327 See Also 328 -------- 329 petsc.DMSwarmRegisterPetscDatatypeField 330 331 """ 332 cdef const char *cfieldname = NULL 333 cdef PetscInt cblocksize = asInt(blocksize) 334 cdef PetscDataType ctype = PETSC_DATATYPE_UNKNOWN 335 if dtype == IntType: ctype = PETSC_INT 336 if dtype == RealType: ctype = PETSC_REAL 337 if dtype == ScalarType: ctype = PETSC_SCALAR 338 if dtype == ComplexType: ctype = PETSC_COMPLEX 339 assert ctype != PETSC_DATATYPE_UNKNOWN 340 fieldname = str2bytes(fieldname, &cfieldname) 341 CHKERR(DMSwarmRegisterPetscDatatypeField(self.dm, cfieldname, cblocksize, ctype)) 342 343 def getField(self, fieldname: str) -> Sequence[int | float | complex]: 344 """Return arrays storing all entries associated with a field. 345 346 Not collective. 347 348 The returned array contains underlying values of the field. 349 350 The array must be returned to the `DMSwarm` using a matching call to 351 `restoreField`. 352 353 Parameters 354 ---------- 355 fieldname 356 The textual name to identify this field. 357 358 Returns 359 ------- 360 `numpy.ndarray` 361 The type of the entries in the array will match the type of the 362 field. The array is two dimensional with shape ``(num_points, blocksize)``. 363 364 See Also 365 -------- 366 restoreField, petsc.DMSwarmGetField 367 368 """ 369 cdef const char *cfieldname = NULL 370 cdef PetscInt blocksize = 0 371 cdef PetscDataType ctype = PETSC_DATATYPE_UNKNOWN 372 cdef PetscReal *data = NULL 373 cdef PetscInt nlocal = 0 374 fieldname = str2bytes(fieldname, &cfieldname) 375 CHKERR(DMSwarmGetField(self.dm, cfieldname, &blocksize, &ctype, <void**> &data)) 376 CHKERR(DMSwarmGetLocalSize(self.dm, &nlocal)) 377 cdef int typenum = -1 378 if ctype == PETSC_INT: typenum = NPY_PETSC_INT 379 if ctype == PETSC_REAL: typenum = NPY_PETSC_REAL 380 if ctype == PETSC_SCALAR: typenum = NPY_PETSC_SCALAR 381 if ctype == PETSC_COMPLEX: typenum = NPY_PETSC_COMPLEX 382 assert typenum != -1 383 cdef npy_intp s[2] 384 s[0] = <npy_intp> nlocal 385 s[1] = <npy_intp> blocksize 386 return <object> PyArray_SimpleNewFromData(2, s, typenum, data) 387 388 def restoreField(self, fieldname: str) -> None: 389 """Restore accesses associated with a registered field. 390 391 Not collective. 392 393 Parameters 394 ---------- 395 fieldname 396 The textual name to identify this field. 397 398 See Also 399 -------- 400 getField, petsc.DMSwarmRestoreField 401 402 """ 403 cdef const char *cfieldname = NULL 404 cdef PetscInt blocksize = 0 405 cdef PetscDataType ctype = PETSC_DATATYPE_UNKNOWN 406 fieldname = str2bytes(fieldname, &cfieldname) 407 CHKERR(DMSwarmRestoreField(self.dm, cfieldname, &blocksize, &ctype, <void**> 0)) 408 409 def vectorDefineField(self, fieldname: str) -> None: 410 """Set the fields from which to define a `Vec` object. 411 412 Collective. 413 414 The field will be used when `DM.createLocalVec`, or 415 `DM.createGlobalVec` is called. 416 417 Parameters 418 ---------- 419 fieldname 420 The textual names given to a registered field. 421 422 See Also 423 -------- 424 petsc.DMSwarmVectorDefineField 425 426 """ 427 cdef const char *cval = NULL 428 fieldname = str2bytes(fieldname, &cval) 429 CHKERR(DMSwarmVectorDefineField(self.dm, cval)) 430 431 def addPoint(self) -> None: 432 """Add space for one new point in the `DMSwarm`. 433 434 Not collective. 435 436 See Also 437 -------- 438 petsc.DMSwarmAddPoint 439 440 """ 441 CHKERR(DMSwarmAddPoint(self.dm)) 442 443 def addNPoints(self, npoints: int) -> None: 444 """Add space for a number of new points in the `DMSwarm`. 445 446 Not collective. 447 448 Parameters 449 ---------- 450 npoints 451 The number of new points to add. 452 453 See Also 454 -------- 455 petsc.DMSwarmAddNPoints 456 457 """ 458 cdef PetscInt cnpoints = asInt(npoints) 459 CHKERR(DMSwarmAddNPoints(self.dm, cnpoints)) 460 461 def removePoint(self) -> None: 462 """Remove the last point from the `DMSwarm`. 463 464 Not collective. 465 466 See Also 467 -------- 468 petsc.DMSwarmRemovePoint 469 470 """ 471 CHKERR(DMSwarmRemovePoint(self.dm)) 472 473 def removePointAtIndex(self, index: int) -> None: 474 """Remove a specific point from the `DMSwarm`. 475 476 Not collective. 477 478 Parameters 479 ---------- 480 index 481 Index of point to remove 482 483 See Also 484 -------- 485 petsc.DMSwarmRemovePointAtIndex 486 487 """ 488 cdef PetscInt cindex = asInt(index) 489 CHKERR(DMSwarmRemovePointAtIndex(self.dm, cindex)) 490 491 def copyPoint(self, pi: int, pj: int) -> None: 492 """Copy point pi to point pj in the `DMSwarm`. 493 494 Not collective. 495 496 Parameters 497 ---------- 498 pi 499 The index of the point to copy (source). 500 pj 501 The point index where the copy should be located (destination). 502 503 See Also 504 -------- 505 petsc.DMSwarmCopyPoint 506 507 """ 508 cdef PetscInt cpi = asInt(pi) 509 cdef PetscInt cpj = asInt(pj) 510 CHKERR(DMSwarmCopyPoint(self.dm, cpi, cpj)) 511 512 def getLocalSize(self) -> int: 513 """Return the local length of fields registered. 514 515 Not collective. 516 517 See Also 518 -------- 519 petsc.DMSwarmGetLocalSize 520 521 """ 522 cdef PetscInt size = asInt(0) 523 CHKERR(DMSwarmGetLocalSize(self.dm, &size)) 524 return toInt(size) 525 526 def getSize(self) -> int: 527 """Return the total length of fields registered. 528 529 Collective. 530 531 See Also 532 -------- 533 petsc.DMSwarmGetSize 534 535 """ 536 cdef PetscInt size = asInt(0) 537 CHKERR(DMSwarmGetSize(self.dm, &size)) 538 return toInt(size) 539 540 def migrate(self, remove_sent_points: bool = False) -> None: 541 """Relocate points defined in the `DMSwarm` to other MPI ranks. 542 543 Collective. 544 545 Parameters 546 ---------- 547 remove_sent_points 548 Flag indicating if sent points should be removed from the current 549 MPI rank. 550 551 See Also 552 -------- 553 petsc.DMSwarmMigrate 554 555 """ 556 cdef PetscBool remove_pts = asBool(remove_sent_points) 557 CHKERR(DMSwarmMigrate(self.dm, remove_pts)) 558 559 def collectViewCreate(self) -> None: 560 """Apply a collection method and gather points in neighbor ranks. 561 562 Collective. 563 564 See Also 565 -------- 566 collectViewDestroy, petsc.DMSwarmCollectViewCreate 567 568 """ 569 CHKERR(DMSwarmCollectViewCreate(self.dm)) 570 571 def collectViewDestroy(self) -> None: 572 """Reset the `DMSwarm` to the size prior to calling `collectViewCreate`. 573 574 Collective. 575 576 See Also 577 -------- 578 collectViewCreate, petsc.DMSwarmCollectViewDestroy 579 580 """ 581 CHKERR(DMSwarmCollectViewDestroy(self.dm)) 582 583 def setCellDM(self, DM dm) -> None: 584 """Attach a `DM` to a `DMSwarm`. 585 586 Collective. 587 588 Parameters 589 ---------- 590 dm 591 The `DM` to attach to the `DMSwarm`. 592 593 See Also 594 -------- 595 getCellDM, petsc.DMSwarmSetCellDM 596 597 """ 598 CHKERR(DMSwarmSetCellDM(self.dm, dm.dm)) 599 600 def getCellDM(self) -> DM: 601 """Return `DM` cell attached to `DMSwarm`. 602 603 Collective. 604 605 See Also 606 -------- 607 setCellDM, petsc.DMSwarmGetCellDM 608 609 """ 610 cdef PetscDM newdm = NULL 611 CHKERR(DMSwarmGetCellDM(self.dm, &newdm)) 612 cdef DM dm = subtype_DM(newdm)() 613 dm.dm = newdm 614 CHKERR(PetscINCREF(dm.obj)) 615 return dm 616 617 def setType(self, dmswarm_type: Type | str) -> None: 618 """Set particular flavor of `DMSwarm`. 619 620 Collective. 621 622 Parameters 623 ---------- 624 dmswarm_type 625 The `DMSwarm` type. 626 627 See Also 628 -------- 629 petsc.DMSwarmSetType 630 631 """ 632 cdef PetscDMSwarmType cval = dmswarm_type 633 CHKERR(DMSwarmSetType(self.dm, cval)) 634 635 def setPointsUniformCoordinates( 636 self, 637 min: Sequence[float], 638 max: Sequence[float], 639 npoints: Sequence[int], 640 mode: InsertMode | None = None) -> Self: 641 """Set point coordinates in a `DMSwarm` on a regular (ijk) grid. 642 643 Collective. 644 645 Parameters 646 ---------- 647 min 648 Minimum coordinate values in the x, y, z directions (array of 649 length ``dim``). 650 max 651 Maximum coordinate values in the x, y, z directions (array of 652 length ``dim``). 653 npoints 654 Number of points in each spatial direction (array of length ``dim``). 655 mode 656 Indicates whether to append points to the swarm (`InsertMode.ADD`), 657 or override existing points (`InsertMode.INSERT`). 658 659 See Also 660 -------- 661 petsc.DMSwarmSetPointsUniformCoordinates 662 663 """ 664 cdef PetscInt dim = asInt(0) 665 CHKERR(DMGetDimension(self.dm, &dim)) 666 cdef PetscReal cmin[3] 667 cmin[0] = cmin[1] = cmin[2] = asReal(0.) 668 for i from 0 <= i < dim: cmin[i] = min[i] 669 cdef PetscReal cmax[3] 670 cmax[0] = cmax[1] = cmax[2] = asReal(0.) 671 for i from 0 <= i < dim: cmax[i] = max[i] 672 cdef PetscInt cnpoints[3] 673 cnpoints[0] = cnpoints[1] = cnpoints[2] = asInt(0) 674 for i from 0 <= i < dim: cnpoints[i] = npoints[i] 675 cdef PetscInsertMode cmode = insertmode(mode) 676 CHKERR(DMSwarmSetPointsUniformCoordinates(self.dm, cmin, cmax, cnpoints, cmode)) 677 return self 678 679 def setPointCoordinates( 680 self, 681 coordinates: Sequence[float], 682 redundant: bool = False, 683 mode: InsertMode | None = None) -> None: 684 """Set point coordinates in a `DMSwarm` from a user-defined list. 685 686 Collective. 687 688 Parameters 689 ---------- 690 coordinates 691 The coordinate values. 692 redundant 693 If set to `True`, it is assumed that coordinates are only valid on 694 rank 0 and should be broadcast to other ranks. 695 mode 696 Indicates whether to append points to the swarm (`InsertMode.ADD`), 697 or override existing points (`InsertMode.INSERT`). 698 699 See Also 700 -------- 701 petsc.DMSwarmSetPointCoordinates 702 703 """ 704 cdef ndarray xyz = iarray(coordinates, NPY_PETSC_REAL) 705 if PyArray_ISFORTRAN(xyz): xyz = PyArray_Copy(xyz) 706 if PyArray_NDIM(xyz) != 2: raise ValueError( 707 ("coordinates must have two dimensions: " 708 "coordinates.ndim=%d") % (PyArray_NDIM(xyz))) 709 cdef PetscInt cnpoints = <PetscInt> PyArray_DIM(xyz, 0) 710 cdef PetscBool credundant = asBool(redundant) 711 cdef PetscInsertMode cmode = insertmode(mode) 712 cdef PetscReal *coords = <PetscReal*> PyArray_DATA(xyz) 713 CHKERR(DMSwarmSetPointCoordinates(self.dm, cnpoints, coords, credundant, cmode)) 714 715 def insertPointUsingCellDM(self, layoutType: PICLayoutType, fill_param: int) -> None: 716 """Insert point coordinates within each cell. 717 718 Not collective. 719 720 Parameters 721 ---------- 722 layout_type 723 Method used to fill each cell with the cell DM. 724 fill_param 725 Parameter controlling how many points per cell are added (the 726 meaning of this parameter is dependent on the layout type). 727 728 See Also 729 -------- 730 petsc.DMSwarmInsertPointsUsingCellDM 731 732 """ 733 cdef PetscDMSwarmPICLayoutType clayoutType = layoutType 734 cdef PetscInt cfill_param = asInt(fill_param) 735 CHKERR(DMSwarmInsertPointsUsingCellDM(self.dm, clayoutType, cfill_param)) 736 737 def setPointCoordinatesCellwise(self, coordinates: Sequence[float]) -> None: 738 """Insert point coordinates within each cell. 739 740 Not collective. 741 742 Point coordinates are defined over the reference cell. 743 744 Parameters 745 ---------- 746 coordinates 747 The coordinates (defined in the local coordinate system for each 748 cell) to insert. 749 750 See Also 751 -------- 752 petsc.DMSwarmSetPointCoordinatesCellwise 753 754 """ 755 cdef ndarray xyz = iarray(coordinates, NPY_PETSC_REAL) 756 if PyArray_ISFORTRAN(xyz): xyz = PyArray_Copy(xyz) 757 if PyArray_NDIM(xyz) != 2: raise ValueError( 758 ("coordinates must have two dimensions: " 759 "coordinates.ndim=%d") % (PyArray_NDIM(xyz))) 760 cdef PetscInt cnpoints = <PetscInt> PyArray_DIM(xyz, 0) 761 cdef PetscReal *coords = <PetscReal*> PyArray_DATA(xyz) 762 CHKERR(DMSwarmSetPointCoordinatesCellwise(self.dm, cnpoints, coords)) 763 764 def viewFieldsXDMF(self, filename: str, fieldnames: Sequence[str]) -> None: 765 """Write a selection of `DMSwarm` fields to an XDMF3 file. 766 767 Collective. 768 769 Parameters 770 ---------- 771 filename 772 The file name of the XDMF file (must have the extension .xmf). 773 fieldnames 774 Array containing the textual names of fields to write. 775 776 See Also 777 -------- 778 petsc.DMSwarmViewFieldsXDMF 779 780 """ 781 cdef const char *cval = NULL 782 cdef const char *cfilename = NULL 783 filename = str2bytes(filename, &cfilename) 784 cdef PetscInt cnfields = <PetscInt> len(fieldnames) 785 cdef const char** cfieldnames = NULL 786 cdef object unused = oarray_p(empty_p(cnfields), NULL, <void**>&cfieldnames) 787 fieldnames = list(fieldnames) 788 for i from 0 <= i < cnfields: 789 fieldnames[i] = str2bytes(fieldnames[i], &cval) 790 cfieldnames[i] = cval 791 CHKERR(DMSwarmViewFieldsXDMF(self.dm, cfilename, cnfields, cfieldnames)) 792 793 def viewXDMF(self, filename: str) -> None: 794 """Write this `DMSwarm` fields to an XDMF3 file. 795 796 Collective. 797 798 Parameters 799 ---------- 800 filename 801 The file name of the XDMF file (must have the extension .xmf). 802 803 See Also 804 -------- 805 petsc.DMSwarmViewXDMF 806 807 """ 808 cdef const char *cval = NULL 809 filename = str2bytes(filename, &cval) 810 CHKERR(DMSwarmViewXDMF(self.dm, cval)) 811 812 def sortGetAccess(self) -> None: 813 """Setup up a `DMSwarm` point sort context. 814 815 Not collective. 816 817 The point sort context is used for efficient traversal of points within 818 a cell. 819 820 You must call `sortRestoreAccess` when you no longer need access to the 821 sort context. 822 823 See Also 824 -------- 825 sortRestoreAccess, petsc.DMSwarmSortGetAccess 826 827 """ 828 CHKERR(DMSwarmSortGetAccess(self.dm)) 829 830 def sortRestoreAccess(self) -> None: 831 """Invalidate the `DMSwarm` point sorting context. 832 833 Not collective. 834 835 See Also 836 -------- 837 sortGetAccess, petsc.DMSwarmSortRestoreAccess 838 839 """ 840 CHKERR(DMSwarmSortRestoreAccess(self.dm)) 841 842 def sortGetPointsPerCell(self, e: int) -> list[int]: 843 """Create an array of point indices for all points in a cell. 844 845 Not collective. 846 847 Parameters 848 ---------- 849 e 850 The index of the cell. 851 852 See Also 853 -------- 854 petsc.DMSwarmSortGetPointsPerCell 855 856 """ 857 cdef PetscInt ce = asInt(e) 858 cdef PetscInt cnpoints = asInt(0) 859 cdef PetscInt *cpidlist = NULL 860 cdef list pidlist = [] 861 CHKERR(DMSwarmSortGetPointsPerCell(self.dm, ce, &cnpoints, &cpidlist)) 862 npoints = asInt(cnpoints) 863 for i from 0 <= i < npoints: pidlist.append(asInt(cpidlist[i])) 864 return pidlist 865 866 def sortGetNumberOfPointsPerCell(self, e: int) -> int: 867 """Return the number of points in a cell. 868 869 Not collective. 870 871 Parameters 872 ---------- 873 e 874 The index of the cell. 875 876 See Also 877 -------- 878 petsc.DMSwarmSortGetNumberOfPointsPerCell 879 880 """ 881 cdef PetscInt ce = asInt(e) 882 cdef PetscInt npoints = asInt(0) 883 CHKERR(DMSwarmSortGetNumberOfPointsPerCell(self.dm, ce, &npoints)) 884 return toInt(npoints) 885 886 def sortGetIsValid(self) -> bool: 887 """Return whether the sort context is up-to-date. 888 889 Not collective. 890 891 Returns the flag associated with a `DMSwarm` point sorting context. 892 893 See Also 894 -------- 895 petsc.DMSwarmSortGetIsValid 896 897 """ 898 cdef PetscBool isValid = asBool(False) 899 CHKERR(DMSwarmSortGetIsValid(self.dm, &isValid)) 900 return toBool(isValid) 901 902 def sortGetSizes(self) -> tuple[int, int]: 903 """Return the sizes associated with a `DMSwarm` point sorting context. 904 905 Not collective. 906 907 Returns 908 ------- 909 ncells : int 910 Number of cells within the sort context. 911 npoints : int 912 Number of points used to create the sort context. 913 914 See Also 915 -------- 916 petsc.DMSwarmSortGetSizes 917 918 """ 919 cdef PetscInt ncells = asInt(0) 920 cdef PetscInt npoints = asInt(0) 921 CHKERR(DMSwarmSortGetSizes(self.dm, &ncells, &npoints)) 922 return (toInt(ncells), toInt(npoints)) 923 924 def projectFields(self, DM dm, fieldnames: Sequence[str], vecs: Sequence[Vec], mode: ScatterModeSpec = None) -> None: 925 """Project a set of `DMSwarm` fields onto the cell `DM`. 926 927 Collective. 928 929 Parameters 930 ---------- 931 dm 932 The continuum `DM`. 933 fieldnames 934 The textual names of the swarm fields to project. 935 936 See Also 937 -------- 938 petsc.DMSwarmProjectFields 939 940 """ 941 cdef const char *cval = NULL 942 cdef PetscInt cnfields = <PetscInt> len(fieldnames) 943 cdef const char** cfieldnames = NULL 944 cdef object unused = oarray_p(empty_p(cnfields), NULL, <void**>&cfieldnames) 945 cdef PetscVec *cfieldvecs = NULL 946 cdef object unused2 = oarray_p(empty_p(cnfields), NULL, <void**>&cfieldvecs) 947 cdef PetscScatterMode cmode = scattermode(mode) 948 fieldnames = list(fieldnames) 949 for i from 0 <= i < cnfields: 950 fieldnames[i] = str2bytes(fieldnames[i], &cval) 951 cfieldnames[i] = cval 952 cfieldvecs[i] = (<Vec?>(vecs[i])).vec 953 CHKERR(DMSwarmProjectFields(self.dm, dm.dm, cnfields, cfieldnames, cfieldvecs, cmode)) 954 return 955 956 def addCellDM(self, CellDM celldm) -> None: 957 """Add a cell DM to the `DMSwarm`. 958 959 Logically collective. 960 961 Parameters 962 ------- 963 cellDM : CellDM 964 The cell DM object 965 966 See Also 967 -------- 968 petsc.DMSwarmAddCellDM 969 970 """ 971 CHKERR(DMSwarmAddCellDM(self.dm, celldm.cdm)) 972 973 def setCellDMActive(self, name : str) -> None: 974 """Activate a cell DM in the `DMSwarm`. 975 976 Logically collective. 977 978 Parameters 979 ------- 980 name : str 981 The name of the cell DM to activate 982 983 See Also 984 -------- 985 petsc.DMSwarmSetCellDMActive 986 987 """ 988 cdef const char *cname = NULL 989 str2bytes(name, &cname) 990 CHKERR(DMSwarmSetCellDMActive(self.dm, cname)) 991 992 def getCellDMActive(self) -> CellDM: 993 """Return the active cell DM in the `DMSwarm`. 994 995 Not collective. 996 997 See Also 998 -------- 999 petsc.DMSwarmGetCellDMActive 1000 1001 """ 1002 cdef PetscDMSwarmCellDM newcdm = NULL 1003 CHKERR(DMSwarmGetCellDMActive(self.dm, &newcdm)) 1004 cdef CellDM cdm = CellDM() 1005 cdm.cdm = newcdm 1006 CHKERR(PetscINCREF(cdm.obj)) 1007 return cdm 1008 1009 def getCellDMByName(self, name: str) -> CellDM: 1010 """Return the cell DM with the given name in the `DMSwarm`. 1011 1012 Not collective. 1013 1014 See Also 1015 -------- 1016 petsc.DMSwarmGetCellDMByName 1017 1018 """ 1019 cdef PetscDMSwarmCellDM newcdm = NULL 1020 cdef const char *cname = NULL 1021 str2bytes(name, &cname) 1022 CHKERR(DMSwarmGetCellDMByName(self.dm, cname, &newcdm)) 1023 cdef CellDM cdm = CellDM() 1024 cdm.cdm = newcdm 1025 CHKERR(PetscINCREF(cdm.obj)) 1026 return cdm 1027 1028 def getCellDMNames(self) -> list[str]: 1029 """Return the names of all cell DMs in the `DMSwarm`. 1030 1031 Not collective. 1032 1033 See Also 1034 -------- 1035 petsc.DMSwarmGetCellDMNames 1036 1037 """ 1038 cdef PetscInt i = 0, nc = 0 1039 cdef const char **c = NULL 1040 CHKERR(DMSwarmGetCellDMNames(self.dm, &nc, &c)) 1041 cdef list names = [] 1042 for i from 0 <= i < nc: 1043 names.append(bytes2str(c[i])) 1044 return names 1045 1046 def computeMoments(self, coord: str, weight: str) -> list[float]: 1047 """Return the moments defined in the active cell DM. 1048 1049 Collective. 1050 1051 Notes 1052 ----- 1053 We integrate the given weight field over the given coordinate 1054 1055 See Also 1056 -------- 1057 petsc.DMSwarmComputeMoments 1058 1059 """ 1060 cdef PetscReal *mom = NULL 1061 cdef PetscInt cbs = 0 1062 cdef const char *ccoord = NULL 1063 cdef const char *cweight = NULL 1064 str2bytes(coord, &ccoord) 1065 str2bytes(weight, &cweight) 1066 CHKERR(DMSwarmGetFieldInfo(self.dm, ccoord, &cbs, NULL)) 1067 cdef object moments = oarray_r(empty_r(asInt(cbs) + 2), NULL, &mom) 1068 CHKERR(DMSwarmComputeMoments(self.dm, ccoord, cweight, mom)) 1069 return moments 1070 1071# -------------------------------------------------------------------- 1072 1073cdef class CellDM(Object): 1074 """CellDM object. 1075 1076 See Also 1077 -------- 1078 petsc.DMSwarmCellDM 1079 1080 """ 1081 # 1082 1083 def __cinit__(self): 1084 self.obj = <PetscObject*> &self.cdm 1085 self.cdm = NULL 1086 1087 # 1088 1089 def view(self, Viewer viewer=None) -> None: 1090 """View the cell DM. 1091 1092 Collective. 1093 1094 Parameters 1095 ---------- 1096 viewer 1097 A `Viewer` instance or `None` for the default viewer. 1098 1099 See Also 1100 -------- 1101 Viewer, petsc.DMSwarmCellDMView 1102 1103 """ 1104 cdef PetscViewer vwr = NULL 1105 if viewer is not None: vwr = viewer.vwr 1106 CHKERR(DMSwarmCellDMView(self.cdm, vwr)) 1107 1108 def destroy(self) -> Self: 1109 """Destroy the cell DM. 1110 1111 Collective. 1112 1113 See Also 1114 -------- 1115 create, petsc.DMSwarmCellDMDestroy 1116 1117 """ 1118 CHKERR(DMSwarmCellDMDestroy(&self.cdm)) 1119 return self 1120 1121 def create( 1122 self, 1123 DM dm, 1124 fields: Sequence[str], 1125 coords: Sequence[str]) -> Self: 1126 """Create the cell DM. 1127 1128 Collective. 1129 1130 Parameters 1131 ---------- 1132 dm 1133 The underlying DM on which to place particles 1134 fields 1135 The swarm fields represented on the DM 1136 coords 1137 The swarm fields to use as coordinates in the DM 1138 1139 See Also 1140 -------- 1141 destroy, petsc.DMSwarmCellDMCreate 1142 1143 """ 1144 cdef const char *cval = NULL 1145 cdef PetscInt i = 0, nf = <PetscInt> len(fields) 1146 cdef const char** fieldnames = NULL 1147 cdef object unuseda = oarray_p(empty_p(nf), NULL, <void**>&fieldnames) 1148 fields = list(fields) 1149 for i from 0 <= i < nf: 1150 str2bytes(fields[i], &cval) 1151 fieldnames[i] = cval 1152 cdef PetscInt nc = <PetscInt> len(coords) 1153 cdef const char** coordnames = NULL 1154 cdef object unusedb = oarray_p(empty_p(nf), NULL, <void**>&coordnames) 1155 coords = list(coords) 1156 for i from 0 <= i < nc: 1157 str2bytes(coords[i], &cval) 1158 coordnames[i] = cval 1159 cdef PetscDMSwarmCellDM newcdm = NULL 1160 CHKERR(DMSwarmCellDMCreate(dm.dm, nf, fieldnames, nc, coordnames, &newcdm)) 1161 CHKERR(PetscCLEAR(self.obj)); self.cdm = newcdm 1162 return self 1163 1164 def getDM(self) -> DM: 1165 """Return the underlying DM. 1166 1167 Not collective. 1168 1169 See Also 1170 -------- 1171 petsc.DMSwarmCellDMGetDM 1172 1173 """ 1174 cdef PetscDM newdm = NULL 1175 CHKERR(DMSwarmCellDMGetDM(self.cdm, &newdm)) 1176 cdef DM dm = subtype_DM(newdm)() 1177 dm.dm = newdm 1178 CHKERR(PetscINCREF(dm.obj)) 1179 return dm 1180 1181 def getCellID(self) -> str: 1182 """Return the cellid field for this cell DM. 1183 1184 Not collective. 1185 1186 See Also 1187 -------- 1188 petsc.DMSwarmCellDMGetCellID 1189 1190 """ 1191 cdef const char *cellid = NULL 1192 CHKERR(DMSwarmCellDMGetCellID(self.cdm, &cellid)) 1193 return bytes2str(cellid) 1194 1195 def getBlockSize(self, DM sw) -> int: 1196 """Return the block size for this cell DM. 1197 1198 Not collective. 1199 1200 Parameters 1201 ---------- 1202 sw: 1203 The `DMSwarm` object 1204 1205 See Also 1206 -------- 1207 petsc.DMSwarmCellDMGetBlockSize 1208 1209 """ 1210 cdef PetscInt bs = 0 1211 CHKERR(DMSwarmCellDMGetBlockSize(self.cdm, sw.dm, &bs)) 1212 return toInt(bs) 1213 1214 def getFields(self) -> list[str]: 1215 """Return the swarm fields defined on this cell DM. 1216 1217 Not collective. 1218 1219 See Also 1220 -------- 1221 petsc.DMSwarmCellDMGetFields 1222 1223 """ 1224 cdef PetscInt i = 0, nf = 0 1225 cdef const char **f = NULL 1226 CHKERR(DMSwarmCellDMGetFields(self.cdm, &nf, &f)) 1227 cdef list fields = [] 1228 for i from 0 <= i < nf: 1229 fields.append(bytes2str(f[i])) 1230 return fields 1231 1232 def getCoordinateFields(self) -> list[str]: 1233 """Return the swarm fields used as coordinates on this cell DM. 1234 1235 Not collective. 1236 1237 See Also 1238 -------- 1239 petsc.DMSwarmCellDMGetCoordinateFields 1240 1241 """ 1242 cdef PetscInt i = 0, ncf = 0 1243 cdef const char **cf = NULL 1244 CHKERR(DMSwarmCellDMGetCoordinateFields(self.cdm, &ncf, &cf)) 1245 cdef list fields = [] 1246 for i from 0 <= i < ncf: 1247 fields.append(bytes2str(cf[i])) 1248 return fields 1249 1250# -------------------------------------------------------------------- 1251 1252del DMSwarmType 1253del DMSwarmMigrateType 1254del DMSwarmCollectType 1255del DMSwarmPICLayoutType 1256