1# -------------------------------------------------------------------- 2 3cdef class DMPlex(DM): 4 """Encapsulate an unstructured mesh. 5 6 DMPlex encapsulates both topology and geometry. 7 It is capable of parallel refinement and coarsening (using Pragmatic or ParMmg) 8 and parallel redistribution for load balancing. 9 10 """ 11 12 # 13 14 def create(self, comm: Comm | None = None) -> Self: 15 """Create a `DMPlex` object, which encapsulates an unstructured mesh. 16 17 Collective. 18 19 Parameters 20 ---------- 21 comm 22 MPI communicator, defaults to `Sys.getDefaultComm`. 23 24 See Also 25 -------- 26 DM, DMPlex, DM.create, DM.setType, petsc.DMPlexCreate 27 28 """ 29 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 30 cdef PetscDM newdm = NULL 31 CHKERR(DMPlexCreate(ccomm, &newdm)) 32 CHKERR(PetscCLEAR(self.obj)); self.dm = newdm 33 return self 34 35 def createFromCellList(self, dim: int, cells: Sequence[int], coords: Sequence[float], interpolate: bool | None = True, comm: Comm | None = None) -> Self: 36 """Create a `DMPlex` from a list of vertices for each cell on process 0. 37 38 Collective. 39 40 Parameters 41 ---------- 42 dim 43 The topological dimension of the mesh. 44 cells 45 An array of number of cells times number of vertices on each cell. 46 coords 47 An array of number of vertices times spatial dimension for coordinates. 48 interpolate 49 Flag to interpolate the mesh. 50 comm 51 MPI communicator, defaults to `Sys.getDefaultComm`. 52 53 See Also 54 -------- 55 DM, DMPlex, DMPlex.create, DMPlex.interpolate, 56 petsc.DMPlexCreateFromCellListPetsc 57 58 """ 59 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 60 cdef PetscBool interp = interpolate 61 cdef PetscDM newdm = NULL 62 cdef PetscInt cdim = asInt(dim) 63 cdef PetscInt numCells = 0 64 cdef PetscInt numCorners = 0 65 cdef PetscInt *cellVertices = NULL 66 cdef PetscInt numVertices = 0 67 cdef PetscInt spaceDim= 0 68 cdef PetscReal *vertexCoords = NULL 69 cdef int npy_flags = NPY_ARRAY_ALIGNED|NPY_ARRAY_NOTSWAPPED|NPY_ARRAY_CARRAY 70 cells = PyArray_FROM_OTF(cells, NPY_PETSC_INT, npy_flags) 71 coords = PyArray_FROM_OTF(coords, NPY_PETSC_REAL, npy_flags) 72 if PyArray_NDIM(cells) != 2: raise ValueError( 73 ("cell indices must have two dimensions: " 74 "cells.ndim=%d") % (PyArray_NDIM(cells))) 75 if PyArray_NDIM(coords) != 2: raise ValueError( 76 ("coords vertices must have two dimensions: " 77 "coords.ndim=%d") % (PyArray_NDIM(coords))) 78 numCells = <PetscInt> PyArray_DIM(cells, 0) 79 numCorners = <PetscInt> PyArray_DIM(cells, 1) 80 numVertices = <PetscInt> PyArray_DIM(coords, 0) 81 spaceDim = <PetscInt> PyArray_DIM(coords, 1) 82 cellVertices = <PetscInt*> PyArray_DATA(cells) 83 vertexCoords = <PetscReal*> PyArray_DATA(coords) 84 CHKERR(DMPlexCreateFromCellListPetsc(ccomm, cdim, numCells, numVertices, 85 numCorners, interp, cellVertices, 86 spaceDim, vertexCoords, &newdm)) 87 CHKERR(PetscCLEAR(self.obj)); self.dm = newdm 88 return self 89 90 def createBoxMesh(self, faces: Sequence[int], lower: Sequence[float] | None = (0, 0, 0), upper: Sequence[float] | None = (1, 1, 1), 91 simplex: bool | None = True, periodic: Sequence | str | int | bool | None = False, interpolate: bool | None = True, localizationHeight: int | None = 0, sparseLocalize: bool | None = True, comm: Comm | None = None) -> Self: 92 """Create a mesh on the tensor product of intervals. 93 94 Collective. 95 96 Parameters 97 ---------- 98 faces 99 Number of faces per dimension, or `None` for the default. 100 lower 101 The lower left corner. 102 upper 103 The upper right corner. 104 simplex 105 `True` for simplices, `False` for tensor cells. 106 periodic 107 The boundary type for the X, Y, Z direction, 108 or `None` for `DM.BoundaryType.NONE`. 109 interpolate 110 Flag to create intermediate mesh entities (edges, faces). 111 localizationHeight 112 Flag to localize edges and faces in addition to cells; 113 only significant for periodic meshes. 114 sparseLocalize 115 Flag to localize coordinates only for cells near the 116 periodic boundary; only significant for periodic meshes. 117 comm 118 MPI communicator, defaults to `Sys.getDefaultComm`. 119 120 See Also 121 -------- 122 DM, DMPlex, DM.setFromOptions, DMPlex.createFromFile, DM.setType 123 DM.create, petsc.DMPlexCreateBoxMesh 124 125 """ 126 cdef Py_ssize_t i = 0 127 cdef PetscInt dim = 0, *cfaces = NULL 128 faces = iarray_i(faces, &dim, &cfaces) 129 assert dim >= 1 and dim <= 3 130 cdef PetscReal clower[3] 131 clower[0] = clower[1] = clower[2] = 0 132 for i from 0 <= i < dim: clower[i] = lower[i] 133 cdef PetscReal cupper[3] 134 cupper[0] = cupper[1] = cupper[2] = 1 135 for i from 0 <= i < dim: cupper[i] = upper[i] 136 cdef PetscDMBoundaryType btype[3] 137 asBoundary(periodic, &btype[0], &btype[1], &btype[2]) 138 cdef PetscBool csimplex = simplex 139 cdef PetscBool cinterp = interpolate 140 cdef PetscInt clocalizationHeight = asInt(localizationHeight) 141 cdef PetscBool csparseLocalize = sparseLocalize 142 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 143 cdef PetscDM newdm = NULL 144 CHKERR(DMPlexCreateBoxMesh(ccomm, dim, csimplex, cfaces, 145 clower, cupper, btype, cinterp, clocalizationHeight, csparseLocalize, &newdm)) 146 CHKERR(PetscCLEAR(self.obj)); self.dm = newdm 147 return self 148 149 def createBoxSurfaceMesh(self, faces: Sequence[int], lower: Sequence[float] | None = (0, 0, 0), upper: Sequence[float] | None = (1, 1, 1), 150 interpolate: bool | None = True, comm: Comm | None = None) -> Self: 151 """Create a mesh on the surface of a box mesh using tensor cells. 152 153 Collective. 154 155 Parameters 156 ---------- 157 faces 158 Number of faces per dimension, or `None` for the default. 159 lower 160 The lower left corner. 161 upper 162 The upper right corner. 163 interpolate 164 Flag to create intermediate mesh pieces (edges, faces). 165 comm 166 MPI communicator, defaults to `Sys.getDefaultComm`. 167 168 See Also 169 -------- 170 DM, DMPlex, DM.setFromOptions, DMPlex.createBoxMesh 171 DMPlex.createFromFile, DM.setType, DM.create 172 petsc.DMPlexCreateBoxSurfaceMesh 173 174 """ 175 cdef Py_ssize_t i = 0 176 cdef PetscInt dim = 0, *cfaces = NULL 177 faces = iarray_i(faces, &dim, &cfaces) 178 assert dim >= 1 and dim <= 3 179 cdef PetscReal clower[3] 180 clower[0] = clower[1] = clower[2] = 0 181 for i from 0 <= i < dim: clower[i] = lower[i] 182 cdef PetscReal cupper[3] 183 cupper[0] = cupper[1] = cupper[2] = 1 184 for i from 0 <= i < dim: cupper[i] = upper[i] 185 cdef PetscBool cinterp = interpolate 186 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 187 cdef PetscDM newdm = NULL 188 CHKERR(DMPlexCreateBoxSurfaceMesh(ccomm, dim, cfaces, clower, cupper, cinterp, &newdm)) 189 CHKERR(PetscCLEAR(self.obj)); self.dm = newdm 190 return self 191 192 def createFromFile(self, filename: str, plexname: str | None = "unnamed", interpolate: bool | None = True, comm: Comm | None = None) -> Self: 193 """Create `DMPlex` from a file. 194 195 Collective. 196 197 Parameters 198 ---------- 199 filename 200 A file name. 201 plexname 202 The name of the resulting `DMPlex`, 203 also used for intra-datafile lookup by some formats. 204 interpolate 205 Flag to create intermediate mesh pieces (edges, faces). 206 comm 207 MPI communicator, defaults to `Sys.getDefaultComm`. 208 209 See Also 210 -------- 211 DM, DMPlex, DMPlex.createFromCellList, DMPlex.create, Object.setName 212 DM.view, DM.load, petsc_options, petsc.DMPlexCreateFromFile 213 214 """ 215 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 216 cdef PetscBool interp = interpolate 217 cdef PetscDM newdm = NULL 218 cdef const char *cfile = NULL 219 cdef const char *pname = NULL 220 filename = str2bytes(filename, &cfile) 221 plexname = str2bytes(plexname, &pname) 222 CHKERR(DMPlexCreateFromFile(ccomm, cfile, pname, interp, &newdm)) 223 CHKERR(PetscCLEAR(self.obj)); self.dm = newdm 224 return self 225 226 def createCGNS(self, cgid: int, interpolate: bool | None = True, comm: Comm | None = None) -> Self: 227 """Create a `DMPlex` mesh from a CGNS file. 228 229 Collective. 230 231 Parameters 232 ---------- 233 cgid 234 The CG id associated with a file and obtained using cg_open. 235 interpolate 236 Create faces and edges in the mesh. 237 comm 238 MPI communicator, defaults to `Sys.getDefaultComm`. 239 240 See Also 241 -------- 242 DM, DMPlex, DMPlex.create, DMPlex.createCGNSFromFile 243 DMPlex.createExodus, petsc.DMPlexCreateCGNS 244 245 """ 246 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 247 cdef PetscBool interp = interpolate 248 cdef PetscDM newdm = NULL 249 cdef PetscInt ccgid = asInt(cgid) 250 CHKERR(DMPlexCreateCGNS(ccomm, ccgid, interp, &newdm)) 251 CHKERR(PetscCLEAR(self.obj)); self.dm = newdm 252 return self 253 254 def createCGNSFromFile(self, filename: str, interpolate: bool | None = True, comm: Comm | None = None) -> Self: 255 """"Create a `DMPlex` mesh from a CGNS file. 256 257 Collective. 258 259 Parameters 260 ---------- 261 filename 262 The name of the CGNS file. 263 interpolate 264 Create faces and edges in the mesh. 265 comm 266 MPI communicator, defaults to `Sys.getDefaultComm`. 267 268 See Also 269 -------- 270 DM, DMPlex, DMPlex.create, DMPlex.createCGNS, DMPlex.createExodus 271 petsc.DMPlexCreateCGNS 272 273 """ 274 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 275 cdef PetscBool interp = interpolate 276 cdef PetscDM newdm = NULL 277 cdef const char *cfile = NULL 278 filename = str2bytes(filename, &cfile) 279 CHKERR(DMPlexCreateCGNSFromFile(ccomm, cfile, interp, &newdm)) 280 CHKERR(PetscCLEAR(self.obj)); self.dm = newdm 281 return self 282 283 def createExodusFromFile(self, filename: str, interpolate: bool | None = True, comm: Comm | None = None) -> Self: 284 """Create a `DMPlex` mesh from an ExodusII file. 285 286 Collective. 287 288 Parameters 289 ---------- 290 filename 291 The name of the ExodusII file. 292 interpolate 293 Create faces and edges in the mesh. 294 comm 295 MPI communicator, defaults to `Sys.getDefaultComm`. 296 297 See Also 298 -------- 299 DM, DMPlex, DM.create, DMPlex.createExodus 300 petsc.DMPlexCreateExodusFromFile 301 302 """ 303 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 304 cdef PetscBool interp = interpolate 305 cdef PetscDM newdm = NULL 306 cdef const char *cfile = NULL 307 filename = str2bytes(filename, &cfile) 308 CHKERR(DMPlexCreateExodusFromFile(ccomm, cfile, interp, &newdm)) 309 CHKERR(PetscCLEAR(self.obj)); self.dm = newdm 310 return self 311 312 def createExodus(self, exoid: int, interpolate: bool | None = True, comm: Comm | None = None) -> Self: 313 """Create a `DMPlex` mesh from an ExodusII file ID. 314 315 Collective. 316 317 Parameters 318 ---------- 319 exoid 320 The ExodusII id associated with a file obtained using ``ex_open``. 321 interpolate 322 Create faces and edges in the mesh, 323 comm 324 MPI communicator, defaults to `Sys.getDefaultComm`. 325 326 See Also 327 -------- 328 DM, DMPlex, DM.create, petsc.DMPlexCreateExodus 329 330 """ 331 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 332 cdef PetscBool interp = interpolate 333 cdef PetscDM newdm = NULL 334 cdef PetscInt cexoid = asInt(exoid) 335 CHKERR(DMPlexCreateExodus(ccomm, <int> cexoid, interp, &newdm)) 336 CHKERR(PetscCLEAR(self.obj)); self.dm = newdm 337 return self 338 339 def createGmsh(self, Viewer viewer, interpolate: bool | None = True, comm: Comm | None = None) -> Self: 340 """Create a `DMPlex` mesh from a Gmsh file viewer. 341 342 Collective. 343 344 Parameters 345 ---------- 346 viewer 347 The `Viewer` associated with a Gmsh file. 348 interpolate 349 Create faces and edges in the mesh. 350 comm 351 MPI communicator, defaults to `Sys.getDefaultComm`. 352 353 Notes 354 ----- 355 ``-dm_plex_gmsh_hybrid`` forces triangular prisms to use tensor order.\n 356 ``-dm_plex_gmsh_periodic`` allows for reading Gmsh periodic section.\n 357 ``-dm_plex_gmsh_highorder`` allows for generating high-order coordinates.\n 358 ``-dm_plex_gmsh_project`` projects high-order coordinates to a different space, 359 use the prefix ``-dm_plex_gmsh_project_`` to define the space.\n 360 ``-dm_plex_gmsh_use_regions`` generates labels with region names.\n 361 ``-dm_plex_gmsh_mark_vertices`` adds vertices to generated labels.\n 362 ``-dm_plex_gmsh_multiple_tags`` allows multiple tags for default labels.\n 363 ``-dm_plex_gmsh_spacedim <d>`` embedding space dimension. 364 365 See Also 366 -------- 367 DM, DMPlex, DM.create, petsc_options, petsc.DMPlexCreateGmsh 368 369 """ 370 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 371 cdef PetscBool interp = interpolate 372 cdef PetscDM newdm = NULL 373 CHKERR(DMPlexCreateGmsh(ccomm, viewer.vwr, interp, &newdm)) 374 CHKERR(PetscCLEAR(self.obj)); self.dm = newdm 375 return self 376 377 def createCoordinateSpace(self, degree: int, localized: bool, project: bool) -> None: 378 """Create a finite element space for the coordinates. 379 380 Collective. 381 382 Parameters 383 ---------- 384 degree 385 The degree of the finite element. 386 localized 387 Flag to create a localized (DG) coordinate space. 388 project 389 Flag to project current coordinates into the space. 390 391 See Also 392 -------- 393 DM, DMPlex, petsc.DMPlexCreateCoordinateSpace 394 395 """ 396 cdef PetscInt cdegree = asInt(degree) 397 cdef PetscBool clocalized = localized 398 cdef PetscBool cproject = project 399 CHKERR(DMPlexCreateCoordinateSpace(self.dm, cdegree, clocalized, cproject)) 400 401 def createCohesiveSubmesh(self, hasLagrange: bool, value: int) -> DMPlex: 402 """Extract the hypersurface defined by one face of the cohesive cells. 403 404 Collective. 405 406 Parameters 407 ---------- 408 hasLagrange 409 Flag indicating whether the mesh has Lagrange dofs in the cohesive cells. 410 value 411 A label value. 412 413 See Also 414 -------- 415 DM, DMPlex, petsc.DMPlexCreateCohesiveSubmesh 416 417 """ 418 cdef PetscBool flag = hasLagrange 419 cdef PetscInt cvalue = asInt(value) 420 cdef DM subdm = DMPlex() 421 CHKERR(DMPlexCreateCohesiveSubmesh(self.dm, flag, NULL, cvalue, &subdm.dm)) 422 return subdm 423 424 def filter(self, label: DMLabel | None = None, value: int | None = None, ignoreHalo: bool = False, 425 sanitizeSubMesh: bool = False, comm: Comm | None = None) -> tuple[DMPlex, SF]: 426 """Extract a subset of mesh cells defined by a label as a separate mesh. 427 428 Collective. 429 430 Parameters 431 ---------- 432 cellLabel 433 label marking cells to be contained in the new mesh 434 value 435 label value to use 436 ignoreHalo 437 Flag indicating if labeled points that are in the halo are ignored 438 sanitizeSubmesh 439 Flag indicating if a subpoint is forced to be owned by a rank that owns 440 a subcell that contains that point in its closure 441 comm 442 The communicator you want the mesh on 443 444 See Also 445 -------- 446 DM, DMPlex, petsc.DMPlexFilter 447 448 """ 449 cdef DM subdm = DMPlex() 450 cdef PetscDMLabel clbl = NULL 451 cdef MPI_Comm ccomm = MPI_COMM_NULL 452 cdef MPI_Comm dmcomm = MPI_COMM_NULL 453 cdef PetscInt cvalue = -1 454 cdef PetscBool cignoreHalo = ignoreHalo 455 cdef PetscBool csanitize = sanitizeSubMesh 456 cdef SF sf = SF() 457 if label is not None: 458 clbl = (<DMLabel?>label).dmlabel 459 if value is not None: 460 cvalue = asInt(value) 461 CHKERR(PetscObjectGetComm(<PetscObject>self.dm, &dmcomm)) 462 ccomm = def_Comm(comm, dmcomm) 463 CHKERR(DMPlexFilter(self.dm, clbl, cvalue, cignoreHalo, csanitize, ccomm, &sf.sf, &subdm.dm)) 464 return subdm, sf 465 466 def getChart(self) -> tuple[int, int]: 467 """Return the interval for all mesh points [``pStart``, ``pEnd``). 468 469 Not collective. 470 471 Returns 472 ------- 473 pStart : int 474 The first mesh point. 475 pEnd : int 476 The upper bound for mesh points. 477 478 See Also 479 -------- 480 DM, DMPlex, DMPlex.create, DMPlex.setChart, petsc.DMPlexGetChart 481 482 """ 483 cdef PetscInt pStart = 0, pEnd = 0 484 CHKERR(DMPlexGetChart(self.dm, &pStart, &pEnd)) 485 return toInt(pStart), toInt(pEnd) 486 487 def setChart(self, pStart: int, pEnd: int) -> None: 488 """Set the interval for all mesh points [``pStart``, ``pEnd``). 489 490 Not collective. 491 492 Parameters 493 ---------- 494 pStart 495 The first mesh point. 496 pEnd 497 The upper bound for mesh points. 498 499 See Also 500 -------- 501 DM, DMPlex, DMPlex.create, DMPlex.getChart, petsc.DMPlexSetChart 502 503 """ 504 cdef PetscInt cStart = asInt(pStart) 505 cdef PetscInt cEnd = asInt(pEnd) 506 CHKERR(DMPlexSetChart(self.dm, cStart, cEnd)) 507 508 def getConeSize(self, p: int) -> int: 509 """Return the number of in-edges for this point in the DAG. 510 511 Not collective. 512 513 Parameters 514 ---------- 515 p 516 The point, which must lie in the chart set with `DMPlex.setChart`. 517 518 See Also 519 -------- 520 DM, DMPlex, DMPlex.create, DMPlex.setConeSize, DMPlex.setChart 521 petsc.DMPlexGetConeSize 522 523 """ 524 cdef PetscInt cp = asInt(p) 525 cdef PetscInt pStart = 0, pEnd = 0 526 CHKERR(DMPlexGetChart(self.dm, &pStart, &pEnd)) 527 assert cp>=pStart and cp<pEnd 528 cdef PetscInt csize = 0 529 CHKERR(DMPlexGetConeSize(self.dm, cp, &csize)) 530 return toInt(csize) 531 532 def setConeSize(self, p: int, size: int) -> None: 533 """Set the number of in-edges for this point in the DAG. 534 535 Not collective. 536 537 Parameters 538 ---------- 539 p 540 The point, which must lie in the chart set with `DMPlex.setChart`. 541 size 542 The cone size for point ``p``. 543 544 See Also 545 -------- 546 DM, DMPlex, DMPlex.create, DMPlex.getConeSize, DMPlex.setChart 547 petsc.DMPlexSetConeSize 548 549 """ 550 cdef PetscInt cp = asInt(p) 551 cdef PetscInt pStart = 0, pEnd = 0 552 CHKERR(DMPlexGetChart(self.dm, &pStart, &pEnd)) 553 assert cp>=pStart and cp<pEnd 554 cdef PetscInt csize = asInt(size) 555 CHKERR(DMPlexSetConeSize(self.dm, cp, csize)) 556 557 def getCone(self, p: int) -> ArrayInt: 558 """Return the points on the in-edges for this point in the DAG. 559 560 Not collective. 561 562 Parameters 563 ---------- 564 p 565 The point, which must lie in the chart set with `DMPlex.setChart`. 566 567 See Also 568 -------- 569 DM, DMPlex, DMPlex.getConeSize, DMPlex.setCone, DMPlex.setChart 570 petsc.DMPlexGetCone 571 572 """ 573 cdef PetscInt cp = asInt(p) 574 cdef PetscInt pStart = 0, pEnd = 0 575 CHKERR(DMPlexGetChart(self.dm, &pStart, &pEnd)) 576 assert cp>=pStart and cp<pEnd 577 cdef PetscInt ncone = 0 578 cdef const PetscInt *icone = NULL 579 CHKERR(DMPlexGetConeSize(self.dm, cp, &ncone)) 580 CHKERR(DMPlexGetCone(self.dm, cp, &icone)) 581 return array_i(ncone, icone) 582 583 def setCone(self, p: int, cone: Sequence[int], orientation: Sequence[int] | None = None) -> None: 584 """Set the points on the in-edges for this point in the DAG. 585 586 Not collective. 587 588 Parameters 589 ---------- 590 p 591 The point, which must lie in the chart set with `DMPlex.setChart`. 592 cone 593 An array of points which are on the in-edges for point ``p``. 594 orientation 595 An array of orientations, defaults to `None`. 596 597 See Also 598 -------- 599 DM, DMPlex, DMPlex.create, DMPlex.getCone, DMPlex.setChart 600 DMPlex.setConeSize, DM.setUp, DMPlex.setSupport 601 DMPlex.setSupportSize, petsc.DMPlexSetCone 602 603 """ 604 cdef PetscInt cp = asInt(p) 605 cdef PetscInt pStart = 0, pEnd = 0 606 CHKERR(DMPlexGetChart(self.dm, &pStart, &pEnd)) 607 assert cp>=pStart and cp<pEnd 608 # 609 cdef PetscInt ncone = 0 610 cdef PetscInt *icone = NULL 611 cone = iarray_i(cone, &ncone, &icone) 612 CHKERR(DMPlexSetConeSize(self.dm, cp, ncone)) 613 CHKERR(DMPlexSetCone(self.dm, cp, icone)) 614 # 615 cdef PetscInt norie = 0 616 cdef PetscInt *iorie = NULL 617 if orientation is not None: 618 orientation = iarray_i(orientation, &norie, &iorie) 619 assert norie == ncone 620 CHKERR(DMPlexSetConeOrientation(self.dm, cp, iorie)) 621 622 def insertCone(self, p: int, conePos: int, conePoint: int) -> None: 623 """DMPlexInsertCone - Insert a point into the in-edges for the point p in the DAG. 624 625 Not collective. 626 627 Parameters 628 ---------- 629 p 630 The point, which must lie in the chart set with `DMPlex.setChart`. 631 conePos 632 The local index in the cone where the point should be put. 633 conePoint 634 The mesh point to insert. 635 636 See Also 637 -------- 638 DM, DMPlex, DMPlex.create, DMPlex.getCone, DMPlex.setChart 639 DMPlex.setConeSize, DM.setUp, petsc.DMPlexInsertCone 640 641 """ 642 cdef PetscInt cp = asInt(p) 643 cdef PetscInt cconePos = asInt(conePos) 644 cdef PetscInt cconePoint = asInt(conePoint) 645 CHKERR(DMPlexInsertCone(self.dm, cp, cconePos, cconePoint)) 646 647 def insertConeOrientation(self, p: int, conePos: int, coneOrientation: int) -> None: 648 """Insert a point orientation for the in-edge for the point p in the DAG. 649 650 Not collective. 651 652 Parameters 653 ---------- 654 p 655 The point, which must lie in the chart set with `DMPlex.setChart` 656 conePos 657 The local index in the cone where the point should be put. 658 coneOrientation 659 The point orientation to insert. 660 661 See Also 662 -------- 663 DM, DMPlex, DMPlex.create, DMPlex.getCone, DMPlex.setChart 664 DMPlex.setConeSize, DM.setUp, petsc.DMPlexInsertConeOrientation 665 666 """ 667 cdef PetscInt cp = asInt(p) 668 cdef PetscInt cconePos = asInt(conePos) 669 cdef PetscInt cconeOrientation = asInt(coneOrientation) 670 CHKERR(DMPlexInsertConeOrientation(self.dm, cp, cconePos, cconeOrientation)) 671 672 def getConeOrientation(self, p: int) -> ArrayInt: 673 """Return the orientations on the in-edges for this point in the DAG. 674 675 Not collective. 676 677 Parameters 678 ---------- 679 p 680 The point, which must lie in the chart set with `DMPlex.setChart`. 681 682 See Also 683 -------- 684 DM, DMPlex, DMPlex.create, DMPlex.getCone, DMPlex.setCone 685 DMPlex.setChart, petsc.DMPlexGetConeOrientation 686 687 """ 688 cdef PetscInt cp = asInt(p) 689 cdef PetscInt pStart = 0, pEnd = 0 690 CHKERR(DMPlexGetChart(self.dm, &pStart, &pEnd)) 691 assert cp>=pStart and cp<pEnd 692 cdef PetscInt norie = 0 693 cdef const PetscInt *iorie = NULL 694 CHKERR(DMPlexGetConeSize(self.dm, cp, &norie)) 695 CHKERR(DMPlexGetConeOrientation(self.dm, cp, &iorie)) 696 return array_i(norie, iorie) 697 698 def setConeOrientation(self, p: int, orientation: Sequence[int]) -> None: 699 """Set the orientations on the in-edges for this point in the DAG. 700 701 Not collective. 702 703 Parameters 704 ---------- 705 p 706 The point, which must lie in the chart set with `DMPlex.setChart`. 707 orientation 708 An array of orientations. 709 710 See Also 711 -------- 712 DM, DMPlex, DMPlex.create, DMPlex.getConeOrientation, DMPlex.setCone 713 DMPlex.setChart, DMPlex.setConeSize, DM.setUp 714 petsc.DMPlexSetConeOrientation 715 716 """ 717 cdef PetscInt cp = asInt(p) 718 cdef PetscInt pStart = 0, pEnd = 0 719 CHKERR(DMPlexGetChart(self.dm, &pStart, &pEnd)) 720 assert cp>=pStart and cp<pEnd 721 cdef PetscInt ncone = 0 722 CHKERR(DMPlexGetConeSize(self.dm, cp, &ncone)) 723 cdef PetscInt norie = 0 724 cdef PetscInt *iorie = NULL 725 orientation = iarray_i(orientation, &norie, &iorie) 726 assert norie == ncone 727 CHKERR(DMPlexSetConeOrientation(self.dm, cp, iorie)) 728 729 def setCellType(self, p: int, ctype: DM.PolytopeType) -> None: 730 """Set the polytope type of a given cell. 731 732 Not collective. 733 734 Parameters 735 ---------- 736 p 737 The cell. 738 ctype 739 The polytope type of the cell. 740 741 See Also 742 -------- 743 DM, DMPlex, DMPlex.getCellTypeLabel, DMPlex.getDepth, DM.createLabel 744 petsc.DMPlexSetCellType 745 746 """ 747 cdef PetscInt cp = asInt(p) 748 cdef PetscDMPolytopeType val = ctype 749 CHKERR(DMPlexSetCellType(self.dm, cp, val)) 750 751 def getCellType(self, p: int) -> DM.PolytopeType: 752 """Return the polytope type of a given cell. 753 754 Not collective. 755 756 Parameters 757 ---------- 758 p 759 The cell. 760 761 See Also 762 -------- 763 DM, DMPlex, DM.PolytopeType, DMPlex.getCellTypeLabel, DMPlex.getDepth 764 petsc.DMPlexGetCellType 765 766 """ 767 cdef PetscInt cp = asInt(p) 768 cdef PetscDMPolytopeType ctype = DM_POLYTOPE_UNKNOWN 769 CHKERR(DMPlexGetCellType(self.dm, cp, &ctype)) 770 return toInt(ctype) 771 772 def getCellTypeLabel(self) -> DMLabel: 773 """Return the `DMLabel` recording the polytope type of each cell. 774 775 Not collective. 776 777 See Also 778 -------- 779 DM, DMPlex, DMPlex.getCellType, DM.createLabel 780 petsc.DMPlexGetCellTypeLabel 781 782 """ 783 cdef DMLabel label = DMLabel() 784 CHKERR(DMPlexGetCellTypeLabel(self.dm, &label.dmlabel)) 785 CHKERR(PetscINCREF(label.obj)) 786 return label 787 788 def getSupportSize(self, p: int) -> int: 789 """Return the number of out-edges for this point in the DAG. 790 791 Not collective. 792 793 Parameters 794 ---------- 795 p 796 The point, which must lie in the chart set with `DMPlex.setChart`. 797 798 See Also 799 -------- 800 DM, DMPlex, DMPlex.create, DMPlex.setConeSize, DMPlex.setChart 801 DMPlex.getConeSize, petsc.DMPlexGetSupportSize 802 803 """ 804 cdef PetscInt cp = asInt(p) 805 cdef PetscInt pStart = 0, pEnd = 0 806 CHKERR(DMPlexGetChart(self.dm, &pStart, &pEnd)) 807 assert cp>=pStart and cp<pEnd 808 cdef PetscInt ssize = 0 809 CHKERR(DMPlexGetSupportSize(self.dm, cp, &ssize)) 810 return toInt(ssize) 811 812 def setSupportSize(self, p: int, size: int) -> None: 813 """Set the number of out-edges for this point in the DAG. 814 815 Not collective. 816 817 Parameters 818 ---------- 819 p 820 The point, which must lie in the chart set with `DMPlex.setChart`. 821 size 822 The support size for point ``p``. 823 824 See Also 825 -------- 826 DM, DMPlex, DMPlex.create, DMPlex.getSupportSize, DMPlex.setChart 827 petsc.DMPlexSetSupportSize 828 829 """ 830 cdef PetscInt cp = asInt(p) 831 cdef PetscInt pStart = 0, pEnd = 0 832 CHKERR(DMPlexGetChart(self.dm, &pStart, &pEnd)) 833 assert cp>=pStart and cp<pEnd 834 cdef PetscInt ssize = asInt(size) 835 CHKERR(DMPlexSetSupportSize(self.dm, cp, ssize)) 836 837 def getSupport(self, p: int) -> ArrayInt: 838 """Return the points on the out-edges for this point in the DAG. 839 840 Not collective. 841 842 Parameters 843 ---------- 844 p 845 The point, which must lie in the chart set with `DMPlex.setChart`. 846 847 See Also 848 -------- 849 DM, DMPlex, DMPlex.getSupportSize, DMPlex.setSupport, DMPlex.getCone 850 DMPlex.setChart, petsc.DMPlexGetSupport 851 852 """ 853 cdef PetscInt cp = asInt(p) 854 cdef PetscInt pStart = 0, pEnd = 0 855 CHKERR(DMPlexGetChart(self.dm, &pStart, &pEnd)) 856 assert cp>=pStart and cp<pEnd 857 cdef PetscInt nsupp = 0 858 cdef const PetscInt *isupp = NULL 859 CHKERR(DMPlexGetSupportSize(self.dm, cp, &nsupp)) 860 CHKERR(DMPlexGetSupport(self.dm, cp, &isupp)) 861 return array_i(nsupp, isupp) 862 863 def setSupport(self, p: int, supp: Sequence[int]) -> None: 864 """Set the points on the out-edges for this point in the DAG. 865 866 Not collective. 867 868 Parameters 869 ---------- 870 p 871 The point, which must lie in the chart set with `DMPlex.setChart`. 872 supp 873 An array of points which are on the out-edges for point ``p``. 874 875 See Also 876 -------- 877 DM, DMPlex, DMPlex.setCone, DMPlex.setConeSize, DMPlex.create 878 DMPlex.getSupport, DMPlex.setChart, DMPlex.setSupportSize, DM.setUp 879 petsc.DMPlexSetSupport 880 881 """ 882 cdef PetscInt cp = asInt(p) 883 cdef PetscInt pStart = 0, pEnd = 0 884 CHKERR(DMPlexGetChart(self.dm, &pStart, &pEnd)) 885 assert cp>=pStart and cp<pEnd 886 cdef PetscInt nsupp = 0 887 cdef PetscInt *isupp = NULL 888 supp = iarray_i(supp, &nsupp, &isupp) 889 CHKERR(DMPlexSetSupportSize(self.dm, cp, nsupp)) 890 CHKERR(DMPlexSetSupport(self.dm, cp, isupp)) 891 892 def getMaxSizes(self) -> tuple[int, int]: 893 """Return the maximum number of in-edges and out-edges of the DAG. 894 895 Not collective. 896 897 Returns 898 ------- 899 maxConeSize : int 900 The maximum number of in-edges. 901 maxSupportSize : int 902 The maximum number of out-edges. 903 904 See Also 905 -------- 906 DM, DMPlex, DMPlex.create, DMPlex.setConeSize, DMPlex.setChart 907 petsc.DMPlexGetMaxSizes 908 909 """ 910 cdef PetscInt maxConeSize = 0, maxSupportSize = 0 911 CHKERR(DMPlexGetMaxSizes(self.dm, &maxConeSize, &maxSupportSize)) 912 return toInt(maxConeSize), toInt(maxSupportSize) 913 914 def symmetrize(self) -> None: 915 """Create support (out-edge) information from cone (in-edge) information. 916 917 Not collective. 918 919 See Also 920 -------- 921 DM, DMPlex, DMPlex.create, DMPlex.setChart, DMPlex.setConeSize 922 DMPlex.setCone, petsc.DMPlexSymmetrize 923 924 """ 925 CHKERR(DMPlexSymmetrize(self.dm)) 926 927 def stratify(self) -> None: 928 """Calculate the strata of DAG. 929 930 Collective. 931 932 See Also 933 -------- 934 DM, DMPlex, DMPlex.create, DMPlex.symmetrize, petsc.DMPlexStratify 935 936 """ 937 CHKERR(DMPlexStratify(self.dm)) 938 939 def orient(self) -> None: 940 """Give a consistent orientation to the input mesh. 941 942 Collective. 943 944 See Also 945 -------- 946 DM, DMPlex, DM.create, petsc.DMPlexOrient 947 948 """ 949 CHKERR(DMPlexOrient(self.dm)) 950 951 def getCellNumbering(self) -> IS: 952 """Return a global cell numbering for all cells on this process. 953 954 Collective the first time it is called. 955 956 See Also 957 -------- 958 DM, DMPlex, DMPlex.getVertexNumbering, petsc.DMPlexGetCellNumbering 959 960 """ 961 cdef IS iset = IS() 962 CHKERR(DMPlexGetCellNumbering(self.dm, &iset.iset)) 963 CHKERR(PetscINCREF(iset.obj)) 964 return iset 965 966 def getVertexNumbering(self) -> IS: 967 """Return a global vertex numbering for all vertices on this process. 968 969 Collective the first time it is called. 970 971 See Also 972 -------- 973 DM, DMPlex, DMPlex.getCellNumbering, petsc.DMPlexGetVertexNumbering 974 975 """ 976 cdef IS iset = IS() 977 CHKERR(DMPlexGetVertexNumbering(self.dm, &iset.iset)) 978 CHKERR(PetscINCREF(iset.obj)) 979 return iset 980 981 def createPointNumbering(self) -> IS: 982 """Create a global numbering for all points. 983 984 Collective. 985 986 See Also 987 -------- 988 DM, DMPlex, DMPlex.getCellNumbering, petsc.DMPlexCreatePointNumbering 989 990 """ 991 cdef IS iset = IS() 992 CHKERR(DMPlexCreatePointNumbering(self.dm, &iset.iset)) 993 return iset 994 995 def getDepth(self) -> int: 996 """Return the depth of the DAG representing this mesh. 997 998 Not collective. 999 1000 See Also 1001 -------- 1002 DM, DMPlex, DMPlex.getDepthStratum, DMPlex.symmetrize 1003 petsc.DMPlexGetDepth 1004 1005 """ 1006 cdef PetscInt depth = 0 1007 CHKERR(DMPlexGetDepth(self.dm, &depth)) 1008 return toInt(depth) 1009 1010 def getDepthStratum(self, svalue: int) -> tuple[int, int]: 1011 """Return the bounds [``start``, ``end``) for all points at a certain depth. 1012 1013 Not collective. 1014 1015 Parameters 1016 ---------- 1017 svalue 1018 The requested depth. 1019 1020 Returns 1021 ------- 1022 pStart : int 1023 The first stratum point. 1024 pEnd : int 1025 The upper bound for stratum points. 1026 1027 See Also 1028 -------- 1029 DM, DMPlex, DMPlex.getHeightStratum, DMPlex.getDepth 1030 DMPlex.symmetrize, DMPlex.interpolate, petsc.DMPlexGetDepthStratum 1031 1032 """ 1033 cdef PetscInt csvalue = asInt(svalue), sStart = 0, sEnd = 0 1034 CHKERR(DMPlexGetDepthStratum(self.dm, csvalue, &sStart, &sEnd)) 1035 return (toInt(sStart), toInt(sEnd)) 1036 1037 def getHeightStratum(self, svalue: int) -> tuple[int, int]: 1038 """Return the bounds [``start``, ``end``) for all points at a certain height. 1039 1040 Not collective. 1041 1042 Parameters 1043 ---------- 1044 svalue 1045 The requested height. 1046 1047 Returns 1048 ------- 1049 pStart : int 1050 The first stratum point. 1051 pEnd : int 1052 The upper bound for stratum points. 1053 1054 See Also 1055 -------- 1056 DM, DMPlex, DMPlex.getDepthStratum, DMPlex.getDepth 1057 petsc.DMPlexGetHeightStratum 1058 1059 """ 1060 cdef PetscInt csvalue = asInt(svalue), sStart = 0, sEnd = 0 1061 CHKERR(DMPlexGetHeightStratum(self.dm, csvalue, &sStart, &sEnd)) 1062 return (toInt(sStart), toInt(sEnd)) 1063 1064 def getPointDepth(self, point: int) -> int: 1065 """Return the *depth* of a given point. 1066 1067 Not collective. 1068 1069 Parameters 1070 ---------- 1071 point 1072 The point. 1073 1074 See Also 1075 -------- 1076 DM, DMPlex, DMPlex.getDepthStratum, DMPlex.getDepth 1077 petsc.DMPlexGetPointDepth 1078 1079 """ 1080 cdef PetscInt cpoint = asInt(point) 1081 cdef PetscInt depth = 0 1082 CHKERR(DMPlexGetPointDepth(self.dm, cpoint, &depth)) 1083 return toInt(depth) 1084 1085 def getPointHeight(self, point: int) -> int: 1086 """Return the *height* of a given point. 1087 1088 Not collective. 1089 1090 Parameters 1091 ---------- 1092 point 1093 The point. 1094 1095 See Also 1096 -------- 1097 DM, DMPlex, DMPlex.getHeightStratum 1098 petsc.DMPlexGetPointHeight 1099 1100 """ 1101 cdef PetscInt cpoint = asInt(point) 1102 cdef PetscInt height = 0 1103 CHKERR(DMPlexGetPointHeight(self.dm, cpoint, &height)) 1104 return toInt(height) 1105 1106 def getMeet(self, points: Sequence[int]) -> ArrayInt: 1107 """Return an array for the meet of the set of points. 1108 1109 Not collective. 1110 1111 Parameters 1112 ---------- 1113 points 1114 The input points. 1115 1116 See Also 1117 -------- 1118 DM, DMPlex, DMPlex.getJoin, petsc.DMPlexGetMeet 1119 1120 """ 1121 cdef PetscInt numPoints = 0 1122 cdef PetscInt *ipoints = NULL 1123 cdef PetscInt numCoveringPoints = 0 1124 cdef const PetscInt *coveringPoints = NULL 1125 points = iarray_i(points, &numPoints, &ipoints) 1126 CHKERR(DMPlexGetMeet(self.dm, numPoints, ipoints, &numCoveringPoints, &coveringPoints)) 1127 try: 1128 return array_i(numCoveringPoints, coveringPoints) 1129 finally: 1130 CHKERR(DMPlexRestoreMeet(self.dm, numPoints, ipoints, &numCoveringPoints, &coveringPoints)) 1131 1132 def getJoin(self, points: Sequence[int]) -> ArrayInt: 1133 """Return an array for the join of the set of points. 1134 1135 Not collective. 1136 1137 Parameters 1138 ---------- 1139 points 1140 The input points. 1141 1142 See Also 1143 -------- 1144 DM, DMPlex, DMPlex.getMeet, petsc.DMPlexGetJoin 1145 1146 """ 1147 cdef PetscInt numPoints = 0 1148 cdef PetscInt *ipoints = NULL 1149 cdef PetscInt numCoveringPoints = 0 1150 cdef const PetscInt *coveringPoints = NULL 1151 points = iarray_i(points, &numPoints, &ipoints) 1152 CHKERR(DMPlexGetJoin(self.dm, numPoints, ipoints, &numCoveringPoints, &coveringPoints)) 1153 try: 1154 return array_i(numCoveringPoints, coveringPoints) 1155 finally: 1156 CHKERR(DMPlexRestoreJoin(self.dm, numPoints, ipoints, &numCoveringPoints, &coveringPoints)) 1157 1158 def getFullJoin(self, points: Sequence[int]) -> ArrayInt: 1159 """Return an array for the join of the set of points. 1160 1161 Not collective. 1162 1163 Parameters 1164 ---------- 1165 points 1166 The input points. 1167 1168 See Also 1169 -------- 1170 DM, DMPlex, DMPlex.getJoin, DMPlex.getMeet, petsc.DMPlexGetFullJoin 1171 1172 """ 1173 cdef PetscInt numPoints = 0 1174 cdef PetscInt *ipoints = NULL 1175 cdef PetscInt numCoveringPoints = 0 1176 cdef const PetscInt *coveringPoints = NULL 1177 points = iarray_i(points, &numPoints, &ipoints) 1178 CHKERR(DMPlexGetFullJoin(self.dm, numPoints, ipoints, &numCoveringPoints, &coveringPoints)) 1179 try: 1180 return array_i(numCoveringPoints, coveringPoints) 1181 finally: 1182 CHKERR(DMPlexRestoreJoin(self.dm, numPoints, ipoints, &numCoveringPoints, &coveringPoints)) 1183 1184 def getTransitiveClosure(self, p: int, useCone: bool | None = True) -> tuple[ArrayInt, ArrayInt]: 1185 """Return the points and orientations on the transitive closure of this point. 1186 1187 Not collective. 1188 1189 Parameters 1190 ---------- 1191 p 1192 The mesh point. 1193 useCone 1194 `True` for the closure, otherwise return the star. 1195 1196 Returns 1197 ------- 1198 points : ArrayInt 1199 The points. 1200 orientations : ArrayInt 1201 The orientations. 1202 1203 See Also 1204 -------- 1205 DM, DMPlex, DMPlex.create, DMPlex.setCone, DMPlex.setChart 1206 DMPlex.getCone, petsc.DMPlexGetTransitiveClosure 1207 1208 """ 1209 cdef PetscInt cp = asInt(p) 1210 cdef PetscInt pStart = 0, pEnd = 0 1211 CHKERR(DMPlexGetChart(self.dm, &pStart, &pEnd)) 1212 assert cp>=pStart and cp<pEnd 1213 cdef PetscBool cuseCone = useCone 1214 cdef PetscInt numPoints = 0 1215 cdef PetscInt *points = NULL 1216 CHKERR(DMPlexGetTransitiveClosure(self.dm, cp, cuseCone, &numPoints, &points)) 1217 try: 1218 out = array_i(2*numPoints, points) 1219 finally: 1220 CHKERR(DMPlexRestoreTransitiveClosure(self.dm, cp, cuseCone, &numPoints, &points)) 1221 return out[::2], out[1::2] 1222 1223 def vecGetClosure(self, Section sec, Vec vec, p: int) -> ArrayScalar: 1224 """Return an array of values on the closure of ``p``. 1225 1226 Not collective. 1227 1228 Parameters 1229 ---------- 1230 sec 1231 The section describing the layout in ``vec``. 1232 vec 1233 The local vector. 1234 p 1235 The point in the `DMPlex`. 1236 1237 See Also 1238 -------- 1239 DM, DMPlex, petsc.DMPlexVecRestoreClosure 1240 1241 """ 1242 cdef PetscInt cp = asInt(p), csize = 0 1243 cdef PetscScalar *cvals = NULL 1244 CHKERR(DMPlexVecGetClosure(self.dm, sec.sec, vec.vec, cp, &csize, &cvals)) 1245 try: 1246 closure = array_s(csize, cvals) 1247 finally: 1248 CHKERR(DMPlexVecRestoreClosure(self.dm, sec.sec, vec.vec, cp, &csize, &cvals)) 1249 return closure 1250 1251 def getVecClosure(self, Section sec or None, Vec vec, point: int) -> ArrayScalar: 1252 """Return an array of the values on the closure of a point. 1253 1254 Not collective. 1255 1256 Parameters 1257 ---------- 1258 sec 1259 The `Section` describing the layout in ``vec`` 1260 or `None` to use the default section. 1261 vec 1262 The local vector. 1263 point 1264 The point in the `DMPlex`. 1265 1266 See Also 1267 -------- 1268 DM, DMPlex, petsc.DMPlexVecRestoreClosure 1269 1270 """ 1271 cdef PetscSection csec = sec.sec if sec is not None else NULL 1272 cdef PetscInt cp = asInt(point), csize = 0 1273 cdef PetscScalar *cvals = NULL 1274 CHKERR(DMPlexVecGetClosure(self.dm, csec, vec.vec, cp, &csize, &cvals)) 1275 try: 1276 closure = array_s(csize, cvals) 1277 finally: 1278 CHKERR(DMPlexVecRestoreClosure(self.dm, csec, vec.vec, cp, &csize, &cvals)) 1279 return closure 1280 1281 def setVecClosure(self, Section sec or None, Vec vec, point: int, values: Sequence[Scalar], addv: InsertModeSpec | None = None) -> None: 1282 """Set an array of the values on the closure of ``point``. 1283 1284 Not collective. 1285 1286 Parameters 1287 ---------- 1288 sec 1289 The section describing the layout in ``vec``, 1290 or `None` to use the default section. 1291 vec 1292 The local vector. 1293 point 1294 The point in the `DMPlex`. 1295 values 1296 The array of values. 1297 mode 1298 The insertion mode. 1299 1300 See Also 1301 -------- 1302 DM, DMPlex, petsc.DMPlexVecSetClosure 1303 1304 """ 1305 cdef PetscSection csec = sec.sec if sec is not None else NULL 1306 cdef PetscInt cp = asInt(point) 1307 cdef PetscInt csize = 0 1308 cdef PetscScalar *cvals = NULL 1309 cdef object unused = iarray_s(values, &csize, &cvals) 1310 cdef PetscInsertMode im = insertmode(addv) 1311 CHKERR(DMPlexVecSetClosure(self.dm, csec, vec.vec, cp, cvals, im)) 1312 1313 def setMatClosure(self, Section sec or None, Section gsec or None, 1314 Mat mat, point: int, values: Sequence[Scalar], addv: InsertModeSpec | None = None) -> None: 1315 """Set an array of the values on the closure of ``point``. 1316 1317 Not collective. 1318 1319 Parameters 1320 ---------- 1321 sec 1322 The section describing the layout in ``mat``, 1323 or `None` to use the default section. 1324 gsec 1325 The section describing the layout in ``mat``, 1326 or `None` to use the default global section. 1327 mat 1328 The matrix. 1329 point 1330 The point in the `DMPlex`. 1331 values 1332 The array of values. 1333 mode 1334 The insertion mode. 1335 1336 See Also 1337 -------- 1338 DM, DMPlex, petsc.DMPlexMatSetClosure 1339 1340 """ 1341 cdef PetscSection csec = sec.sec if sec is not None else NULL 1342 cdef PetscSection cgsec = gsec.sec if gsec is not None else NULL 1343 cdef PetscInt cp = asInt(point) 1344 cdef PetscInt csize = 0 1345 cdef PetscScalar *cvals = NULL 1346 cdef object unused = iarray_s(values, &csize, &cvals) 1347 cdef PetscInsertMode im = insertmode(addv) 1348 CHKERR(DMPlexMatSetClosure(self.dm, csec, cgsec, mat.mat, cp, cvals, im)) 1349 1350 def generate(self, DMPlex boundary, name: str | None = None, interpolate: bool | None = True) -> Self: 1351 """Generate a mesh. 1352 1353 Not collective. 1354 1355 Parameters 1356 ---------- 1357 boundary 1358 The `DMPlex` boundary object. 1359 name 1360 The mesh generation package name. 1361 interpolate 1362 Flag to create intermediate mesh elements. 1363 1364 See Also 1365 -------- 1366 DM, DMPlex, DMPlex.create, DM.refine, petsc_options 1367 petsc.DMPlexGenerate 1368 1369 """ 1370 cdef PetscBool interp = interpolate 1371 cdef const char *cname = NULL 1372 if name: name = str2bytes(name, &cname) 1373 cdef PetscDM newdm = NULL 1374 CHKERR(DMPlexGenerate(boundary.dm, cname, interp, &newdm)) 1375 CHKERR(PetscCLEAR(self.obj)); self.dm = newdm 1376 return self 1377 1378 def setTriangleOptions(self, opts: str) -> None: 1379 """Set the options used for the Triangle mesh generator. 1380 1381 Not collective. 1382 1383 Parameters 1384 ---------- 1385 opts 1386 The command line options. 1387 1388 See Also 1389 -------- 1390 petsc_options, DM, DMPlex, DMPlex.setTetGenOptions, DMPlex.generate 1391 petsc.DMPlexTriangleSetOptions 1392 1393 """ 1394 cdef const char *copts = NULL 1395 opts = str2bytes(opts, &copts) 1396 CHKERR(DMPlexTriangleSetOptions(self.dm, copts)) 1397 1398 def setTetGenOptions(self, opts: str) -> None: 1399 """Set the options used for the Tetgen mesh generator. 1400 1401 Not collective. 1402 1403 Parameters 1404 ---------- 1405 opts 1406 The command line options. 1407 1408 See Also 1409 -------- 1410 petsc_options, DM, DMPlex, DMPlex.setTriangleOptions, DMPlex.generate 1411 petsc.DMPlexTetgenSetOptions 1412 1413 """ 1414 cdef const char *copts = NULL 1415 opts = str2bytes(opts, &copts) 1416 CHKERR(DMPlexTetgenSetOptions(self.dm, copts)) 1417 1418 def markBoundaryFaces(self, label: str, value: int | None = None) -> DMLabel: 1419 """Mark all faces on the boundary. 1420 1421 Not collective. 1422 1423 Parameters 1424 ---------- 1425 value 1426 The marker value, or `DETERMINE` or `None` to use some 1427 value in the closure (or 1 if none are found). 1428 1429 See Also 1430 -------- 1431 DM, DMPlex, DMLabel.create, DM.createLabel 1432 petsc.DMPlexMarkBoundaryFaces 1433 1434 """ 1435 cdef PetscInt ival = PETSC_DETERMINE 1436 if value is not None: ival = asInt(value) 1437 if not self.hasLabel(label): 1438 self.createLabel(label) 1439 cdef const char *cval = NULL 1440 label = str2bytes(label, &cval) 1441 cdef PetscDMLabel clbl = NULL 1442 CHKERR(DMGetLabel(self.dm, cval, &clbl)) 1443 CHKERR(DMPlexMarkBoundaryFaces(self.dm, ival, clbl)) 1444 1445 def labelComplete(self, DMLabel label) -> None: 1446 """Add the transitive closure to the surface. 1447 1448 Not collective. 1449 1450 Parameters 1451 ---------- 1452 label 1453 A `DMLabel` marking the surface points. 1454 1455 See Also 1456 -------- 1457 DM, DMPlex, DMPlex.labelCohesiveComplete, petsc.DMPlexLabelComplete 1458 1459 """ 1460 CHKERR(DMPlexLabelComplete(self.dm, label.dmlabel)) 1461 1462 def labelCohesiveComplete(self, DMLabel label, DMLabel bdlabel, bdvalue: int, 1463 flip: bool, split: bool, DMPlex subdm) -> None: 1464 """Add all other mesh pieces to complete the surface. 1465 1466 Not collective. 1467 1468 Parameters 1469 ---------- 1470 label 1471 A `DMLabel` marking the surface. 1472 bdlabel 1473 A `DMLabel` marking the vertices on the boundary 1474 which will not be duplicated. 1475 bdvalue 1476 Value of `DMLabel` marking the vertices on the boundary. 1477 flip 1478 Flag to flip the submesh normal and replace points 1479 on the other side. 1480 split 1481 Flag to split faces incident on the surface boundary, 1482 rather than clamping those faces to the boundary 1483 subdm 1484 The `DMPlex` associated with the label. 1485 1486 See Also 1487 -------- 1488 DM, DMPlex, DMPlex.labelComplete 1489 petsc.DMPlexLabelCohesiveComplete 1490 1491 """ 1492 cdef PetscBool flg = flip 1493 cdef PetscBool flg2 = split 1494 cdef PetscInt val = asInt(bdvalue) 1495 CHKERR(DMPlexLabelCohesiveComplete(self.dm, label.dmlabel, bdlabel.dmlabel, val, flg, flg2, subdm.dm)) 1496 1497 def setAdjacencyUseAnchors(self, useAnchors: bool = True) -> None: 1498 """Define adjacency in the mesh using the point-to-point constraints. 1499 1500 Logically collective. 1501 1502 Parameters 1503 ---------- 1504 useAnchors 1505 Flag to use the constraints. 1506 If `True`, then constrained points are omitted from `DMPlex.getAdjacency`, 1507 and their anchor points appear in their place. 1508 1509 See Also 1510 -------- 1511 DMPlex, DMPlex.getAdjacency, DMPlex.distribute 1512 petsc.DMPlexSetAdjacencyUseAnchors 1513 1514 """ 1515 cdef PetscBool flag = useAnchors 1516 CHKERR(DMPlexSetAdjacencyUseAnchors(self.dm, flag)) 1517 1518 def getAdjacencyUseAnchors(self) -> bool: 1519 """Query whether adjacency in the mesh uses the point-to-point constraints. 1520 1521 Not collective. 1522 1523 See Also 1524 -------- 1525 DMPlex, DMPlex.getAdjacency, DMPlex.distribute 1526 petsc.DMPlexGetAdjacencyUseAnchors 1527 1528 """ 1529 cdef PetscBool flag = PETSC_FALSE 1530 CHKERR(DMPlexGetAdjacencyUseAnchors(self.dm, &flag)) 1531 return toBool(flag) 1532 1533 def getAdjacency(self, p: int) -> ArrayInt: 1534 """Return all points adjacent to the given point. 1535 1536 Not collective. 1537 1538 Parameters 1539 ---------- 1540 p 1541 The point. 1542 1543 See Also 1544 -------- 1545 DMPlex, DMPlex.distribute, petsc.DMPlexGetAdjacency 1546 1547 """ 1548 cdef PetscInt cp = asInt(p) 1549 cdef PetscInt nadj = PETSC_DETERMINE 1550 cdef PetscInt *iadj = NULL 1551 CHKERR(DMPlexGetAdjacency(self.dm, cp, &nadj, &iadj)) 1552 try: 1553 adjacency = array_i(nadj, iadj) 1554 finally: 1555 CHKERR(PetscFree(iadj)) 1556 return adjacency 1557 1558 def setPartitioner(self, Partitioner part) -> None: 1559 """Set the mesh partitioner. 1560 1561 Logically collective. 1562 1563 Parameters 1564 ---------- 1565 part 1566 The partitioner. 1567 1568 See Also 1569 -------- 1570 DM, DMPlex, Partitioner, DMPlex.distribute, DMPlex.getPartitioner 1571 Partitioner.create, petsc.DMPlexSetPartitioner 1572 1573 """ 1574 CHKERR(DMPlexSetPartitioner(self.dm, part.part)) 1575 1576 def getPartitioner(self) -> Partitioner: 1577 """Return the mesh partitioner. 1578 1579 Not collective. 1580 1581 See Also 1582 -------- 1583 DM, DMPlex, Partitioner, Section, DMPlex.distribute 1584 DMPlex.setPartitioner, Partitioner.create 1585 petsc.PetscPartitionerDMPlexPartition, petsc.DMPlexGetPartitioner 1586 1587 """ 1588 cdef Partitioner part = Partitioner() 1589 CHKERR(DMPlexGetPartitioner(self.dm, &part.part)) 1590 CHKERR(PetscINCREF(part.obj)) 1591 return part 1592 1593 def rebalanceSharedPoints(self, entityDepth: int | None = 0, useInitialGuess: bool | None = True, parallel: bool | None = True) -> bool: 1594 """Redistribute shared points in order to achieve better balancing. 1595 1596 Collective. 1597 1598 Parameters 1599 ---------- 1600 entityDepth 1601 Depth of the entity to balance (e.g., 0 -> balance vertices). 1602 useInitialGuess 1603 Whether to use the current distribution as initial guess. 1604 parallel 1605 Whether to use ParMETIS and do the partition in parallel 1606 or gather the graph onto a single process. 1607 1608 Returns 1609 ------- 1610 success : bool 1611 Whether the graph partitioning was successful or not. 1612 Unsuccessful simply means no change to the partitioning. 1613 1614 See Also 1615 -------- 1616 DM, DMPlex, DMPlex.distribute, petsc_options 1617 petsc.DMPlexRebalanceSharedPoints 1618 1619 """ 1620 cdef PetscInt centityDepth = asInt(entityDepth) 1621 cdef PetscBool cuseInitialGuess = asBool(useInitialGuess) 1622 cdef PetscBool cparallel = asBool(parallel) 1623 cdef PetscBool csuccess = PETSC_FALSE 1624 CHKERR(DMPlexRebalanceSharedPoints(self.dm, centityDepth, cuseInitialGuess, cparallel, &csuccess)) 1625 return toBool(csuccess) 1626 1627 def distribute(self, overlap: int | None = 0) -> SF | None: 1628 """Distribute the mesh and any associated sections. 1629 1630 Collective. 1631 1632 Parameters 1633 ---------- 1634 overlap 1635 The overlap of partitions. 1636 1637 Returns 1638 ------- 1639 sf : SF or None 1640 The `SF` used for point distribution, or `None` if not distributed. 1641 1642 See Also 1643 -------- 1644 DM, DMPlex, DMPlex.create, petsc.DMPlexDistribute 1645 1646 """ 1647 cdef PetscDM dmParallel = NULL 1648 cdef PetscInt coverlap = asInt(overlap) 1649 cdef SF sf = SF() 1650 CHKERR(DMPlexDistribute(self.dm, coverlap, &sf.sf, &dmParallel)) 1651 if dmParallel != NULL: 1652 CHKERR(PetscCLEAR(self.obj)); self.dm = dmParallel 1653 return sf 1654 1655 def distributeOverlap(self, overlap: int | None = 0) -> SF: 1656 """Add partition overlap to a distributed non-overlapping `DMPlex`. 1657 1658 Collective. 1659 1660 Parameters 1661 ---------- 1662 overlap 1663 The overlap of partitions (the same on all ranks). 1664 1665 Returns 1666 ------- 1667 sf : SF 1668 The `SF` used for point distribution. 1669 1670 See Also 1671 -------- 1672 DM, DMPlex, SF, DMPlex.create, DMPlex.distribute, 1673 petsc.DMPlexDistributeOverlap 1674 1675 """ 1676 cdef PetscInt coverlap = asInt(overlap) 1677 cdef SF sf = SF() 1678 cdef PetscDM dmOverlap = NULL 1679 CHKERR(DMPlexDistributeOverlap(self.dm, coverlap, 1680 &sf.sf, &dmOverlap)) 1681 CHKERR(PetscCLEAR(self.obj)); self.dm = dmOverlap 1682 return sf 1683 1684 def isDistributed(self) -> bool: 1685 """Return the flag indicating if the mesh is distributed. 1686 1687 Collective. 1688 1689 See Also 1690 -------- 1691 DM, DMPlex, DMPlex.distribute, petsc.DMPlexIsDistributed 1692 1693 """ 1694 cdef PetscBool flag = PETSC_FALSE 1695 CHKERR(DMPlexIsDistributed(self.dm, &flag)) 1696 return toBool(flag) 1697 1698 def isSimplex(self) -> bool: 1699 """Return the flag indicating if the first cell is a simplex. 1700 1701 Not collective. 1702 1703 See Also 1704 -------- 1705 DM, DMPlex, DMPlex.getCellType, DMPlex.getHeightStratum 1706 petsc.DMPlexIsSimplex 1707 1708 """ 1709 cdef PetscBool flag = PETSC_FALSE 1710 CHKERR(DMPlexIsSimplex(self.dm, &flag)) 1711 return toBool(flag) 1712 1713 def distributeGetDefault(self) -> bool: 1714 """Return a flag indicating whether the `DM` should be distributed by default. 1715 1716 Not collective. 1717 1718 Returns 1719 ------- 1720 dist : bool 1721 Flag indicating whether the `DMPlex` should be distributed by default. 1722 1723 See Also 1724 -------- 1725 DM, DMPlex, DMPlex.distributeSetDefault, DMPlex.distribute 1726 petsc.DMPlexDistributeGetDefault 1727 1728 """ 1729 cdef PetscBool dist = PETSC_FALSE 1730 CHKERR(DMPlexDistributeGetDefault(self.dm, &dist)) 1731 return toBool(dist) 1732 1733 def distributeSetDefault(self, flag: bool) -> None: 1734 """Set flag indicating whether the `DMPlex` should be distributed by default. 1735 1736 Logically collective. 1737 1738 Parameters 1739 ---------- 1740 flag 1741 Flag indicating whether the `DMPlex` should be distributed by default. 1742 1743 See Also 1744 -------- 1745 DMPlex, DMPlex.distributeGetDefault, DMPlex.distribute 1746 petsc.DMPlexDistributeSetDefault 1747 1748 """ 1749 cdef PetscBool dist = asBool(flag) 1750 CHKERR(DMPlexDistributeSetDefault(self.dm, dist)) 1751 return 1752 1753 def distributionSetName(self, name: str) -> None: 1754 """Set the name of the specific parallel distribution. 1755 1756 Logically collective. 1757 1758 Parameters 1759 ---------- 1760 name 1761 The name of the specific parallel distribution. 1762 1763 See Also 1764 -------- 1765 DMPlex, DMPlex.distributionGetName, DMPlex.topologyView 1766 DMPlex.topologyLoad, petsc.DMPlexDistributionSetName 1767 1768 """ 1769 cdef const char *cname = NULL 1770 if name is not None: 1771 name = str2bytes(name, &cname) 1772 CHKERR(DMPlexDistributionSetName(self.dm, cname)) 1773 1774 def distributionGetName(self) -> str: 1775 """Retrieve the name of the specific parallel distribution. 1776 1777 Not collective. 1778 1779 Returns 1780 ------- 1781 name : str 1782 The name of the specific parallel distribution. 1783 1784 See Also 1785 -------- 1786 DMPlex, DMPlex.distributionSetName, DMPlex.topologyView 1787 DMPlex.topologyLoad, petsc.DMPlexDistributionGetName 1788 1789 """ 1790 cdef const char *cname = NULL 1791 CHKERR(DMPlexDistributionGetName(self.dm, &cname)) 1792 return bytes2str(cname) 1793 1794 def interpolate(self) -> None: 1795 """Convert to a mesh with all intermediate faces, edges, etc. 1796 1797 Collective. 1798 1799 See Also 1800 -------- 1801 DMPlex, DMPlex.uninterpolate, DMPlex.createFromCellList 1802 petsc.DMPlexInterpolate 1803 1804 """ 1805 cdef PetscDM newdm = NULL 1806 CHKERR(DMPlexInterpolate(self.dm, &newdm)) 1807 CHKERR(PetscCLEAR(self.obj)); self.dm = newdm 1808 1809 def uninterpolate(self) -> None: 1810 """Convert to a mesh with only cells and vertices. 1811 1812 Collective. 1813 1814 See Also 1815 -------- 1816 DMPlex, DMPlex.interpolate, DMPlex.createFromCellList 1817 petsc.DMPlexUninterpolate 1818 1819 """ 1820 cdef PetscDM newdm = NULL 1821 CHKERR(DMPlexUninterpolate(self.dm, &newdm)) 1822 CHKERR(PetscCLEAR(self.obj)); self.dm = newdm 1823 1824 def distributeField(self, SF sf, Section sec, Vec vec, 1825 Section newsec=None, Vec newvec=None) -> tuple[Section, Vec]: 1826 """Distribute field data with a with a given `SF`. 1827 1828 Collective. 1829 1830 Parameters 1831 ---------- 1832 sf 1833 The `SF` describing the communication pattern. 1834 sec 1835 The `Section` for existing data layout. 1836 vec 1837 The existing data in a local vector. 1838 newsec 1839 The `SF` describing the new data layout. 1840 newvec 1841 The new data in a local vector. 1842 1843 Returns 1844 ------- 1845 newSection : Section 1846 The `SF` describing the new data layout. 1847 newVec : Vec 1848 The new data in a local vector. 1849 1850 See Also 1851 -------- 1852 DMPlex, DMPlex.distribute, petsc.DMPlexDistributeField 1853 1854 """ 1855 cdef MPI_Comm ccomm = MPI_COMM_NULL 1856 if newsec is None: newsec = Section() 1857 if newvec is None: newvec = Vec() 1858 if newsec.sec == NULL: 1859 CHKERR(PetscObjectGetComm(<PetscObject>sec.sec, &ccomm)) 1860 CHKERR(PetscSectionCreate(ccomm, &newsec.sec)) 1861 if newvec.vec == NULL: 1862 CHKERR(PetscObjectGetComm(<PetscObject>vec.vec, &ccomm)) 1863 CHKERR(VecCreate(ccomm, &newvec.vec)) 1864 CHKERR(DMPlexDistributeField(self.dm, sf.sf, 1865 sec.sec, vec.vec, 1866 newsec.sec, newvec.vec)) 1867 return (newsec, newvec) 1868 1869 def getMinRadius(self) -> float: 1870 """Return the minimum distance from any cell centroid to a face. 1871 1872 Not collective. 1873 1874 See Also 1875 -------- 1876 DMPlex, DM.getCoordinates, petsc.DMPlexGetMinRadius 1877 1878 """ 1879 cdef PetscReal cminradius = 0. 1880 CHKERR(DMPlexGetMinRadius(self.dm, &cminradius)) 1881 return asReal(cminradius) 1882 1883 def createCoarsePointIS(self) -> IS: 1884 """Create an `IS` covering the coarse `DMPlex` chart with the fine points as data. 1885 1886 Collective. 1887 1888 Returns 1889 ------- 1890 fpointIS : IS 1891 The `IS` of all the fine points which exist in the original coarse mesh. 1892 1893 See Also 1894 -------- 1895 DM, DMPlex, IS, DM.refine, DMPlex.setRefinementUniform 1896 petsc.DMPlexCreateCoarsePointIS 1897 1898 """ 1899 cdef IS fpoint = IS() 1900 CHKERR(DMPlexCreateCoarsePointIS(self.dm, &fpoint.iset)) 1901 return fpoint 1902 1903 def createSection(self, numComp: Sequence[int], numDof: Sequence[int], 1904 bcField: Sequence[int] | None = None, bcComps: Sequence[IS] | None = None, bcPoints: Sequence[IS] | None = None, 1905 IS perm=None) -> Section: 1906 """Create a `Section` based upon the DOF layout specification provided. 1907 1908 Not collective. 1909 1910 Parameters 1911 ---------- 1912 numComp 1913 An array of size ``numFields`` holding the number of components per field. 1914 numDof 1915 An array of size ``numFields*(dim+1)`` holding the number of DOFs 1916 per field on a mesh piece of dimension ``dim``. 1917 bcField 1918 An array of size ``numBC`` giving the field number for each boundary 1919 condition, where ``numBC`` is the number of boundary conditions. 1920 bcComps 1921 An array of size ``numBC`` giving an `IS` holding the field 1922 components to which each boundary condition applies. 1923 bcPoints 1924 An array of size ``numBC`` giving an `IS` holding the `DMPlex` points 1925 to which each boundary condition applies. 1926 perm 1927 Permutation of the chart. 1928 1929 See Also 1930 -------- 1931 DM, DMPlex, DMPlex.create, Section.create, Section.setPermutation 1932 petsc.DMPlexCreateSection 1933 1934 """ 1935 # topological dimension 1936 cdef PetscInt dim = 0 1937 CHKERR(DMGetDimension(self.dm, &dim)) 1938 # components and DOFs 1939 cdef PetscInt ncomp = 0, ndof = 0 1940 cdef PetscInt *icomp = NULL, *idof = NULL 1941 numComp = iarray_i(numComp, &ncomp, &icomp) 1942 numDof = iarray_i(numDof, &ndof, &idof) 1943 assert ndof == ncomp*(dim+1) 1944 # boundary conditions 1945 cdef PetscInt nbc = 0, i = 0 1946 cdef PetscInt *bcfield = NULL 1947 cdef PetscIS *bccomps = NULL 1948 cdef PetscIS *bcpoints = NULL 1949 cdef object unused1, unused2 1950 if bcField is not None: 1951 bcField = iarray_i(bcField, &nbc, &bcfield) 1952 if bcComps is not None: 1953 bcComps = list(bcComps) 1954 assert len(bcComps) == nbc 1955 unused1 = oarray_p(empty_p(nbc), NULL, <void**>&bccomps) 1956 for i from 0 <= i < nbc: 1957 bccomps[i] = (<IS?>bcComps[<Py_ssize_t>i]).iset 1958 if bcPoints is not None: 1959 bcPoints = list(bcPoints) 1960 assert len(bcPoints) == nbc 1961 unused2 = oarray_p(empty_p(nbc), NULL, <void**>&bcpoints) 1962 for i from 0 <= i < nbc: 1963 bcpoints[i] = (<IS?>bcPoints[<Py_ssize_t>i]).iset 1964 else: 1965 raise ValueError("bcPoints is a required argument") 1966 else: 1967 assert bcComps is None 1968 assert bcPoints is None 1969 # optional chart permutations 1970 cdef PetscIS cperm = NULL 1971 if perm is not None: cperm = perm.iset 1972 # create section 1973 cdef Section sec = Section() 1974 CHKERR(DMPlexCreateSection(self.dm, NULL, icomp, idof, 1975 nbc, bcfield, bccomps, bcpoints, 1976 cperm, &sec.sec)) 1977 return sec 1978 1979 def getPointLocal(self, point: int) -> tuple[int, int]: 1980 """Return location of point data in local `Vec`. 1981 1982 Not collective. 1983 1984 Parameters 1985 ---------- 1986 point 1987 The topological point. 1988 1989 Returns 1990 ------- 1991 start : int 1992 Start of point data. 1993 end : int 1994 End of point data. 1995 1996 See Also 1997 -------- 1998 DM, DMPlex, DMPlex.getPointLocalField, Section.getOffset 1999 Section.getDof, petsc.DMPlexGetPointLocal 2000 2001 """ 2002 cdef PetscInt start = 0, end = 0 2003 cdef PetscInt cpoint = asInt(point) 2004 CHKERR(DMPlexGetPointLocal(self.dm, cpoint, &start, &end)) 2005 return toInt(start), toInt(end) 2006 2007 def getPointLocalField(self, point: int, field: int) -> tuple[int, int]: 2008 """Return location of point field data in local `Vec`. 2009 2010 Not collective. 2011 2012 Parameters 2013 ---------- 2014 point 2015 The topological point. 2016 field 2017 The field number. 2018 2019 Returns 2020 ------- 2021 start : int 2022 Start of point data. 2023 end : int 2024 End of point data. 2025 2026 See Also 2027 -------- 2028 DM, DMPlex, DMPlex.getPointLocal, Section.getOffset 2029 petsc.DMPlexGetPointLocalField 2030 2031 """ 2032 cdef PetscInt start = 0, end = 0 2033 cdef PetscInt cpoint = asInt(point) 2034 cdef PetscInt cfield = asInt(field) 2035 CHKERR(DMPlexGetPointLocalField(self.dm, cpoint, cfield, &start, &end)) 2036 return toInt(start), toInt(end) 2037 2038 def getPointGlobal(self, point: int) -> tuple[int, int]: 2039 """Return location of point data in global `Vec`. 2040 2041 Not collective. 2042 2043 Parameters 2044 ---------- 2045 point 2046 The topological point. 2047 2048 Returns 2049 ------- 2050 start : int 2051 Start of point data; returns ``-(globalStart+1)`` if point is not owned. 2052 end : int 2053 End of point data; returns ``-(globalEnd+1)`` if point is not owned. 2054 2055 See Also 2056 -------- 2057 DM, DMPlex, DMPlex.getPointGlobalField, Section.getOffset 2058 Section.getDof, DMPlex.getPointLocal, petsc.DMPlexGetPointGlobal 2059 2060 """ 2061 cdef PetscInt start = 0, end = 0 2062 cdef PetscInt cpoint = asInt(point) 2063 CHKERR(DMPlexGetPointGlobal(self.dm, cpoint, &start, &end)) 2064 return toInt(start), toInt(end) 2065 2066 def getPointGlobalField(self, point: int, field: int) -> tuple[int, int]: 2067 """Return location of point field data in global `Vec`. 2068 2069 Not collective. 2070 2071 Parameters 2072 ---------- 2073 point 2074 The topological point. 2075 field 2076 The field number. 2077 2078 Returns 2079 ------- 2080 start : int 2081 Start of point data; returns ``-(globalStart+1)`` if point is not owned. 2082 end : int 2083 End of point data; returns ``-(globalEnd+1)`` if point is not owned. 2084 2085 See Also 2086 -------- 2087 DM, DMPlex, DMPlex.getPointGlobal, Section.getOffset, Section.getDof 2088 DMPlex.getPointLocal, petsc.DMPlexGetPointGlobalField 2089 2090 """ 2091 cdef PetscInt start = 0, end = 0 2092 cdef PetscInt cpoint = asInt(point) 2093 cdef PetscInt cfield = asInt(field) 2094 CHKERR(DMPlexGetPointGlobalField(self.dm, cpoint, cfield, &start, &end)) 2095 return toInt(start), toInt(end) 2096 2097 def createClosureIndex(self, Section sec or None) -> None: 2098 """Calculate an index for ``sec`` for the closure operation. 2099 2100 Not collective. 2101 2102 Parameters 2103 ---------- 2104 sec 2105 The `Section` describing the layout in the local vector, 2106 or `None` to use the default section. 2107 2108 See Also 2109 -------- 2110 DM, DMPlex, Section, DMPlex.vecGetClosure 2111 petsc.DMPlexCreateClosureIndex 2112 2113 """ 2114 cdef PetscSection csec = sec.sec if sec is not None else NULL 2115 CHKERR(DMPlexCreateClosureIndex(self.dm, csec)) 2116 2117 # 2118 2119 def setRefinementUniform(self, refinementUniform: bool | None = True) -> None: 2120 """Set the flag for uniform refinement. 2121 2122 Logically collective. 2123 2124 Parameters 2125 ---------- 2126 refinementUniform 2127 The flag for uniform refinement. 2128 2129 See Also 2130 -------- 2131 DM, DMPlex, DM.refine, DMPlex.getRefinementUniform 2132 DMPlex.getRefinementLimit, DMPlex.setRefinementLimit 2133 petsc.DMPlexSetRefinementUniform 2134 2135 """ 2136 cdef PetscBool flag = refinementUniform 2137 CHKERR(DMPlexSetRefinementUniform(self.dm, flag)) 2138 2139 def getRefinementUniform(self) -> bool: 2140 """Retrieve the flag for uniform refinement. 2141 2142 Not collective. 2143 2144 Returns 2145 ------- 2146 refinementUniform : bool 2147 The flag for uniform refinement. 2148 2149 See Also 2150 -------- 2151 DM, DMPlex, DM.refine, DMPlex.setRefinementUniform 2152 DMPlex.getRefinementLimit, DMPlex.setRefinementLimit 2153 petsc.DMPlexGetRefinementUniform 2154 2155 """ 2156 cdef PetscBool flag = PETSC_FALSE 2157 CHKERR(DMPlexGetRefinementUniform(self.dm, &flag)) 2158 return toBool(flag) 2159 2160 def setRefinementLimit(self, refinementLimit: float) -> None: 2161 """Set the maximum cell volume for refinement. 2162 2163 Logically collective. 2164 2165 Parameters 2166 ---------- 2167 refinementLimit 2168 The maximum cell volume in the refined mesh. 2169 2170 See Also 2171 -------- 2172 DM, DMPlex, DM.refine, DMPlex.getRefinementLimit 2173 DMPlex.getRefinementUniform, DMPlex.setRefinementUniform 2174 petsc.DMPlexSetRefinementLimit 2175 2176 """ 2177 cdef PetscReal rval = asReal(refinementLimit) 2178 CHKERR(DMPlexSetRefinementLimit(self.dm, rval)) 2179 2180 def getRefinementLimit(self) -> float: 2181 """Retrieve the maximum cell volume for refinement. 2182 2183 Not collective. 2184 2185 See Also 2186 -------- 2187 DM, DMPlex, DM.refine, DMPlex.setRefinementLimit 2188 DMPlex.getRefinementUniform, DMPlex.setRefinementUniform 2189 petsc.DMPlexGetRefinementLimit 2190 2191 """ 2192 cdef PetscReal rval = 0.0 2193 CHKERR(DMPlexGetRefinementLimit(self.dm, &rval)) 2194 return toReal(rval) 2195 2196 def getOrdering(self, otype: Mat.OrderingType) -> IS: 2197 """Calculate a reordering of the mesh. 2198 2199 Collective. 2200 2201 Parameters 2202 ---------- 2203 otype 2204 Type of reordering, see `Mat.OrderingType`. 2205 2206 Returns 2207 ------- 2208 perm : IS 2209 The point permutation. 2210 2211 See Also 2212 -------- 2213 DMPlex, DMPlex.permute, Mat.OrderingType, Mat.getOrdering 2214 petsc.DMPlexGetOrdering 2215 2216 """ 2217 cdef PetscMatOrderingType cval = NULL 2218 cdef PetscDMLabel label = NULL 2219 otype = str2bytes(otype, &cval) 2220 cdef IS perm = IS() 2221 CHKERR(DMPlexGetOrdering(self.dm, cval, label, &perm.iset)) 2222 return perm 2223 2224 def permute(self, IS perm) -> DMPlex: 2225 """Reorder the mesh according to the input permutation. 2226 2227 Collective. 2228 2229 Parameters 2230 ---------- 2231 perm 2232 The point permutation, ``perm[old point number] = new point number``. 2233 2234 Returns 2235 ------- 2236 pdm : DMPlex 2237 The permuted `DMPlex`. 2238 2239 See Also 2240 -------- 2241 DMPlex, Mat.permute, petsc.DMPlexPermute 2242 2243 """ 2244 cdef DMPlex dm = <DMPlex>type(self)() 2245 CHKERR(DMPlexPermute(self.dm, perm.iset, &dm.dm)) 2246 return dm 2247 2248 def reorderGetDefault(self) -> DM.ReorderDefaultFlag: 2249 """Return flag indicating whether the `DMPlex` should be reordered by default. 2250 2251 Not collective. 2252 2253 See Also 2254 -------- 2255 `DMPlex.reorderSetDefault`, petsc.DMPlexReorderGetDefault 2256 2257 """ 2258 cdef PetscDMReorderDefaultFlag reorder = DM_REORDER_DEFAULT_NOTSET 2259 CHKERR(DMPlexReorderGetDefault(self.dm, &reorder)) 2260 return reorder 2261 2262 def reorderSetDefault(self, flag: DM.ReorderDefaultFlag) -> None: 2263 """Set flag indicating whether the DM should be reordered by default. 2264 2265 Logically collective. 2266 2267 Parameters 2268 ---------- 2269 reorder 2270 Flag for reordering. 2271 2272 See Also 2273 -------- 2274 DMPlex.reorderGetDefault, petsc.DMPlexReorderSetDefault 2275 2276 """ 2277 cdef PetscDMReorderDefaultFlag reorder = flag 2278 CHKERR(DMPlexReorderSetDefault(self.dm, reorder)) 2279 return 2280 2281 # 2282 2283 def computeCellGeometryFVM(self, cell: int) -> tuple[float, ArrayReal, ArrayReal]: 2284 """Compute the volume for a given cell. 2285 2286 Not collective. 2287 2288 Parameters 2289 ---------- 2290 cell 2291 The cell. 2292 2293 Returns 2294 ------- 2295 volume : float 2296 The cell volume. 2297 centroid : ArrayReal 2298 The cell centroid. 2299 normal : ArrayReal 2300 The cell normal, if appropriate. 2301 2302 See Also 2303 -------- 2304 DMPlex, DM.getCoordinateSection, DM.getCoordinates 2305 petsc.DMPlexComputeCellGeometryFVM 2306 2307 """ 2308 cdef PetscInt cdim = 0 2309 cdef PetscInt ccell = asInt(cell) 2310 CHKERR(DMGetCoordinateDim(self.dm, &cdim)) 2311 cdef PetscReal vol = 0, centroid[3], normal[3] 2312 CHKERR(DMPlexComputeCellGeometryFVM(self.dm, ccell, &vol, centroid, normal)) 2313 return (toReal(vol), array_r(cdim, centroid), array_r(cdim, normal)) 2314 2315 def constructGhostCells(self, labelName: str | None = None) -> int: 2316 """Construct ghost cells which connect to every boundary face. 2317 2318 Collective. 2319 2320 Parameters 2321 ---------- 2322 labelName 2323 The name of the label specifying the boundary faces. 2324 Defaults to ``"Face Sets"``. 2325 2326 Returns 2327 ------- 2328 numGhostCells : int 2329 The number of ghost cells added to the `DMPlex`. 2330 2331 See Also 2332 -------- 2333 DM, DMPlex, DM.create, petsc.DMPlexConstructGhostCells 2334 2335 """ 2336 cdef const char *cname = NULL 2337 labelName = str2bytes(labelName, &cname) 2338 cdef PetscInt numGhostCells = 0 2339 cdef PetscDM dmGhosted = NULL 2340 CHKERR(DMPlexConstructGhostCells(self.dm, cname, &numGhostCells, &dmGhosted)) 2341 CHKERR(PetscCLEAR(self.obj)); self.dm = dmGhosted 2342 return toInt(numGhostCells) 2343 2344 def getSubpointIS(self) -> IS: 2345 """Return an `IS` covering the entire subdm chart. 2346 2347 Not collective. 2348 2349 Returns 2350 ------- 2351 iset : IS 2352 The `IS` containing subdm's parent's points. 2353 2354 See Also 2355 -------- 2356 DM, DMPlex, petsc.DMPlexGetSubpointIS 2357 2358 """ 2359 cdef IS iset = IS() 2360 CHKERR(DMPlexGetSubpointIS(self.dm, &iset.iset)) 2361 PetscINCREF(iset.obj) 2362 return iset 2363 2364 def getSubpointMap(self) -> DMLabel: 2365 """Return a `DMLabel` with point dimension as values. 2366 2367 Not collective. 2368 2369 Returns 2370 ------- 2371 label : DMLabel 2372 The `DMLabel` whose values are subdm's point dimensions. 2373 2374 See Also 2375 -------- 2376 DM, DMPlex, petsc.DMPlexGetSubpointMap 2377 2378 """ 2379 cdef DMLabel label = DMLabel() 2380 CHKERR(DMPlexGetSubpointMap(self.dm, &label.dmlabel)) 2381 PetscINCREF(label.obj) 2382 return label 2383 2384 # Metric 2385 2386 def metricSetFromOptions(self) -> None: 2387 """Configure the object from the options database. 2388 2389 Collective. 2390 2391 See Also 2392 -------- 2393 petsc_options 2394 2395 """ 2396 # FIXME petsc.DMPlexMetricSetFromOptions 2397 CHKERR(DMPlexMetricSetFromOptions(self.dm)) 2398 2399 def metricSetUniform(self, uniform: bool) -> None: 2400 """Record whether the metric is uniform or not. 2401 2402 Logically collective. 2403 2404 Parameters 2405 ---------- 2406 uniform 2407 Flag indicating whether the metric is uniform or not. 2408 2409 See Also 2410 -------- 2411 DMPlex.metricIsUniform, DMPlex.metricSetIsotropic 2412 DMPlex.metricSetRestrictAnisotropyFirst, petsc.DMPlexMetricSetUniform 2413 2414 """ 2415 cdef PetscBool bval = asBool(uniform) 2416 CHKERR(DMPlexMetricSetUniform(self.dm, bval)) 2417 2418 def metricIsUniform(self) -> bool: 2419 """Return the flag indicating whether the metric is uniform or not. 2420 2421 Not collective. 2422 2423 See Also 2424 -------- 2425 DMPlex.metricSetUniform, DMPlex.metricRestrictAnisotropyFirst 2426 petsc.DMPlexMetricIsUniform 2427 2428 """ 2429 cdef PetscBool uniform = PETSC_FALSE 2430 CHKERR(DMPlexMetricIsUniform(self.dm, &uniform)) 2431 return toBool(uniform) 2432 2433 def metricSetIsotropic(self, isotropic: bool) -> None: 2434 """Record whether the metric is isotropic or not. 2435 2436 Logically collective. 2437 2438 Parameters 2439 ---------- 2440 isotropic 2441 Flag indicating whether the metric is isotropic or not. 2442 2443 See Also 2444 -------- 2445 DMPlex.metricIsIsotropic, DMPlex.metricSetUniform 2446 DMPlex.metricSetRestrictAnisotropyFirst, petsc.DMPlexMetricSetIsotropic 2447 2448 """ 2449 cdef PetscBool bval = asBool(isotropic) 2450 CHKERR(DMPlexMetricSetIsotropic(self.dm, bval)) 2451 2452 def metricIsIsotropic(self) -> bool: 2453 """Return the flag indicating whether the metric is isotropic or not. 2454 2455 Not collective. 2456 2457 See Also 2458 -------- 2459 DMPlex.metricSetIsotropic, DMPlex.metricIsUniform 2460 DMPlex.metricRestrictAnisotropyFirst, petsc.DMPlexMetricIsIsotropic 2461 2462 """ 2463 cdef PetscBool isotropic = PETSC_FALSE 2464 CHKERR(DMPlexMetricIsIsotropic(self.dm, &isotropic)) 2465 return toBool(isotropic) 2466 2467 def metricSetRestrictAnisotropyFirst(self, restrictAnisotropyFirst: bool) -> None: 2468 """Record whether anisotropy is be restricted before normalization or after. 2469 2470 Logically collective. 2471 2472 Parameters 2473 ---------- 2474 restrictAnisotropyFirst 2475 Flag indicating if anisotropy is restricted before normalization or after. 2476 2477 See Also 2478 -------- 2479 DMPlex.metricSetIsotropic, DMPlex.metricRestrictAnisotropyFirst 2480 petsc.DMPlexMetricSetRestrictAnisotropyFirst 2481 2482 """ 2483 cdef PetscBool bval = asBool(restrictAnisotropyFirst) 2484 CHKERR(DMPlexMetricSetRestrictAnisotropyFirst(self.dm, bval)) 2485 2486 def metricRestrictAnisotropyFirst(self) -> bool: 2487 """Return ``true`` if anisotropy is restricted before normalization. 2488 2489 Not collective. 2490 2491 See Also 2492 -------- 2493 DMPlex.metricIsIsotropic, DMPlex.metricSetRestrictAnisotropyFirst 2494 petsc.DMPlexMetricRestrictAnisotropyFirst 2495 2496 """ 2497 cdef PetscBool restrictAnisotropyFirst = PETSC_FALSE 2498 CHKERR(DMPlexMetricRestrictAnisotropyFirst(self.dm, &restrictAnisotropyFirst)) 2499 return toBool(restrictAnisotropyFirst) 2500 2501 def metricSetNoInsertion(self, noInsert: bool) -> None: 2502 """Set the flag indicating whether node insertion should be turned off. 2503 2504 Logically collective. 2505 2506 Parameters 2507 ---------- 2508 noInsert 2509 Flag indicating whether node insertion and deletion should be turned off. 2510 2511 See Also 2512 -------- 2513 DMPlex.metricNoInsertion, DMPlex.metricSetNoSwapping 2514 DMPlex.metricSetNoMovement, DMPlex.metricSetNoSurf 2515 petsc.DMPlexMetricSetNoInsertion 2516 2517 """ 2518 cdef PetscBool bval = asBool(noInsert) 2519 CHKERR(DMPlexMetricSetNoInsertion(self.dm, bval)) 2520 2521 def metricNoInsertion(self) -> bool: 2522 """Return the flag indicating whether node insertion and deletion are turned off. 2523 2524 Not collective. 2525 2526 See Also 2527 -------- 2528 DMPlex.metricSetNoInsertion, DMPlex.metricNoSwapping 2529 DMPlex.metricNoMovement, DMPlex.metricNoSurf 2530 petsc.DMPlexMetricNoInsertion 2531 2532 """ 2533 cdef PetscBool noInsert = PETSC_FALSE 2534 CHKERR(DMPlexMetricNoInsertion(self.dm, &noInsert)) 2535 return toBool(noInsert) 2536 2537 def metricSetNoSwapping(self, noSwap: bool) -> None: 2538 """Set the flag indicating whether facet swapping should be turned off. 2539 2540 Logically collective. 2541 2542 Parameters 2543 ---------- 2544 noSwap 2545 Flag indicating whether facet swapping should be turned off. 2546 2547 See Also 2548 -------- 2549 DMPlex.metricNoSwapping, DMPlex.metricSetNoInsertion 2550 DMPlex.metricSetNoMovement, DMPlex.metricSetNoSurf 2551 petsc.DMPlexMetricSetNoSwapping 2552 2553 """ 2554 cdef PetscBool bval = asBool(noSwap) 2555 CHKERR(DMPlexMetricSetNoSwapping(self.dm, bval)) 2556 2557 def metricNoSwapping(self) -> bool: 2558 """Return the flag indicating whether facet swapping is turned off. 2559 2560 Not collective. 2561 2562 See Also 2563 -------- 2564 DMPlex.metricSetNoSwapping, DMPlex.metricNoInsertion 2565 DMPlex.metricNoMovement, DMPlex.metricNoSurf 2566 petsc.DMPlexMetricNoSwapping 2567 2568 """ 2569 cdef PetscBool noSwap = PETSC_FALSE 2570 CHKERR(DMPlexMetricNoSwapping(self.dm, &noSwap)) 2571 return toBool(noSwap) 2572 2573 def metricSetNoMovement(self, noMove: bool) -> None: 2574 """Set the flag indicating whether node movement should be turned off. 2575 2576 Logically collective. 2577 2578 Parameters 2579 ---------- 2580 noMove 2581 Flag indicating whether node movement should be turned off. 2582 2583 See Also 2584 -------- 2585 DMPlex.metricNoMovement, DMPlex.metricSetNoInsertion 2586 DMPlex.metricSetNoSwapping, DMPlex.metricSetNoSurf 2587 petsc.DMPlexMetricSetNoMovement 2588 2589 """ 2590 cdef PetscBool bval = asBool(noMove) 2591 CHKERR(DMPlexMetricSetNoMovement(self.dm, bval)) 2592 2593 def metricNoMovement(self) -> bool: 2594 """Return the flag indicating whether node movement is turned off. 2595 2596 Not collective. 2597 2598 See Also 2599 -------- 2600 DMPlex.metricSetNoMovement, DMPlex.metricNoInsertion 2601 DMPlex.metricNoSwapping, DMPlex.metricNoSurf 2602 petsc.DMPlexMetricNoMovement 2603 2604 """ 2605 cdef PetscBool noMove = PETSC_FALSE 2606 CHKERR(DMPlexMetricNoMovement(self.dm, &noMove)) 2607 return toBool(noMove) 2608 2609 def metricSetNoSurf(self, noSurf: bool) -> None: 2610 """Set the flag indicating whether surface modification should be turned off. 2611 2612 Logically collective. 2613 2614 Parameters 2615 ---------- 2616 noSurf 2617 Flag indicating whether surface modification should be turned off. 2618 2619 See Also 2620 -------- 2621 DMPlex.metricNoSurf, DMPlex.metricSetNoMovement 2622 DMPlex.metricSetNoInsertion, DMPlex.metricSetNoSwapping 2623 petsc.DMPlexMetricSetNoSurf 2624 2625 """ 2626 cdef PetscBool bval = asBool(noSurf) 2627 CHKERR(DMPlexMetricSetNoSurf(self.dm, bval)) 2628 2629 def metricNoSurf(self) -> bool: 2630 """Return the flag indicating whether surface modification is turned off. 2631 2632 Not collective. 2633 2634 See Also 2635 -------- 2636 DMPlex.metricSetNoSurf, DMPlex.metricNoMovement 2637 DMPlex.metricNoInsertion, DMPlex.metricNoSwapping 2638 petsc.DMPlexMetricNoSurf 2639 2640 """ 2641 cdef PetscBool noSurf = PETSC_FALSE 2642 CHKERR(DMPlexMetricNoSurf(self.dm, &noSurf)) 2643 return toBool(noSurf) 2644 2645 def metricSetVerbosity(self, verbosity: int) -> None: 2646 """Set the verbosity of the mesh adaptation package. 2647 2648 Logically collective. 2649 2650 Parameters 2651 ---------- 2652 verbosity 2653 The verbosity, where -1 is silent and 10 is maximum. 2654 2655 See Also 2656 -------- 2657 DMPlex.metricGetVerbosity, DMPlex.metricSetNumIterations 2658 petsc.DMPlexMetricSetVerbosity 2659 2660 """ 2661 cdef PetscInt ival = asInt(verbosity) 2662 CHKERR(DMPlexMetricSetVerbosity(self.dm, ival)) 2663 2664 def metricGetVerbosity(self) -> int: 2665 """Return the verbosity of the mesh adaptation package. 2666 2667 Not collective. 2668 2669 Returns 2670 ------- 2671 verbosity : int 2672 The verbosity, where -1 is silent and 10 is maximum. 2673 2674 See Also 2675 -------- 2676 DMPlex.metricSetVerbosity, DMPlex.metricGetNumIterations 2677 petsc.DMPlexMetricGetVerbosity 2678 2679 """ 2680 cdef PetscInt verbosity = 0 2681 CHKERR(DMPlexMetricGetVerbosity(self.dm, &verbosity)) 2682 return toInt(verbosity) 2683 2684 def metricSetNumIterations(self, numIter: int) -> None: 2685 """Set the number of parallel adaptation iterations. 2686 2687 Logically collective. 2688 2689 Parameters 2690 ---------- 2691 numIter 2692 The number of parallel adaptation iterations. 2693 2694 See Also 2695 -------- 2696 DMPlex.metricSetVerbosity, DMPlex.metricGetNumIterations 2697 petsc.DMPlexMetricSetNumIterations 2698 2699 """ 2700 cdef PetscInt ival = asInt(numIter) 2701 CHKERR(DMPlexMetricSetNumIterations(self.dm, ival)) 2702 2703 def metricGetNumIterations(self) -> int: 2704 """Return the number of parallel adaptation iterations. 2705 2706 Not collective. 2707 2708 See Also 2709 -------- 2710 DMPlex.metricSetNumIterations, DMPlex.metricGetVerbosity 2711 petsc.DMPlexMetricGetNumIterations 2712 2713 """ 2714 cdef PetscInt numIter = 0 2715 CHKERR(DMPlexMetricGetNumIterations(self.dm, &numIter)) 2716 return toInt(numIter) 2717 2718 def metricSetMinimumMagnitude(self, h_min: float) -> None: 2719 """Set the minimum tolerated metric magnitude. 2720 2721 Logically collective. 2722 2723 Parameters 2724 ---------- 2725 h_min 2726 The minimum tolerated metric magnitude. 2727 2728 See Also 2729 -------- 2730 DMPlex.metricGetMinimumMagnitude, DMPlex.metricSetMaximumMagnitude 2731 petsc.DMPlexMetricSetMinimumMagnitude 2732 2733 """ 2734 cdef PetscReal rval = asReal(h_min) 2735 CHKERR(DMPlexMetricSetMinimumMagnitude(self.dm, rval)) 2736 2737 def metricGetMinimumMagnitude(self) -> float: 2738 """Return the minimum tolerated metric magnitude. 2739 2740 Not collective. 2741 2742 See Also 2743 -------- 2744 DMPlex.metricSetMinimumMagnitude, DMPlex.metricGetMaximumMagnitude 2745 petsc.DMPlexMetricGetMinimumMagnitude 2746 2747 """ 2748 cdef PetscReal h_min = 0 2749 CHKERR(DMPlexMetricGetMinimumMagnitude(self.dm, &h_min)) 2750 return toReal(h_min) 2751 2752 def metricSetMaximumMagnitude(self, h_max: float) -> None: 2753 """Set the maximum tolerated metric magnitude. 2754 2755 Logically collective. 2756 2757 Parameters 2758 ---------- 2759 h_max 2760 The maximum tolerated metric magnitude. 2761 2762 See Also 2763 -------- 2764 DMPlex.metricGetMaximumMagnitude, DMPlex.metricSetMinimumMagnitude 2765 petsc.DMPlexMetricSetMaximumMagnitude 2766 2767 """ 2768 cdef PetscReal rval = asReal(h_max) 2769 CHKERR(DMPlexMetricSetMaximumMagnitude(self.dm, rval)) 2770 2771 def metricGetMaximumMagnitude(self) -> float: 2772 """Return the maximum tolerated metric magnitude. 2773 2774 Not collective. 2775 2776 See Also 2777 -------- 2778 DMPlex.metricSetMaximumMagnitude, DMPlex.metricGetMinimumMagnitude 2779 petsc.DMPlexMetricGetMaximumMagnitude 2780 2781 """ 2782 cdef PetscReal h_max = 0 2783 CHKERR(DMPlexMetricGetMaximumMagnitude(self.dm, &h_max)) 2784 return toReal(h_max) 2785 2786 def metricSetMaximumAnisotropy(self, a_max: float) -> None: 2787 """Set the maximum tolerated metric anisotropy. 2788 2789 Logically collective. 2790 2791 Parameters 2792 ---------- 2793 a_max 2794 The maximum tolerated metric anisotropy. 2795 2796 See Also 2797 -------- 2798 DMPlex.metricGetMaximumAnisotropy, DMPlex.metricSetMaximumMagnitude 2799 petsc.DMPlexMetricSetMaximumAnisotropy 2800 2801 """ 2802 cdef PetscReal rval = asReal(a_max) 2803 CHKERR(DMPlexMetricSetMaximumAnisotropy(self.dm, rval)) 2804 2805 def metricGetMaximumAnisotropy(self) -> float: 2806 """Return the maximum tolerated metric anisotropy. 2807 2808 Not collective. 2809 2810 See Also 2811 -------- 2812 DMPlex.metricSetMaximumAnisotropy, DMPlex.metricGetMaximumMagnitude 2813 petsc.DMPlexMetricGetMaximumAnisotropy 2814 2815 """ 2816 cdef PetscReal a_max = 0 2817 CHKERR(DMPlexMetricGetMaximumAnisotropy(self.dm, &a_max)) 2818 return toReal(a_max) 2819 2820 def metricSetTargetComplexity(self, targetComplexity: float) -> None: 2821 """Set the target metric complexity. 2822 2823 Logically collective. 2824 2825 Parameters 2826 ---------- 2827 targetComplexity 2828 The target metric complexity. 2829 2830 See Also 2831 -------- 2832 DMPlex.metricGetTargetComplexity, DMPlex.metricSetNormalizationOrder 2833 petsc.DMPlexMetricSetTargetComplexity 2834 2835 """ 2836 cdef PetscReal rval = asReal(targetComplexity) 2837 CHKERR(DMPlexMetricSetTargetComplexity(self.dm, rval)) 2838 2839 def metricGetTargetComplexity(self) -> float: 2840 """Return the target metric complexity. 2841 2842 Not collective. 2843 2844 See Also 2845 -------- 2846 DMPlex.metricSetTargetComplexity, DMPlex.metricGetNormalizationOrder 2847 petsc.DMPlexMetricGetTargetComplexity 2848 2849 """ 2850 cdef PetscReal targetComplexity = 0 2851 CHKERR(DMPlexMetricGetTargetComplexity(self.dm, &targetComplexity)) 2852 return toReal(targetComplexity) 2853 2854 def metricSetNormalizationOrder(self, p: float) -> None: 2855 """Set the order p for L-p normalization. 2856 2857 Logically collective. 2858 2859 Parameters 2860 ---------- 2861 p 2862 The normalization order. 2863 2864 See Also 2865 -------- 2866 DMPlex.metricGetNormalizationOrder, DMPlex.metricSetTargetComplexity 2867 petsc.DMPlexMetricSetNormalizationOrder 2868 2869 """ 2870 cdef PetscReal rval = asReal(p) 2871 CHKERR(DMPlexMetricSetNormalizationOrder(self.dm, rval)) 2872 2873 def metricGetNormalizationOrder(self) -> float: 2874 """Return the order p for L-p normalization. 2875 2876 Not collective. 2877 2878 See Also 2879 -------- 2880 DMPlex.metricSetNormalizationOrder, DMPlex.metricGetTargetComplexity 2881 petsc.DMPlexMetricGetNormalizationOrder 2882 2883 """ 2884 cdef PetscReal p = 0 2885 CHKERR(DMPlexMetricGetNormalizationOrder(self.dm, &p)) 2886 return toReal(p) 2887 2888 def metricSetGradationFactor(self, beta: float) -> None: 2889 """Set the metric gradation factor. 2890 2891 Logically collective. 2892 2893 Parameters 2894 ---------- 2895 beta 2896 The metric gradation factor. 2897 2898 See Also 2899 -------- 2900 DMPlex.metricGetGradationFactor, DMPlex.metricSetHausdorffNumber 2901 petsc.DMPlexMetricSetGradationFactor 2902 2903 """ 2904 cdef PetscReal rval = asReal(beta) 2905 CHKERR(DMPlexMetricSetGradationFactor(self.dm, rval)) 2906 2907 def metricGetGradationFactor(self) -> float: 2908 """Return the metric gradation factor. 2909 2910 Not collective. 2911 2912 See Also 2913 -------- 2914 DMPlex.metricSetGradationFactor, DMPlex.metricGetHausdorffNumber 2915 petsc.DMPlexMetricGetGradationFactor 2916 2917 """ 2918 cdef PetscReal beta = 0 2919 CHKERR(DMPlexMetricGetGradationFactor(self.dm, &beta)) 2920 return toReal(beta) 2921 2922 def metricSetHausdorffNumber(self, hausd: float) -> None: 2923 """Set the metric Hausdorff number. 2924 2925 Logically collective. 2926 2927 Parameters 2928 ---------- 2929 hausd 2930 The metric Hausdorff number. 2931 2932 See Also 2933 -------- 2934 DMPlex.metricSetGradationFactor, DMPlex.metricGetHausdorffNumber 2935 petsc.DMPlexMetricSetHausdorffNumber 2936 2937 """ 2938 cdef PetscReal rval = asReal(hausd) 2939 CHKERR(DMPlexMetricSetHausdorffNumber(self.dm, rval)) 2940 2941 def metricGetHausdorffNumber(self) -> float: 2942 """Return the metric Hausdorff number. 2943 2944 Not collective. 2945 2946 See Also 2947 -------- 2948 DMPlex.metricGetGradationFactor, DMPlex.metricSetHausdorffNumber 2949 petsc.DMPlexMetricGetHausdorffNumber 2950 2951 """ 2952 cdef PetscReal hausd = 0 2953 CHKERR(DMPlexMetricGetHausdorffNumber(self.dm, &hausd)) 2954 return toReal(hausd) 2955 2956 def metricCreate(self, field: int | None = 0) -> Vec: 2957 """Create a Riemannian metric field. 2958 2959 Collective. 2960 2961 Parameters 2962 ---------- 2963 field 2964 The field number to use. 2965 2966 See Also 2967 -------- 2968 DMPlex.metricCreateUniform, DMPlex.metricCreateIsotropic 2969 petsc_options, petsc.DMPlexMetricCreate 2970 2971 """ 2972 cdef PetscInt ival = asInt(field) 2973 cdef Vec metric = Vec() 2974 CHKERR(DMPlexMetricCreate(self.dm, ival, &metric.vec)) 2975 return metric 2976 2977 def metricCreateUniform(self, alpha: float, field: int | None = 0) -> Vec: 2978 """Construct a uniform isotropic metric. 2979 2980 Collective. 2981 2982 Parameters 2983 ---------- 2984 alpha 2985 Scaling parameter for the diagonal. 2986 field 2987 The field number to use. 2988 2989 See Also 2990 -------- 2991 DMPlex.metricCreate, DMPlex.metricCreateIsotropic 2992 petsc.DMPlexMetricCreateUniform 2993 2994 """ 2995 cdef PetscInt ival = asInt(field) 2996 cdef PetscReal rval = asReal(alpha) 2997 cdef Vec metric = Vec() 2998 CHKERR(DMPlexMetricCreateUniform(self.dm, ival, rval, &metric.vec)) 2999 return metric 3000 3001 def metricCreateIsotropic(self, Vec indicator, field: int | None = 0) -> Vec: 3002 """Construct an isotropic metric from an error indicator. 3003 3004 Collective. 3005 3006 Parameters 3007 ---------- 3008 indicator 3009 The error indicator. 3010 field 3011 The field number to use. 3012 3013 See Also 3014 -------- 3015 DMPlex.metricCreate, DMPlex.metricCreateUniform 3016 petsc.DMPlexMetricCreateIsotropic 3017 3018 """ 3019 cdef PetscInt ival = asInt(field) 3020 cdef Vec metric = Vec() 3021 CHKERR(DMPlexMetricCreateIsotropic(self.dm, ival, indicator.vec, &metric.vec)) 3022 return metric 3023 3024 def metricDeterminantCreate(self, field: int | None = 0) -> tuple[Vec, DM]: 3025 """Create the determinant field for a Riemannian metric. 3026 3027 Collective. 3028 3029 Parameters 3030 ---------- 3031 field 3032 The field number to use. 3033 3034 Returns 3035 ------- 3036 determinant : Vec 3037 The determinant field. 3038 dmDet : DM 3039 The corresponding DM 3040 3041 See Also 3042 -------- 3043 DMPlex.metricCreateUniform, DMPlex.metricCreateIsotropic 3044 DMPlex.metricCreate, petsc.DMPlexMetricDeterminantCreate 3045 3046 """ 3047 cdef PetscInt ival = asInt(field) 3048 cdef Vec determinant = Vec() 3049 cdef DM dmDet = DM() 3050 CHKERR(DMPlexMetricDeterminantCreate(self.dm, ival, &determinant.vec, &dmDet.dm)) 3051 return (determinant, dmDet) 3052 3053 def metricEnforceSPD(self, Vec metric, Vec ometric, Vec determinant, restrictSizes: bool | None = False, restrictAnisotropy: bool | None = False) -> tuple[Vec, Vec]: 3054 """Enforce symmetric positive-definiteness of a metric. 3055 3056 Collective. 3057 3058 Parameters 3059 ---------- 3060 metric 3061 The metric. 3062 ometric 3063 The output metric. 3064 determinant 3065 The output determinant. 3066 restrictSizes 3067 Flag indicating whether maximum/minimum magnitudes should be enforced. 3068 restrictAnisotropy 3069 Flag indicating whether maximum anisotropy should be enforced. 3070 3071 Returns 3072 ------- 3073 ometric : Vec 3074 The output metric. 3075 determinant : Vec 3076 The output determinant. 3077 3078 See Also 3079 -------- 3080 DMPlex.metricNormalize, DMPlex.metricIntersection2 3081 DMPlex.metricIntersection3, petsc_options, petsc.DMPlexMetricEnforceSPD 3082 3083 """ 3084 cdef PetscBool bval_rs = asBool(restrictSizes) 3085 cdef PetscBool bval_ra = asBool(restrictAnisotropy) 3086 CHKERR(DMPlexMetricEnforceSPD(self.dm, metric.vec, bval_rs, bval_ra, ometric.vec, determinant.vec)) 3087 return (ometric, determinant) 3088 3089 def metricNormalize(self, Vec metric, Vec ometric, Vec determinant, restrictSizes: bool | None = True, restrictAnisotropy: bool | None = True) -> tuple[Vec, Vec]: 3090 """Apply L-p normalization to a metric. 3091 3092 Collective. 3093 3094 Parameters 3095 ---------- 3096 metric 3097 The metric. 3098 ometric 3099 The output metric. 3100 determinant 3101 The output determinant. 3102 restrictSizes 3103 Flag indicating whether maximum/minimum magnitudes should be enforced. 3104 restrictAnisotropy 3105 Flag indicating whether maximum anisotropy should be enforced. 3106 3107 Returns 3108 ------- 3109 ometric : Vec 3110 The output normalized metric. 3111 determinant : Vec 3112 The output determinant. 3113 3114 See Also 3115 -------- 3116 DMPlex.metricEnforceSPD, DMPlex.metricIntersection2 3117 DMPlex.metricIntersection3, petsc_options, petsc.DMPlexMetricNormalize 3118 3119 """ 3120 cdef PetscBool bval_rs = asBool(restrictSizes) 3121 cdef PetscBool bval_ra = asBool(restrictAnisotropy) 3122 CHKERR(DMPlexMetricNormalize(self.dm, metric.vec, bval_rs, bval_ra, ometric.vec, determinant.vec)) 3123 return (ometric, determinant) 3124 3125 def metricAverage2(self, Vec metric1, Vec metric2, Vec metricAvg) -> Vec: 3126 """Compute and return the unweighted average of two metrics. 3127 3128 Collective. 3129 3130 Parameters 3131 ---------- 3132 metric1 3133 The first metric to be averaged. 3134 metric2 3135 The second metric to be averaged. 3136 metricAvg 3137 The output averaged metric. 3138 3139 See Also 3140 -------- 3141 DMPlex.metricAverage3, petsc.DMPlexMetricAverage2 3142 3143 """ 3144 CHKERR(DMPlexMetricAverage2(self.dm, metric1.vec, metric2.vec, metricAvg.vec)) 3145 return metricAvg 3146 3147 def metricAverage3(self, Vec metric1, Vec metric2, Vec metric3, Vec metricAvg) -> Vec: 3148 """Compute and return the unweighted average of three metrics. 3149 3150 Collective. 3151 3152 Parameters 3153 ---------- 3154 metric1 3155 The first metric to be averaged. 3156 metric2 3157 The second metric to be averaged. 3158 metric3 3159 The third metric to be averaged. 3160 metricAvg 3161 The output averaged metric. 3162 3163 See Also 3164 -------- 3165 DMPlex.metricAverage2, petsc.DMPlexMetricAverage3 3166 3167 """ 3168 CHKERR(DMPlexMetricAverage3(self.dm, metric1.vec, metric2.vec, metric3.vec, metricAvg.vec)) 3169 return metricAvg 3170 3171 def metricIntersection2(self, Vec metric1, Vec metric2, Vec metricInt) -> Vec: 3172 """Compute and return the intersection of two metrics. 3173 3174 Collective. 3175 3176 Parameters 3177 ---------- 3178 metric1 3179 The first metric to be intersected. 3180 metric2 3181 The second metric to be intersected. 3182 metricInt 3183 The output intersected metric. 3184 3185 See Also 3186 -------- 3187 DMPlex.metricIntersection3, petsc.DMPlexMetricIntersection2 3188 3189 """ 3190 CHKERR(DMPlexMetricIntersection2(self.dm, metric1.vec, metric2.vec, metricInt.vec)) 3191 return metricInt 3192 3193 def metricIntersection3(self, Vec metric1, Vec metric2, Vec metric3, Vec metricInt) -> Vec: 3194 """Compute the intersection of three metrics. 3195 3196 Collective. 3197 3198 Parameters 3199 ---------- 3200 metric1 3201 The first metric to be intersected. 3202 metric2 3203 The second metric to be intersected. 3204 metric3 3205 The third metric to be intersected. 3206 metricInt 3207 The output intersected metric. 3208 3209 See Also 3210 -------- 3211 DMPlex.metricIntersection2, petsc.DMPlexMetricIntersection3 3212 3213 """ 3214 CHKERR(DMPlexMetricIntersection3(self.dm, metric1.vec, metric2.vec, metric3.vec, metricInt.vec)) 3215 return metricInt 3216 3217 def computeGradientClementInterpolant(self, Vec locX, Vec locC) -> Vec: 3218 """Return the L2 projection of the cellwise gradient of a function onto P1. 3219 3220 Collective. 3221 3222 Parameters 3223 ---------- 3224 locX 3225 The coefficient vector of the function. 3226 locC 3227 The output `Vec` which holds the Clement interpolant of the gradient. 3228 3229 See Also 3230 -------- 3231 DM, DMPlex, petsc.DMPlexComputeGradientClementInterpolant 3232 3233 """ 3234 CHKERR(DMPlexComputeGradientClementInterpolant(self.dm, locX.vec, locC.vec)) 3235 return locC 3236 3237 # View 3238 3239 def topologyView(self, Viewer viewer) -> None: 3240 """Save a `DMPlex` topology into a file. 3241 3242 Collective. 3243 3244 Parameters 3245 ---------- 3246 viewer 3247 The `Viewer` for saving. 3248 3249 See Also 3250 -------- 3251 DM, DMPlex, DM.view, DMPlex.coordinatesView, DMPlex.labelsView 3252 DMPlex.topologyLoad, Viewer, petsc.DMPlexTopologyView 3253 3254 """ 3255 CHKERR(DMPlexTopologyView(self.dm, viewer.vwr)) 3256 3257 def coordinatesView(self, Viewer viewer) -> None: 3258 """Save `DMPlex` coordinates into a file. 3259 3260 Collective. 3261 3262 Parameters 3263 ---------- 3264 viewer 3265 The `Viewer` for saving. 3266 3267 See Also 3268 -------- 3269 DM, DMPlex, DM.view, DMPlex.topologyView, DMPlex.labelsView 3270 DMPlex.coordinatesLoad, Viewer, petsc.DMPlexCoordinatesView 3271 3272 """ 3273 CHKERR(DMPlexCoordinatesView(self.dm, viewer.vwr)) 3274 3275 def labelsView(self, Viewer viewer) -> None: 3276 """Save `DMPlex` labels into a file. 3277 3278 Collective. 3279 3280 Parameters 3281 ---------- 3282 viewer 3283 The `Viewer` for saving. 3284 3285 See Also 3286 -------- 3287 DM, DMPlex, DM.view, DMPlex.topologyView, DMPlex.coordinatesView 3288 DMPlex.labelsLoad, Viewer, petsc.DMPlexLabelsView 3289 3290 """ 3291 CHKERR(DMPlexLabelsView(self.dm, viewer.vwr)) 3292 3293 def sectionView(self, Viewer viewer, DM sectiondm) -> None: 3294 """Save a section associated with a `DMPlex`. 3295 3296 Collective. 3297 3298 Parameters 3299 ---------- 3300 viewer 3301 The `Viewer` for saving. 3302 sectiondm 3303 The `DM` that contains the section to be saved. 3304 3305 See Also 3306 -------- 3307 DM, DMPlex, DM.view, DMPlex.topologyView, DMPlex.coordinatesView 3308 DMPlex.labelsView, DMPlex.globalVectorView, DMPlex.localVectorView 3309 DMPlex.sectionLoad, Viewer, petsc.DMPlexSectionView 3310 3311 """ 3312 CHKERR(DMPlexSectionView(self.dm, viewer.vwr, sectiondm.dm)) 3313 3314 def globalVectorView(self, Viewer viewer, DM sectiondm, Vec vec) -> None: 3315 """Save a global vector. 3316 3317 Collective. 3318 3319 Parameters 3320 ---------- 3321 viewer 3322 The `Viewer` to save data with. 3323 sectiondm 3324 The `DM` containing the global section on which ``vec`` 3325 is defined; may be the same as this `DMPlex` object. 3326 vec 3327 The global vector to be saved. 3328 3329 See Also 3330 -------- 3331 DM, DMPlex, DMPlex.topologyView, DMPlex.sectionView 3332 DMPlex.localVectorView, DMPlex.globalVectorLoad 3333 DMPlex.localVectorLoad, petsc.DMPlexGlobalVectorView 3334 3335 """ 3336 CHKERR(DMPlexGlobalVectorView(self.dm, viewer.vwr, sectiondm.dm, vec.vec)) 3337 3338 def localVectorView(self, Viewer viewer, DM sectiondm, Vec vec) -> None: 3339 """Save a local vector. 3340 3341 Collective. 3342 3343 Parameters 3344 ---------- 3345 viewer 3346 The `Viewer` to save data with. 3347 sectiondm 3348 The `DM` that contains the local section on which ``vec`` is 3349 defined; may be the same as this `DMPlex` object. 3350 vec 3351 The local vector to be saved. 3352 3353 See Also 3354 -------- 3355 DM, DMPlex, DMPlex.topologyView, DMPlex.sectionView 3356 DMPlex.globalVectorView, DMPlex.globalVectorLoad 3357 DMPlex.localVectorLoad, petsc.DMPlexLocalVectorView 3358 3359 """ 3360 CHKERR(DMPlexLocalVectorView(self.dm, viewer.vwr, sectiondm.dm, vec.vec)) 3361 3362 # Load 3363 3364 def topologyLoad(self, Viewer viewer) -> SF: 3365 """Load a topology into this `DMPlex` object. 3366 3367 Collective. 3368 3369 Parameters 3370 ---------- 3371 viewer 3372 The `Viewer` for the saved topology 3373 3374 Returns 3375 ------- 3376 sfxc : SF 3377 The `SF` that pushes points in ``[0, N)`` to the associated points 3378 in the loaded `DMPlex`, where ``N`` is the global number of points. 3379 3380 See Also 3381 -------- 3382 DM, DMPlex, DM.load, DMPlex.coordinatesLoad, DMPlex.labelsLoad 3383 DM.view, SF, Viewer, petsc.DMPlexTopologyLoad 3384 3385 """ 3386 cdef SF sf = SF() 3387 CHKERR(DMPlexTopologyLoad(self.dm, viewer.vwr, &sf.sf)) 3388 return sf 3389 3390 def coordinatesLoad(self, Viewer viewer, SF sfxc) -> None: 3391 """Load coordinates into this `DMPlex` object. 3392 3393 Collective. 3394 3395 Parameters 3396 ---------- 3397 viewer 3398 The `Viewer` for the saved coordinates. 3399 sfxc 3400 The `SF` returned by `topologyLoad`. 3401 3402 See Also 3403 -------- 3404 DM, DMPlex, DM.load, DMPlex.topologyLoad, DMPlex.labelsLoad, DM.view 3405 SF, Viewer, petsc.DMPlexCoordinatesLoad 3406 3407 """ 3408 CHKERR(DMPlexCoordinatesLoad(self.dm, viewer.vwr, sfxc.sf)) 3409 3410 def labelsLoad(self, Viewer viewer, SF sfxc) -> None: 3411 """Load labels into this `DMPlex` object. 3412 3413 Collective. 3414 3415 Parameters 3416 ---------- 3417 viewer 3418 The `Viewer` for the saved labels. 3419 sfxc 3420 The `SF` returned by `topologyLoad`. 3421 3422 See Also 3423 -------- 3424 DM, DMPlex, DM.load, DMPlex.topologyLoad, DMPlex.coordinatesLoad 3425 DM.view, SF, Viewer, petsc.DMPlexLabelsLoad 3426 3427 """ 3428 CHKERR(DMPlexLabelsLoad(self.dm, viewer.vwr, sfxc.sf)) 3429 3430 def sectionLoad(self, Viewer viewer, DM sectiondm, SF sfxc) -> tuple[SF, SF]: 3431 """Load section into a `DM`. 3432 3433 Collective. 3434 3435 Parameters 3436 ---------- 3437 viewer 3438 The `Viewer` that represents the on-disk section (``sectionA``). 3439 sectiondm 3440 The `DM` into which the on-disk section (``sectionA``) is migrated. 3441 sfxc 3442 The `SF` returned by `topologyLoad`. 3443 3444 Returns 3445 ------- 3446 gsf : SF 3447 The `SF` that migrates any on-disk `Vec` data associated with 3448 ``sectionA`` into a global `Vec` associated with the 3449 ``sectiondm``'s global section (`None` if not needed). 3450 lsf : SF 3451 The `SF` that migrates any on-disk `Vec` data associated with 3452 ``sectionA`` into a local `Vec` associated with the ``sectiondm``'s 3453 local section (`None` if not needed). 3454 3455 See Also 3456 -------- 3457 DM, DMPlex, DM.load, DMPlex.topologyLoad, DMPlex.coordinatesLoad 3458 DMPlex.labelsLoad, DMPlex.globalVectorLoad, DMPlex.localVectorLoad 3459 DMPlex.sectionView, SF, Viewer, petsc.DMPlexSectionLoad 3460 3461 """ 3462 cdef SF gsf = SF() 3463 cdef SF lsf = SF() 3464 CHKERR(DMPlexSectionLoad(self.dm, viewer.vwr, sectiondm.dm, sfxc.sf, &gsf.sf, &lsf.sf)) 3465 return gsf, lsf 3466 3467 def globalVectorLoad(self, Viewer viewer, DM sectiondm, SF sf, Vec vec) -> None: 3468 """Load on-disk vector data into a global vector. 3469 3470 Collective. 3471 3472 Parameters 3473 ---------- 3474 viewer 3475 The `Viewer` that represents the on-disk vector data. 3476 sectiondm 3477 The `DM` that contains the global section on which vec is defined. 3478 sf 3479 The `SF` that migrates the on-disk vector data into vec. 3480 vec 3481 The global vector to set values of. 3482 3483 See Also 3484 -------- 3485 DM, DMPlex, DMPlex.topologyLoad, DMPlex.sectionLoad 3486 DMPlex.localVectorLoad, DMPlex.globalVectorView 3487 DMPlex.localVectorView, SF, Viewer, petsc.DMPlexGlobalVectorLoad 3488 3489 """ 3490 CHKERR(DMPlexGlobalVectorLoad(self.dm, viewer.vwr, sectiondm.dm, sf.sf, vec.vec)) 3491 3492 def localVectorLoad(self, Viewer viewer, DM sectiondm, SF sf, Vec vec) -> None: 3493 """Load on-disk vector data into a local vector. 3494 3495 Collective. 3496 3497 Parameters 3498 ---------- 3499 viewer 3500 The `Viewer` that represents the on-disk vector data. 3501 sectiondm 3502 The `DM` that contains the local section on which vec is defined. 3503 sf 3504 The `SF` that migrates the on-disk vector data into vec. 3505 vec 3506 The local vector to set values of. 3507 3508 See Also 3509 -------- 3510 DM, DMPlex, DMPlex.topologyLoad, DMPlex.sectionLoad 3511 DMPlex.globalVectorLoad, DMPlex.globalVectorView 3512 DMPlex.localVectorView, SF, Viewer, petsc.DMPlexLocalVectorLoad 3513 3514 """ 3515 CHKERR(DMPlexLocalVectorLoad(self.dm, viewer.vwr, sectiondm.dm, sf.sf, vec.vec)) 3516 3517 def createNaturalVec(self) -> Vec: 3518 """Return a natural vector. 3519 3520 Collective. 3521 3522 See Also 3523 -------- 3524 petsc.DMPlexCreateNaturalVector 3525 3526 """ 3527 cdef Vec nv = Vec() 3528 CHKERR(DMPlexCreateNaturalVector(self.dm, &nv.vec)) 3529 return nv 3530 3531 def naturalToGlobalBegin(self, Vec nv, Vec gv) -> None: 3532 """Rearrange a `Vec` in the natural order to the Global order. 3533 3534 Collective. 3535 3536 See Also 3537 -------- 3538 petsc.DMPlexNaturalToGlobalBegin 3539 3540 """ 3541 CHKERR(DMPlexNaturalToGlobalBegin(self.dm, nv.vec, gv.vec)) 3542 3543 def naturalToGlobalEnd(self, Vec nv, Vec gv) -> None: 3544 """Rearrange a `Vec` in the natural order to the Global order. 3545 3546 Collective. 3547 3548 See Also 3549 -------- 3550 petsc.DMPlexNaturalToGlobalEnd 3551 3552 """ 3553 CHKERR(DMPlexNaturalToGlobalEnd(self.dm, nv.vec, gv.vec)) 3554 3555 def globalToNaturalBegin(self, Vec gv, Vec nv) -> None: 3556 """Rearrange a `Vec` in the Global order to the natural order. 3557 3558 Collective. 3559 3560 See Also 3561 -------- 3562 petsc.DMPlexGlobalToNaturalBegin 3563 3564 """ 3565 CHKERR(DMPlexGlobalToNaturalBegin(self.dm, gv.vec, nv.vec)) 3566 3567 def globalToNaturalEnd(self, Vec gv, Vec nv) -> None: 3568 """Rearrange a `Vec` in the Global order to the natural order. 3569 3570 Collective. 3571 3572 See Also 3573 -------- 3574 petsc.DMPlexGlobalToNaturalEnd 3575 3576 """ 3577 CHKERR(DMPlexGlobalToNaturalEnd(self.dm, gv.vec, nv.vec)) 3578 3579 def setMigrationSF(self, SF sf) -> None: 3580 """Set the `SF` for migrating from a parent `DM` into this `DM`. 3581 3582 Not collective. 3583 3584 See Also 3585 -------- 3586 petsc.DMPlexSetMigrationSF 3587 3588 """ 3589 CHKERR(DMPlexSetMigrationSF(self.dm, sf.sf)) 3590 3591 def getMigrationSF(self) -> SF: 3592 """Get the `SF` for migrating from a parent `DM` into this `DM`. 3593 3594 Not collective. 3595 3596 See Also 3597 -------- 3598 petsc.DMPlexGetMigrationSF 3599 3600 """ 3601 cdef SF sf = SF() 3602 CHKERR(DMPlexGetMigrationSF(self.dm, &sf.sf)) 3603 CHKERR(PetscINCREF(sf.obj)) 3604 return sf 3605 3606 def createGlobalToNaturalSF(self, Section section, SF sfMigration) -> SF: 3607 """Create the `SF` for mapping Global `Vec` to the Natural `Vec`. 3608 3609 Collective. 3610 3611 See Also 3612 -------- 3613 petsc.DMPlexCreateGlobalToNaturalSF 3614 3615 """ 3616 cdef SF sf = SF() 3617 CHKERR(DMPlexCreateGlobalToNaturalSF(self.dm, section.sec, sfMigration.sf, &sf.sf)) 3618 return sf 3619 3620 def migrateGlobalToNaturalSF(self, DM dmOld, SF sfNaturalOld, SF sfMigration) -> SF: 3621 """Create the `SF` for mapping Global `Vec` to the Natural `Vec` in the new `DM`. 3622 3623 Collective. 3624 3625 See Also 3626 -------- 3627 petsc.DMPlexMigrateGlobalToNaturalSF 3628 3629 """ 3630 cdef SF sf = SF() 3631 CHKERR(DMPlexMigrateGlobalToNaturalSF(dmOld.dm, self.dm, sfNaturalOld.sf, sfMigration.sf, &sf.sf)) 3632 return sf 3633 3634# -------------------------------------------------------------------- 3635 3636 3637class DMPlexTransformType(object): 3638 """Transformation types.""" 3639 REFINEREGULAR = S_(DMPLEXREFINEREGULAR) 3640 REFINEALFELD = S_(DMPLEXREFINEALFELD) 3641 REFINEPOWELLSABIN = S_(DMPLEXREFINEPOWELLSABIN) 3642 REFINEBOUNDARYLAYER = S_(DMPLEXREFINEBOUNDARYLAYER) 3643 REFINESBR = S_(DMPLEXREFINESBR) 3644 REFINETOBOX = S_(DMPLEXREFINETOBOX) 3645 REFINETOSIMPLEX = S_(DMPLEXREFINETOSIMPLEX) 3646 REFINE1D = S_(DMPLEXREFINE1D) 3647 EXTRUDE = S_(DMPLEXEXTRUDETYPE) 3648 TRANSFORMFILTER = S_(DMPLEXTRANSFORMFILTER) 3649 3650 3651cdef class DMPlexTransform(Object): 3652 """Mesh transformations.""" 3653 3654 def __cinit__(self): 3655 self.obj = <PetscObject*> &self.tr 3656 self.tr = NULL 3657 3658 def apply(self, DM dm) -> DM: 3659 """Apply a mesh transformation. 3660 3661 Collective. 3662 3663 """ 3664 # FIXME petsc.DMPlexTransformApply 3665 cdef DMPlex newdm = DMPlex() 3666 CHKERR(DMPlexTransformApply(self.tr, dm.dm, &newdm.dm)) 3667 return newdm 3668 3669 def create(self, comm: Comm | None = None) -> Self: 3670 """Create a mesh transformation. 3671 3672 Collective. 3673 3674 See Also 3675 -------- 3676 petsc.DMPlexTransformCreate 3677 3678 """ 3679 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 3680 cdef PetscDMPlexTransform newtr = NULL 3681 CHKERR(DMPlexTransformCreate(ccomm, &newtr)) 3682 CHKERR(PetscCLEAR(self.obj)) 3683 self.tr = newtr 3684 return self 3685 3686 def destroy(self) -> Self: 3687 """Destroy a mesh transformation. 3688 3689 Collective. 3690 3691 See Also 3692 -------- 3693 petsc.DMPlexTransformDestroy 3694 3695 """ 3696 CHKERR(DMPlexTransformDestroy(&self.tr)) 3697 return self 3698 3699 def getType(self) -> str: 3700 """Return the transformation type name. 3701 3702 Not collective. 3703 3704 See Also 3705 -------- 3706 petsc.DMPlexTransformGetType 3707 3708 """ 3709 cdef PetscDMPlexTransformType cval = NULL 3710 CHKERR(DMPlexTransformGetType(self.tr, &cval)) 3711 return bytes2str(cval) 3712 3713 def setUp(self) -> Self: 3714 """Setup a mesh transformation. 3715 3716 Collective. 3717 3718 """ 3719 # FIXME petsc.DMPlexTransformSetUp 3720 CHKERR(DMPlexTransformSetUp(self.tr)) 3721 return self 3722 3723 def setType(self, tr_type : DMPlexTransformType | str) -> None: 3724 """Set the transformation type. 3725 3726 Collective. 3727 3728 See Also 3729 -------- 3730 petsc.DMPlexTransformSetType 3731 3732 """ 3733 cdef PetscDMPlexTransformType cval = NULL 3734 tr_type = str2bytes(tr_type, &cval) 3735 CHKERR(DMPlexTransformSetType(self.tr, cval)) 3736 3737 def setDM(self, DM dm) -> None: 3738 """Set the `DM` for the transformation. 3739 3740 Logically collective. 3741 3742 """ 3743 # FIXME petsc.DMPlexTransformSetDM 3744 CHKERR(DMPlexTransformSetDM(self.tr, dm.dm)) 3745 3746 def setFromOptions(self) -> None: 3747 """Configure the transformation from the options database. 3748 3749 Collective. 3750 3751 See Also 3752 -------- 3753 petsc_options, petsc.DMPlexTransformSetFromOptions 3754 3755 """ 3756 CHKERR(DMPlexTransformSetFromOptions(self.tr)) 3757 3758 def view(self, Viewer viewer=None) -> None: 3759 """View the mesh transformation. 3760 3761 Collective. 3762 3763 Parameters 3764 ---------- 3765 viewer 3766 A `Viewer` instance or `None` for the default viewer. 3767 3768 See Also 3769 -------- 3770 Viewer, petsc.DMPlexTransformView 3771 3772 """ 3773 cdef PetscViewer vwr = NULL 3774 if viewer is not None: vwr = viewer.vwr 3775 CHKERR(DMPlexTransformView(self.tr, vwr)) 3776