1# -------------------------------------------------------------------- 2 3class PCType(object): 4 """The preconditioner method.""" 5 NONE = S_(PCNONE) 6 JACOBI = S_(PCJACOBI) 7 SOR = S_(PCSOR) 8 LU = S_(PCLU) 9 QR = S_(PCQR) 10 SHELL = S_(PCSHELL) 11 BJACOBI = S_(PCBJACOBI) 12 VPBJACOBI = S_(PCVPBJACOBI) 13 MG = S_(PCMG) 14 EISENSTAT = S_(PCEISENSTAT) 15 ILU = S_(PCILU) 16 ICC = S_(PCICC) 17 ASM = S_(PCASM) 18 GASM = S_(PCGASM) 19 KSP = S_(PCKSP) 20 COMPOSITE = S_(PCCOMPOSITE) 21 REDUNDANT = S_(PCREDUNDANT) 22 SPAI = S_(PCSPAI) 23 NN = S_(PCNN) 24 CHOLESKY = S_(PCCHOLESKY) 25 PBJACOBI = S_(PCPBJACOBI) 26 MAT = S_(PCMAT) 27 HYPRE = S_(PCHYPRE) 28 PARMS = S_(PCPARMS) 29 FIELDSPLIT = S_(PCFIELDSPLIT) 30 TFS = S_(PCTFS) 31 ML = S_(PCML) 32 GALERKIN = S_(PCGALERKIN) 33 EXOTIC = S_(PCEXOTIC) 34 CP = S_(PCCP) 35 LSC = S_(PCLSC) 36 PYTHON = S_(PCPYTHON) 37 PFMG = S_(PCPFMG) 38 SYSPFMG = S_(PCSYSPFMG) 39 REDISTRIBUTE = S_(PCREDISTRIBUTE) 40 SVD = S_(PCSVD) 41 GAMG = S_(PCGAMG) 42 CHOWILUVIENNACL = S_(PCCHOWILUVIENNACL) 43 ROWSCALINGVIENNACL = S_(PCROWSCALINGVIENNACL) 44 SAVIENNACL = S_(PCSAVIENNACL) 45 BDDC = S_(PCBDDC) 46 KACZMARZ = S_(PCKACZMARZ) 47 TELESCOPE = S_(PCTELESCOPE) 48 PATCH = S_(PCPATCH) 49 LMVM = S_(PCLMVM) 50 HMG = S_(PCHMG) 51 DEFLATION = S_(PCDEFLATION) 52 HPDDM = S_(PCHPDDM) 53 H2OPUS = S_(PCH2OPUS) 54 55 56class PCSide(object): 57 """The manner in which the preconditioner is applied.""" 58 # native 59 LEFT = PC_LEFT 60 RIGHT = PC_RIGHT 61 SYMMETRIC = PC_SYMMETRIC 62 # aliases 63 L = LEFT 64 R = RIGHT 65 S = SYMMETRIC 66 67 68class PCASMType(object): 69 """The *ASM* subtype.""" 70 NONE = PC_ASM_NONE 71 BASIC = PC_ASM_BASIC 72 RESTRICT = PC_ASM_RESTRICT 73 INTERPOLATE = PC_ASM_INTERPOLATE 74 75 76class PCGASMType(object): 77 """The *GASM* subtype.""" 78 NONE = PC_GASM_NONE 79 BASIC = PC_GASM_BASIC 80 RESTRICT = PC_GASM_RESTRICT 81 INTERPOLATE = PC_GASM_INTERPOLATE 82 83 84class PCMGType(object): 85 """The *MG* subtype.""" 86 MULTIPLICATIVE = PC_MG_MULTIPLICATIVE 87 ADDITIVE = PC_MG_ADDITIVE 88 FULL = PC_MG_FULL 89 KASKADE = PC_MG_KASKADE 90 91 92class PCMGCycleType(object): 93 """The *MG* cycle type.""" 94 V = PC_MG_CYCLE_V 95 W = PC_MG_CYCLE_W 96 97 98class PCGAMGType(object): 99 """The *GAMG* subtype.""" 100 AGG = S_(PCGAMGAGG) 101 GEO = S_(PCGAMGGEO) 102 CLASSICAL = S_(PCGAMGCLASSICAL) 103 104 105class PCCompositeType(object): 106 """The composite type.""" 107 ADDITIVE = PC_COMPOSITE_ADDITIVE 108 MULTIPLICATIVE = PC_COMPOSITE_MULTIPLICATIVE 109 SYMMETRIC_MULTIPLICATIVE = PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE 110 SPECIAL = PC_COMPOSITE_SPECIAL 111 SCHUR = PC_COMPOSITE_SCHUR 112 113 114class PCFieldSplitSchurPreType(object): 115 """The field split Schur subtype.""" 116 SELF = PC_FIELDSPLIT_SCHUR_PRE_SELF 117 SELFP = PC_FIELDSPLIT_SCHUR_PRE_SELFP 118 A11 = PC_FIELDSPLIT_SCHUR_PRE_A11 119 USER = PC_FIELDSPLIT_SCHUR_PRE_USER 120 FULL = PC_FIELDSPLIT_SCHUR_PRE_FULL 121 122 123class PCFieldSplitSchurFactType(object): 124 """The field split Schur factorization type.""" 125 DIAG = PC_FIELDSPLIT_SCHUR_FACT_DIAG 126 LOWER = PC_FIELDSPLIT_SCHUR_FACT_LOWER 127 UPPER = PC_FIELDSPLIT_SCHUR_FACT_UPPER 128 FULL = PC_FIELDSPLIT_SCHUR_FACT_FULL 129 130 131class PCPatchConstructType(object): 132 """The patch construction type.""" 133 STAR = PC_PATCH_STAR 134 VANKA = PC_PATCH_VANKA 135 PARDECOMP = PC_PATCH_PARDECOMP 136 USER = PC_PATCH_USER 137 PYTHON = PC_PATCH_PYTHON 138 139 140class PCHPDDMCoarseCorrectionType(object): 141 """The *HPDDM* coarse correction type.""" 142 DEFLATED = PC_HPDDM_COARSE_CORRECTION_DEFLATED 143 ADDITIVE = PC_HPDDM_COARSE_CORRECTION_ADDITIVE 144 BALANCED = PC_HPDDM_COARSE_CORRECTION_BALANCED 145 NONE = PC_HPDDM_COARSE_CORRECTION_NONE 146 147 148class PCDeflationSpaceType(object): 149 """The deflation space subtype.""" 150 HAAR = PC_DEFLATION_SPACE_HAAR 151 DB2 = PC_DEFLATION_SPACE_DB2 152 DB4 = PC_DEFLATION_SPACE_DB4 153 DB8 = PC_DEFLATION_SPACE_DB8 154 DB16 = PC_DEFLATION_SPACE_DB16 155 BIORTH22 = PC_DEFLATION_SPACE_BIORTH22 156 MEYER = PC_DEFLATION_SPACE_MEYER 157 AGGREGATION = PC_DEFLATION_SPACE_AGGREGATION 158 USER = PC_DEFLATION_SPACE_USER 159 160 161class PCFailedReason(object): 162 """The reason the preconditioner has failed.""" 163 SETUP_ERROR = PC_SETUP_ERROR 164 NOERROR = PC_NOERROR 165 FACTOR_STRUCT_ZEROPIVOT = PC_FACTOR_STRUCT_ZEROPIVOT 166 FACTOR_NUMERIC_ZEROPIVOT = PC_FACTOR_NUMERIC_ZEROPIVOT 167 FACTOR_OUTMEMORY = PC_FACTOR_OUTMEMORY 168 FACTOR_OTHER = PC_FACTOR_OTHER 169 SUBPC_ERROR = PC_SUBPC_ERROR 170 171# -------------------------------------------------------------------- 172 173 174cdef class PC(Object): 175 """Preconditioners. 176 177 `PC` is described in the `PETSc manual <petsc:sec_ksppc>`. 178 Calling the `PC` with a vector as an argument will `apply` the 179 preconditioner as shown in the example below. 180 181 Examples 182 -------- 183 >>> from petsc4py import PETSc 184 >>> v = PETSc.Vec().createWithArray([1, 2]) 185 >>> m = PETSc.Mat().createDense(2, array=[[1, 0], [0, 1]]) 186 >>> pc = PETSc.PC().create() 187 >>> pc.setOperators(m) 188 >>> u = pc(v) # u is created internally 189 >>> pc.apply(v, u) # u can also be passed as second argument 190 191 See Also 192 -------- 193 petsc.PC 194 195 """ 196 197 Type = PCType 198 Side = PCSide 199 200 ASMType = PCASMType 201 GASMType = PCGASMType 202 MGType = PCMGType 203 MGCycleType = PCMGCycleType 204 GAMGType = PCGAMGType 205 CompositeType = PCCompositeType 206 FieldSplitSchurFactType = PCFieldSplitSchurFactType 207 FieldSplitSchurPreType = PCFieldSplitSchurPreType 208 PatchConstructType = PCPatchConstructType 209 HPDDMCoarseCorrectionType = PCHPDDMCoarseCorrectionType 210 DeflationSpaceType = PCDeflationSpaceType 211 FailedReason = PCFailedReason 212 # Backward compatibility 213 SchurFactType = PCFieldSplitSchurFactType 214 SchurPreType = PCFieldSplitSchurPreType 215 216 # --- xxx --- 217 218 def __cinit__(self): 219 self.obj = <PetscObject*> &self.pc 220 self.pc = NULL 221 222 def __call__(self, x, y=None): 223 if y is None: # XXX do this better 224 y = self.getOperators()[0].createVecLeft() 225 self.apply(x, y) 226 return y 227 228 # --- xxx --- 229 230 def view(self, Viewer viewer=None) -> None: 231 """View the `PC` object. 232 233 Collective. 234 235 Parameters 236 ---------- 237 viewer 238 The visualization context. 239 240 See Also 241 -------- 242 petsc.PCView 243 244 """ 245 cdef PetscViewer vwr = NULL 246 if viewer is not None: vwr = viewer.vwr 247 CHKERR(PCView(self.pc, vwr)) 248 249 def destroy(self) -> Self: 250 """Destroy the `PC` that was created with `create`. 251 252 Collective. 253 254 See Also 255 -------- 256 petsc.PCDestroy 257 258 """ 259 CHKERR(PCDestroy(&self.pc)) 260 self.pc = NULL 261 return self 262 263 def create(self, comm: Comm | None = None) -> Self: 264 """Create an empty `PC`. 265 266 Collective. 267 268 The default preconditioner for sparse matrices is `ILU` or 269 `ICC` with 0 fill on one process and block Jacobi (`BJACOBI`) with `ILU` 270 or `ICC` in parallel. For dense matrices it is always `None`. 271 272 Parameters 273 ---------- 274 comm 275 MPI communicator, defaults to `Sys.getDefaultComm`. 276 277 See Also 278 -------- 279 destroy, petsc.PCCreate 280 281 """ 282 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 283 cdef PetscPC newpc = NULL 284 CHKERR(PCCreate(ccomm, &newpc)) 285 CHKERR(PetscCLEAR(self.obj)); self.pc = newpc 286 return self 287 288 def setType(self, pc_type: Type | str) -> None: 289 """Set the preconditioner type. 290 291 Collective. 292 293 Parameters 294 ---------- 295 pc_type 296 The preconditioner type. 297 298 See Also 299 -------- 300 petsc_options, getType, petsc.TSSetType 301 302 """ 303 cdef PetscPCType cval = NULL 304 pc_type = str2bytes(pc_type, &cval) 305 CHKERR(PCSetType(self.pc, cval)) 306 307 def getType(self) -> str: 308 """Return the preconditioner type. 309 310 Not collective. 311 312 See Also 313 -------- 314 setType, petsc.PCGetType 315 316 """ 317 cdef PetscPCType cval = NULL 318 CHKERR(PCGetType(self.pc, &cval)) 319 return bytes2str(cval) 320 321 def setOptionsPrefix(self, prefix: str | None) -> None: 322 """Set the prefix used for all the `PC` options. 323 324 Logically collective. 325 326 Parameters 327 ---------- 328 prefix 329 The prefix to prepend to all option names. 330 331 See Also 332 -------- 333 petsc_options, petsc.PCSetOptionsPrefix 334 335 """ 336 cdef const char *cval = NULL 337 prefix = str2bytes(prefix, &cval) 338 CHKERR(PCSetOptionsPrefix(self.pc, cval)) 339 340 def getOptionsPrefix(self) -> str: 341 """Return the prefix used for all the `PC` options. 342 343 Not collective. 344 345 See Also 346 -------- 347 petsc_options, petsc.PCGetOptionsPrefix 348 349 """ 350 cdef const char *cval = NULL 351 CHKERR(PCGetOptionsPrefix(self.pc, &cval)) 352 return bytes2str(cval) 353 354 def appendOptionsPrefix(self, prefix: str | None) -> None: 355 """Append to the prefix used for all the `PC` options. 356 357 Logically collective. 358 359 Parameters 360 ---------- 361 prefix 362 The prefix to append to the current prefix. 363 364 See Also 365 -------- 366 petsc_options, petsc.PCAppendOptionsPrefix 367 368 """ 369 cdef const char *cval = NULL 370 prefix = str2bytes(prefix, &cval) 371 CHKERR(PCAppendOptionsPrefix(self.pc, cval)) 372 373 def setFromOptions(self) -> None: 374 """Set various `PC` parameters from user options. 375 376 Collective. 377 378 See Also 379 -------- 380 petsc_options, petsc.PCSetFromOptions 381 382 """ 383 CHKERR(PCSetFromOptions(self.pc)) 384 385 def setOperators(self, Mat A=None, Mat P=None) -> None: 386 """Set the matrices associated with the linear system. 387 388 Logically collective. 389 390 Passing `None` for ``A`` or ``P`` removes the 391 matrix that is currently used. PETSc does not reset the matrix entries 392 of either ``A`` or ``P`` to zero after a linear solve; the user is 393 completely responsible for matrix assembly. See `Mat.zeroEntries` to 394 zero all elements of a matrix. 395 396 Parameters 397 ---------- 398 A 399 the matrix which defines the linear system 400 P 401 the matrix to be used in constructing the preconditioner, usually 402 the same as ``A`` 403 404 Notes 405 ----- 406 Using this directly is rarely needed, the preferred, and equivalent, usage 407 is to call `KSP.setOperators`. 408 409 See Also 410 -------- 411 petsc.PCSetOperators 412 413 """ 414 cdef PetscMat amat=NULL 415 if A is not None: amat = A.mat 416 cdef PetscMat pmat=amat 417 if P is not None: pmat = P.mat 418 CHKERR(PCSetOperators(self.pc, amat, pmat)) 419 420 def getOperators(self) -> tuple[Mat, Mat]: 421 """Return the matrices associated with a linear system. 422 423 Not collective. 424 425 See Also 426 -------- 427 setOperators, petsc.PCGetOperators 428 429 """ 430 cdef Mat A = Mat(), P = Mat() 431 CHKERR(PCGetOperators(self.pc, &A.mat, &P.mat)) 432 CHKERR(PetscINCREF(A.obj)) 433 CHKERR(PetscINCREF(P.obj)) 434 return (A, P) 435 436 def setUseAmat(self, flag: bool) -> None: 437 """Set to indicate to apply `PC` to ``A`` and not ``P``. 438 439 Logically collective. 440 441 Sets a flag to indicate that when the 442 preconditioner needs to apply (part of) the operator during the 443 preconditioning process, it applies to ``A`` provided to 444 `TS.setRHSJacobian`, `TS.setIJacobian`, `SNES.setJacobian`, 445 `KSP.setOperators` or `PC.setOperators` not the ``P``. 446 447 Parameters 448 ---------- 449 flag 450 Set True to use ``A`` and False to use ``P``. 451 452 See Also 453 -------- 454 setOperators, petsc.PCSetUseAmat 455 456 """ 457 cdef PetscBool cflag = PETSC_FALSE 458 if flag: 459 cflag = PETSC_TRUE 460 CHKERR(PCSetUseAmat(self.pc, cflag)) 461 462 def getUseAmat(self) -> bool: 463 """Return the flag to indicate if `PC` is applied to ``A`` or ``P``. 464 465 Logically collective. 466 467 Returns 468 ------- 469 flag : bool 470 True if ``A`` is used and False if ``P``. 471 472 See Also 473 -------- 474 setUseAmat, petsc.PCGetUseAmat 475 476 """ 477 cdef PetscBool cflag = PETSC_FALSE 478 CHKERR(PCGetUseAmat(self.pc, &cflag)) 479 return toBool(cflag) 480 481 def setReusePreconditioner(self, flag: bool) -> None: 482 """Set to indicate the preconditioner is to be reused. 483 484 Logically collective. 485 486 Normally if the ``A`` matrix inside a `PC` 487 changes, the `PC` automatically updates itself using information from 488 the changed matrix. Enable this option prevents this. 489 490 Parameters 491 ---------- 492 flag 493 Set to `True` to use the reuse the current preconditioner and 494 `False` to recompute on changes to the matrix. 495 496 See Also 497 -------- 498 setOperators, petsc.PCSetReusePreconditioner 499 500 """ 501 cdef PetscBool cflag = PETSC_FALSE 502 if flag: 503 cflag = PETSC_TRUE 504 CHKERR(PCSetReusePreconditioner(self.pc, cflag)) 505 506 def setFailedReason(self, reason: FailedReason | str) -> None: 507 """Set the reason the `PC` terminated. 508 509 Logically collective. 510 511 Parameters 512 ---------- 513 reason 514 the reason the `PC` terminated 515 516 See Also 517 -------- 518 petsc.PCSetFailedReason 519 520 """ 521 cdef PetscPCFailedReason val = reason 522 CHKERR(PCSetFailedReason(self.pc, val)) 523 524 def getFailedReason(self) -> FailedReason: 525 """Return the reason the `PC` terminated. 526 527 Not collective. 528 529 After a call to KSPCheckDot() or 530 KSPCheckNorm() inside a KSPSolve(), or after 531 a call to PCReduceFailedReason() this is the maximum 532 reason over all ranks in the 533 PC communicator and hence logically collective. 534 Otherwise it is the local value. 535 536 See Also 537 -------- 538 petsc.PCGetFailedReason 539 540 """ 541 cdef PetscPCFailedReason reason = PC_NOERROR 542 CHKERR(PCGetFailedReason(self.pc, &reason)) 543 return reason 544 545 def setUp(self) -> None: 546 """Set up the internal data structures for the `PC`. 547 548 Collective. 549 550 See Also 551 -------- 552 petsc.PCSetUp 553 554 """ 555 CHKERR(PCSetUp(self.pc)) 556 557 def reset(self) -> None: 558 """Reset the `PC`, removing any allocated vectors and matrices. 559 560 Collective. 561 562 See Also 563 -------- 564 petsc.PCReset 565 566 """ 567 CHKERR(PCReset(self.pc)) 568 569 def setUpOnBlocks(self) -> None: 570 """Set up the `PC` for each block. 571 572 Collective. 573 574 For nested preconditioners such as `BJACOBI`, `setUp` is not 575 called on each sub-`KSP` when `setUp` is called on the outer `PC`. This 576 routine ensures it is called. 577 578 See Also 579 -------- 580 setUp, petsc.PCSetUpOnBlocks 581 582 """ 583 CHKERR(PCSetUpOnBlocks(self.pc)) 584 585 def apply(self, Vec x, Vec y) -> None: 586 """Apply the `PC` to a vector. 587 588 Collective. 589 590 Parameters 591 ---------- 592 x 593 The input vector. 594 y 595 The output vector, cannot be the same as ``x``. 596 597 See Also 598 -------- 599 petsc.PCApply 600 601 """ 602 CHKERR(PCApply(self.pc, x.vec, y.vec)) 603 604 def matApply(self, Mat x, Mat y) -> None: 605 """Apply the `PC` to many vectors stored as `Mat.Type.DENSE`. 606 607 Collective. 608 609 Parameters 610 ---------- 611 x 612 The input matrix. 613 y 614 The output matrix, cannot be the same as ``x``. 615 616 See Also 617 -------- 618 petsc.PCMatApply, petsc.PCApply 619 620 """ 621 CHKERR(PCMatApply(self.pc, x.mat, y.mat)) 622 623 def applyTranspose(self, Vec x, Vec y) -> None: 624 """Apply the transpose of the `PC` to a vector. 625 626 Collective. 627 628 For complex numbers this applies the non-Hermitian 629 transpose. 630 631 Parameters 632 ---------- 633 x 634 The input vector. 635 y 636 The output vector, cannot be the same as ``x``. 637 638 See Also 639 -------- 640 petsc.PCApply 641 642 """ 643 CHKERR(PCApplyTranspose(self.pc, x.vec, y.vec)) 644 645 def matApplyTranspose(self, Mat x, Mat y) -> None: 646 """Apply the transpose of the `PC` to many vectors stored as `Mat.Type.DENSE`. 647 648 Collective. 649 650 Parameters 651 ---------- 652 x 653 The input matrix. 654 y 655 The output matrix, cannot be the same as ``x``. 656 657 See Also 658 -------- 659 petsc.PCMatApply, petsc.PCMatApplyTranspose 660 661 """ 662 CHKERR(PCMatApplyTranspose(self.pc, x.mat, y.mat)) 663 664 def applySymmetricLeft(self, Vec x, Vec y) -> None: 665 """Apply the left part of a symmetric `PC` to a vector. 666 667 Collective. 668 669 Parameters 670 ---------- 671 x 672 The input vector. 673 y 674 The output vector, cannot be the same as ``x``. 675 676 See Also 677 -------- 678 petsc.PCApplySymmetricLeft 679 680 """ 681 CHKERR(PCApplySymmetricLeft(self.pc, x.vec, y.vec)) 682 683 def applySymmetricRight(self, Vec x, Vec y) -> None: 684 """Apply the right part of a symmetric `PC` to a vector. 685 686 Collective. 687 688 Parameters 689 ---------- 690 x 691 The input vector. 692 y 693 The output vector, cannot be the same as ``x``. 694 695 See Also 696 -------- 697 petsc.PCApplySymmetricRight 698 699 """ 700 CHKERR(PCApplySymmetricRight(self.pc, x.vec, y.vec)) 701 702 # --- discretization space --- 703 704 def getDM(self) -> DM: 705 """Return the `DM` associated with the `PC`. 706 707 Not collective. 708 709 See Also 710 -------- 711 petsc.PCGetDM 712 713 """ 714 cdef PetscDM newdm = NULL 715 CHKERR(PCGetDM(self.pc, &newdm)) 716 cdef DM dm = subtype_DM(newdm)() 717 dm.dm = newdm 718 CHKERR(PetscINCREF(dm.obj)) 719 return dm 720 721 def setDM(self, DM dm) -> None: 722 """Set the `DM` that may be used by some preconditioners. 723 724 Logically collective. 725 726 Parameters 727 ---------- 728 dm 729 The `DM` object. 730 731 See Also 732 -------- 733 petsc.PCSetDM 734 735 """ 736 CHKERR(PCSetDM(self.pc, dm.dm)) 737 738 def setCoordinates(self, coordinates: Sequence[Sequence[float]]) -> None: 739 """Set the coordinates for the nodes on the local process. 740 741 Collective. 742 743 Parameters 744 ---------- 745 coordinates 746 The two dimensional coordinate array. 747 748 See Also 749 -------- 750 petsc.PCSetCoordinates 751 752 """ 753 cdef ndarray xyz = iarray(coordinates, NPY_PETSC_REAL) 754 if PyArray_ISFORTRAN(xyz): xyz = PyArray_Copy(xyz) 755 if PyArray_NDIM(xyz) != 2: raise ValueError( 756 ("coordinates must have two dimensions: " 757 "coordinates.ndim=%d") % (PyArray_NDIM(xyz))) 758 cdef PetscInt nvtx = <PetscInt> PyArray_DIM(xyz, 0) 759 cdef PetscInt ndim = <PetscInt> PyArray_DIM(xyz, 1) 760 cdef PetscReal *coords = <PetscReal*> PyArray_DATA(xyz) 761 CHKERR(PCSetCoordinates(self.pc, ndim, nvtx, coords)) 762 763 # --- Python --- 764 765 def createPython(self, context: Any = None, comm: Comm | None = None) -> Self: 766 """Create a preconditioner of Python type. 767 768 Collective. 769 770 Parameters 771 ---------- 772 context 773 An instance of the Python class implementing the required methods. 774 comm 775 MPI communicator, defaults to `Sys.getDefaultComm`. 776 777 See Also 778 -------- 779 petsc_python_pc, setType, setPythonContext, PC.Type.PYTHON 780 781 """ 782 cdef MPI_Comm ccomm = def_Comm(comm, PETSC_COMM_DEFAULT) 783 cdef PetscPC newpc = NULL 784 CHKERR(PCCreate(ccomm, &newpc)) 785 CHKERR(PetscCLEAR(self.obj)); self.pc = newpc 786 CHKERR(PCSetType(self.pc, PCPYTHON)) 787 CHKERR(PCPythonSetContext(self.pc, <void*>context)) 788 return self 789 790 def setPythonContext(self, context: Any) -> None: 791 """Set the instance of the class implementing the required Python methods. 792 793 Not collective. 794 795 See Also 796 -------- 797 petsc_python_pc, getPythonContext 798 799 """ 800 CHKERR(PCPythonSetContext(self.pc, <void*>context)) 801 802 def getPythonContext(self) -> Any: 803 """Return the instance of the class implementing the required Python methods. 804 805 Not collective. 806 807 See Also 808 -------- 809 petsc_python_pc, setPythonContext 810 811 """ 812 cdef void *context = NULL 813 CHKERR(PCPythonGetContext(self.pc, &context)) 814 if context == NULL: return None 815 else: return <object> context 816 817 def setPythonType(self, py_type: str) -> None: 818 """Set the fully qualified Python name of the class to be used. 819 820 Collective. 821 822 See Also 823 -------- 824 petsc_python_pc, setPythonContext, getPythonType, petsc.PCPythonSetType 825 826 """ 827 cdef const char *cval = NULL 828 py_type = str2bytes(py_type, &cval) 829 CHKERR(PCPythonSetType(self.pc, cval)) 830 831 def getPythonType(self) -> str: 832 """Return the fully qualified Python name of the class used by the preconditioner. 833 834 Not collective. 835 836 See Also 837 -------- 838 petsc_python_pc, setPythonContext, setPythonType, petsc.PCPythonGetType 839 840 """ 841 cdef const char *cval = NULL 842 CHKERR(PCPythonGetType(self.pc, &cval)) 843 return bytes2str(cval) 844 845 # --- Block Jacobi --- 846 847 def getBJacobiSubKSP(self) -> list[KSP]: 848 """Return the local `KSP` object for all blocks on this process. 849 850 Not collective. 851 852 See Also 853 -------- 854 petsc.PCBJacobiGetSubKSP 855 856 """ 857 cdef PetscInt n = 0 858 cdef PetscKSP *p = NULL 859 CHKERR(PCBJacobiGetSubKSP(self.pc, &n, NULL, &p)) 860 return [ref_KSP(p[i]) for i from 0 <= i <n] 861 862 # --- ASM --- 863 864 def setASMType(self, asmtype: ASMType) -> None: 865 """Set the type of restriction and interpolation. 866 867 Logically collective. 868 869 Parameters 870 ---------- 871 asmtype 872 The type of ASM you wish to use. 873 874 See Also 875 -------- 876 petsc.PCASMSetType 877 878 """ 879 cdef PetscPCASMType cval = asmtype 880 CHKERR(PCASMSetType(self.pc, cval)) 881 882 def setASMOverlap(self, overlap: int) -> None: 883 """Set the overlap between a pair of subdomains. 884 885 Logically collective. 886 887 Parameters 888 ---------- 889 overlap 890 The amount of overlap between subdomains. 891 892 See Also 893 -------- 894 petsc.PCASMSetOverlap 895 896 """ 897 cdef PetscInt ival = asInt(overlap) 898 CHKERR(PCASMSetOverlap(self.pc, ival)) 899 900 def setASMLocalSubdomains( 901 self, 902 nsd: int, 903 is_sub: Sequence[IS] | None = None, 904 is_local: Sequence[IS] | None = None) -> None: 905 """Set the local subdomains. 906 907 Collective. 908 909 Parameters 910 ---------- 911 nsd 912 The number of subdomains for this process. 913 is_sub 914 Defines the subdomains for this process or `None` to determine 915 internally. 916 is_local 917 Defines the local part of the subdomains for this process, only used 918 for `PC.ASMType.RESTRICT`. 919 920 See Also 921 -------- 922 setASMTotalSubdomains, petsc.PCASMSetLocalSubdomains 923 924 """ 925 cdef PetscInt n = asInt(nsd) 926 cdef PetscInt i = 0 927 cdef PetscIS *isets = NULL 928 cdef PetscIS *isets_local = NULL 929 if is_sub is not None: 930 assert len(is_sub) == nsd 931 CHKERR(PetscMalloc(<size_t>n*sizeof(PetscIS), &isets)) 932 for i in range(n): 933 isets[i] = (<IS?>is_sub[i]).iset 934 if is_local is not None: 935 assert len(is_local) == nsd 936 CHKERR(PetscMalloc(<size_t>n*sizeof(PetscIS), &isets_local)) 937 for i in range(n): 938 isets_local[i] = (<IS?>is_local[i]).iset 939 CHKERR(PCASMSetLocalSubdomains(self.pc, n, isets, isets_local)) 940 CHKERR(PetscFree(isets)) 941 CHKERR(PetscFree(isets_local)) 942 943 def setASMTotalSubdomains( 944 self, 945 nsd: int, 946 is_sub: Sequence[IS] | None = None, 947 is_local: Sequence[IS] | None = None) -> None: 948 """Set the subdomains for all processes. 949 950 Collective. 951 952 Parameters 953 ---------- 954 nsd 955 The number of subdomains for all processes. 956 is_sub 957 Defines the subdomains for all processes or `None` to determine 958 internally. 959 is_local 960 Defines the local part of the subdomains for this process, only used 961 for `PC.ASMType.RESTRICT`. 962 963 See Also 964 -------- 965 setASMLocalSubdomains, petsc.PCASMSetTotalSubdomains 966 967 """ 968 cdef PetscInt n = asInt(nsd) 969 cdef PetscInt i = 0 970 cdef PetscIS *isets = NULL 971 cdef PetscIS *isets_local = NULL 972 if is_sub is not None: 973 assert len(is_sub) == nsd 974 CHKERR(PetscMalloc(<size_t>n*sizeof(PetscIS), &isets)) 975 for i in range(n): 976 isets[i] = (<IS?>is_sub[i]).iset 977 if is_local is not None: 978 assert len(is_local) == nsd 979 CHKERR(PetscMalloc(<size_t>n*sizeof(PetscIS), &isets_local)) 980 for i in range(n): 981 isets_local[i] = (<IS?>is_local[i]).iset 982 CHKERR(PCASMSetTotalSubdomains(self.pc, n, isets, isets_local)) 983 CHKERR(PetscFree(isets)) 984 CHKERR(PetscFree(isets_local)) 985 986 def getASMSubKSP(self) -> list[KSP]: 987 """Return the local `KSP` object for all blocks on this process. 988 989 Not collective. 990 991 See Also 992 -------- 993 petsc.PCASMGetSubKSP 994 995 """ 996 cdef PetscInt n = 0 997 cdef PetscKSP *p = NULL 998 CHKERR(PCASMGetSubKSP(self.pc, &n, NULL, &p)) 999 return [ref_KSP(p[i]) for i from 0 <= i <n] 1000 1001 def setASMSortIndices(self, dosort: bool) -> None: 1002 """Set to sort subdomain indices. 1003 1004 Logically collective. 1005 1006 Parameters 1007 ---------- 1008 dosort 1009 Set to `True` to sort indices 1010 1011 See Also 1012 -------- 1013 petsc.PCASMSetSortIndices 1014 1015 """ 1016 cdef PetscBool cdosort = asBool(dosort) 1017 CHKERR(PCASMSetSortIndices(self.pc, cdosort)) 1018 1019 # --- GASM --- 1020 1021 def setGASMType(self, gasmtype: GASMType) -> None: 1022 """Set the type of restriction and interpolation. 1023 1024 Logically collective. 1025 1026 Parameters 1027 ---------- 1028 gasmtype 1029 The type of `GASM`. 1030 1031 See Also 1032 -------- 1033 petsc.PCGASMSetType 1034 1035 """ 1036 cdef PetscPCGASMType cval = gasmtype 1037 CHKERR(PCGASMSetType(self.pc, cval)) 1038 1039 def setGASMOverlap(self, overlap: int) -> None: 1040 """Set the overlap between a pair of subdomains. 1041 1042 Logically collective. 1043 1044 Parameters 1045 ---------- 1046 overlap 1047 The amount of overlap between subdomains. 1048 1049 See Also 1050 -------- 1051 petsc.PCGASMSetOverlap 1052 1053 """ 1054 cdef PetscInt ival = asInt(overlap) 1055 CHKERR(PCGASMSetOverlap(self.pc, ival)) 1056 1057 # --- GAMG --- 1058 1059 def setGAMGType(self, gamgtype: GAMGType | str) -> None: 1060 """Set the type of algorithm. 1061 1062 Collective. 1063 1064 Parameters 1065 ---------- 1066 gamgtype 1067 The type of `GAMG` 1068 1069 See Also 1070 -------- 1071 petsc.PCGAMGSetType 1072 1073 """ 1074 cdef PetscPCGAMGType cval = NULL 1075 gamgtype = str2bytes(gamgtype, &cval) 1076 CHKERR(PCGAMGSetType(self.pc, cval)) 1077 1078 def setGAMGLevels(self, levels: int) -> None: 1079 """Set the maximum number of levels. 1080 1081 Not collective. 1082 1083 Parameters 1084 ---------- 1085 levels 1086 The maximum number of levels to use. 1087 1088 See Also 1089 -------- 1090 petsc.PCGAMGSetNlevels 1091 1092 """ 1093 cdef PetscInt ival = asInt(levels) 1094 CHKERR(PCGAMGSetNlevels(self.pc, ival)) 1095 1096 def setGAMGSmooths(self, smooths: int) -> None: 1097 """Set the number of smoothing steps used on all levels. 1098 1099 Logically collective. 1100 1101 Parameters 1102 ---------- 1103 smooths 1104 The maximum number of smooths. 1105 1106 See Also 1107 -------- 1108 petsc.PCGAMGSetNSmooths 1109 1110 """ 1111 cdef PetscInt ival = asInt(smooths) 1112 CHKERR(PCGAMGSetNSmooths(self.pc, ival)) 1113 1114 # --- Hypre --- 1115 1116 def getHYPREType(self) -> str: 1117 """Return the `Type.HYPRE` type. 1118 1119 Not collective. 1120 1121 See Also 1122 -------- 1123 petsc.PCHYPREGetType 1124 1125 """ 1126 cdef PetscPCHYPREType cval = NULL 1127 CHKERR(PCHYPREGetType(self.pc, &cval)) 1128 return bytes2str(cval) 1129 1130 def setHYPREType(self, hypretype: str) -> None: 1131 """Set the `Type.HYPRE` type. 1132 1133 Collective. 1134 1135 Parameters 1136 ---------- 1137 hypretype 1138 The name of the type, one of ``"euclid"``, ``"pilut"``, 1139 ``"parasails"``, ``"boomeramg"``, ``"ams"``, ``"ads"`` 1140 1141 See Also 1142 -------- 1143 petsc.PCHYPRESetType 1144 1145 """ 1146 cdef PetscPCHYPREType cval = NULL 1147 hypretype = str2bytes(hypretype, &cval) 1148 CHKERR(PCHYPRESetType(self.pc, cval)) 1149 1150 def setHYPREDiscreteCurl(self, Mat mat) -> None: 1151 """Set the discrete curl matrix. 1152 1153 Collective. 1154 1155 Parameters 1156 ---------- 1157 mat 1158 The discrete curl. 1159 1160 See Also 1161 -------- 1162 petsc.PCHYPRESetDiscreteCurl 1163 1164 """ 1165 CHKERR(PCHYPRESetDiscreteCurl(self.pc, mat.mat)) 1166 1167 def setHYPREDiscreteGradient(self, Mat mat) -> None: 1168 """Set the discrete gradient matrix. 1169 1170 Collective. 1171 1172 Parameters 1173 ---------- 1174 mat 1175 The discrete gradient. 1176 1177 See Also 1178 -------- 1179 petsc.PCHYPRESetDiscreteGradient 1180 1181 """ 1182 CHKERR(PCHYPRESetDiscreteGradient(self.pc, mat.mat)) 1183 1184 def setHYPRESetAlphaPoissonMatrix(self, Mat mat) -> None: 1185 """Set the vector Poisson matrix. 1186 1187 Collective. 1188 1189 Parameters 1190 ---------- 1191 mat 1192 The vector Poisson matrix. 1193 1194 See Also 1195 -------- 1196 petsc.PCHYPRESetAlphaPoissonMatrix 1197 1198 """ 1199 CHKERR(PCHYPRESetAlphaPoissonMatrix(self.pc, mat.mat)) 1200 1201 def setHYPRESetBetaPoissonMatrix(self, Mat mat=None) -> None: 1202 """Set the Posson matrix. 1203 1204 Collective. 1205 1206 Parameters 1207 ---------- 1208 mat 1209 The Poisson matrix or `None` to turn off. 1210 1211 See Also 1212 -------- 1213 petsc.PCHYPRESetBetaPoissonMatrix 1214 1215 """ 1216 cdef PetscMat pmat = NULL 1217 if mat is not None: pmat = mat.mat 1218 CHKERR(PCHYPRESetBetaPoissonMatrix(self.pc, pmat)) 1219 1220 def setHYPRESetInterpolations(self, dim: int, Mat RT_Pi_Full=None, RT_Pi=None, 1221 Mat ND_Pi_Full=None, ND_Pi=None) -> None: 1222 """Set the interpolation matrices. 1223 1224 Collective. 1225 1226 Parameters 1227 ---------- 1228 dim 1229 The dimension of the problem. 1230 RT_Pi_Full 1231 The Raviart-Thomas interpolation matrix or `None` to omit. 1232 RT_Pi 1233 The xyz components of the Raviart-Thomas interpolation matrix, 1234 or `None` to omit. 1235 ND_Pi_Full 1236 The Nedelec interpolation matrix or `None` to omit. 1237 ND_Pi 1238 The xyz components of the Nedelec interpolation matrix, 1239 or `None` to omit. 1240 1241 See Also 1242 -------- 1243 petsc.PCHYPRESetInterpolations 1244 1245 """ 1246 cdef PetscMat RT_full_mat = NULL 1247 if RT_Pi_Full is not None: RT_full_mat = RT_Pi_Full.mat 1248 cdef PetscMat ND_full_mat = NULL 1249 if ND_Pi_Full is not None: ND_full_mat = ND_Pi_Full.mat 1250 cdef PetscInt idim = asInt(dim) 1251 cdef PetscMat *RT_Pi_mat = NULL 1252 if RT_Pi is not None: 1253 PetscMalloc(<size_t>dim*sizeof(PetscMat), &RT_Pi_mat) 1254 assert len(RT_Pi) == dim 1255 for i in range(dim): 1256 RT_Pi_mat[i] = (<Mat?>RT_Pi[i]).mat 1257 cdef PetscMat *ND_Pi_mat = NULL 1258 if ND_Pi is not None: 1259 PetscMalloc(<size_t>dim*sizeof(PetscMat), &ND_Pi_mat) 1260 assert len(ND_Pi) == dim 1261 for i in range(dim): 1262 ND_Pi_mat[dim] = (<Mat?>ND_Pi[i]).mat 1263 CHKERR (PCHYPRESetInterpolations(self.pc, idim, RT_full_mat, RT_Pi_mat, 1264 ND_full_mat, ND_Pi_mat)) 1265 CHKERR (PetscFree(RT_Pi_mat)) 1266 CHKERR (PetscFree(ND_Pi_mat)) 1267 1268 def setHYPRESetEdgeConstantVectors(self, Vec ozz, Vec zoz, Vec zzo=None) -> None: 1269 """Set the representation of the constant vector fields in the edge element basis. 1270 1271 Collective. 1272 1273 Parameters 1274 ---------- 1275 ozz 1276 A vector representing ``[1, 0, 0]`` or ``[1, 0]`` in 2D. 1277 zoz 1278 A vector representing ``[0, 1, 0]`` or ``[0, 1]`` in 2D. 1279 zzo 1280 A vector representing ``[0, 0, 1]`` or `None` in 2D. 1281 1282 See Also 1283 -------- 1284 petsc.PCHYPRESetEdgeConstantVectors 1285 1286 """ 1287 cdef PetscVec zzo_vec = NULL 1288 if zzo is not None: zzo_vec = zzo.vec 1289 CHKERR(PCHYPRESetEdgeConstantVectors(self.pc, ozz.vec, zoz.vec, 1290 zzo_vec)) 1291 1292 def setHYPREAMSSetInteriorNodes(self, Vec interior) -> None: 1293 """Set the list of interior nodes to a zero conductivity region. 1294 1295 Collective. 1296 1297 Parameters 1298 ---------- 1299 interior 1300 A vector where a value of 1.0 indicates an interior node. 1301 1302 See Also 1303 -------- 1304 petsc.PCHYPREAMSSetInteriorNodes 1305 1306 """ 1307 CHKERR(PCHYPREAMSSetInteriorNodes(self.pc, interior.vec)) 1308 1309 # --- Factor --- 1310 1311 def setFactorSolverType(self, solver: Mat.SolverType | str) -> None: 1312 """Set the solver package used to perform the factorization. 1313 1314 Logically collective. 1315 1316 Parameters 1317 ---------- 1318 solver 1319 The solver package used to factorize. 1320 1321 See Also 1322 -------- 1323 petsc.PCFactorSetMatSolverType 1324 1325 """ 1326 cdef PetscMatSolverType cval = NULL 1327 solver = str2bytes(solver, &cval) 1328 CHKERR(PCFactorSetMatSolverType(self.pc, cval)) 1329 1330 def getFactorSolverType(self) -> str: 1331 """Return the solver package used to perform the factorization. 1332 1333 Not collective. 1334 1335 See Also 1336 -------- 1337 petsc.PCFactorGetMatSolverType 1338 1339 """ 1340 cdef PetscMatSolverType cval = NULL 1341 CHKERR(PCFactorGetMatSolverType(self.pc, &cval)) 1342 return bytes2str(cval) 1343 1344 def setFactorSetUpSolverType(self) -> None: 1345 """Set up the factorization solver. 1346 1347 Collective. 1348 1349 This can be called after `KSP.setOperators` or `PC.setOperators`, causes 1350 `petsc.MatGetFactor` to be called so then one may set the options for 1351 that particular factorization object. 1352 1353 See Also 1354 -------- 1355 petsc_options, petsc.PCFactorSetUpMatSolverType 1356 1357 """ 1358 CHKERR(PCFactorSetUpMatSolverType(self.pc)) 1359 1360 def setFactorOrdering( 1361 self, 1362 ord_type: str | None = None, 1363 nzdiag: float | None = None, 1364 reuse: bool | None = None) -> None: 1365 """Set options for the matrix factorization reordering. 1366 1367 Logically collective. 1368 1369 Parameters 1370 ---------- 1371 ord_type 1372 The name of the matrix ordering or `None` to leave unchanged. 1373 nzdiag 1374 Threshold to consider diagonal entries in the matrix as zero. 1375 reuse 1376 Enable to reuse the ordering of a factored matrix. 1377 1378 See Also 1379 -------- 1380 petsc.PCFactorSetMatOrderingType 1381 petsc.PCFactorReorderForNonzeroDiagonal, petsc.PCFactorSetReuseOrdering 1382 1383 """ 1384 cdef PetscMatOrderingType cval = NULL 1385 if ord_type is not None: 1386 ord_type = str2bytes(ord_type, &cval) 1387 CHKERR(PCFactorSetMatOrderingType(self.pc, cval)) 1388 cdef PetscReal rval = 0 1389 if nzdiag is not None: 1390 rval = asReal(nzdiag) 1391 CHKERR(PCFactorReorderForNonzeroDiagonal(self.pc, rval)) 1392 cdef PetscBool bval = PETSC_FALSE 1393 if reuse is not None: 1394 bval = PETSC_TRUE if reuse else PETSC_FALSE 1395 CHKERR(PCFactorSetReuseOrdering(self.pc, bval)) 1396 1397 def setFactorPivot( 1398 self, 1399 zeropivot: float | None = None, 1400 inblocks: bool | None = None) -> None: 1401 """Set options for matrix factorization pivoting. 1402 1403 Logically collective. 1404 1405 Parameters 1406 ---------- 1407 zeropivot 1408 The size at which smaller pivots are treated as zero. 1409 inblocks 1410 Enable to allow pivoting while factoring in blocks. 1411 1412 See Also 1413 -------- 1414 petsc.PCFactorSetZeroPivot, petsc.PCFactorSetPivotInBlocks 1415 1416 """ 1417 cdef PetscReal rval = 0 1418 if zeropivot is not None: 1419 rval = asReal(zeropivot) 1420 CHKERR(PCFactorSetZeroPivot(self.pc, rval)) 1421 cdef PetscBool bval = PETSC_FALSE 1422 if inblocks is not None: 1423 bval = PETSC_TRUE if inblocks else PETSC_FALSE 1424 CHKERR(PCFactorSetPivotInBlocks(self.pc, bval)) 1425 1426 def setFactorShift( 1427 self, 1428 shift_type: Mat.FactorShiftType | None = None, 1429 amount: float | None = None) -> None: 1430 """Set options for shifting diagonal entries of a matrix. 1431 1432 Logically collective. 1433 1434 Parameters 1435 ---------- 1436 shift_type 1437 The type of shift, or `None` to leave unchanged. 1438 amount 1439 The amount of shift. Specify `DEFAULT` to determine internally or 1440 `None` to leave unchanged. 1441 1442 See Also 1443 -------- 1444 petsc.PCFactorSetShiftType, petsc.PCFactorSetShiftAmount 1445 1446 """ 1447 cdef PetscMatFactorShiftType cval = MAT_SHIFT_NONE 1448 if shift_type is not None: 1449 cval = matfactorshifttype(shift_type) 1450 CHKERR(PCFactorSetShiftType(self.pc, cval)) 1451 cdef PetscReal rval = 0 1452 if amount is not None: 1453 rval = asReal(amount) 1454 CHKERR(PCFactorSetShiftAmount(self.pc, rval)) 1455 1456 def setFactorLevels(self, levels: int) -> None: 1457 """Set the number of levels of fill. 1458 1459 Logically collective. 1460 1461 Parameters 1462 ---------- 1463 levels 1464 The number of levels to fill. 1465 1466 See Also 1467 -------- 1468 petsc.PCFactorSetLevels 1469 1470 """ 1471 cdef PetscInt ival = asInt(levels) 1472 CHKERR(PCFactorSetLevels(self.pc, ival)) 1473 1474 def getFactorMatrix(self) -> Mat: 1475 """Return the factored matrix. 1476 1477 Not collective. 1478 1479 See Also 1480 -------- 1481 petsc.PCFactorGetMatrix 1482 1483 """ 1484 cdef Mat mat = Mat() 1485 CHKERR(PCFactorGetMatrix(self.pc, &mat.mat)) 1486 CHKERR(PetscINCREF(mat.obj)) 1487 return mat 1488 1489 # --- FieldSplit --- 1490 1491 def setFieldSplitType(self, ctype: CompositeType) -> None: 1492 """Set the type of composition of a field split preconditioner. 1493 1494 Collective. 1495 1496 Parameters 1497 ---------- 1498 ctype 1499 The type of composition. 1500 1501 See Also 1502 -------- 1503 petsc.PCFieldSplitSetType 1504 1505 """ 1506 cdef PetscPCCompositeType cval = ctype 1507 CHKERR(PCFieldSplitSetType(self.pc, cval)) 1508 1509 def setFieldSplitIS(self, *fields: tuple[str, IS]) -> None: 1510 """Set the elements for the field split by `IS`. 1511 1512 Logically collective. 1513 1514 Solve options for this split will be available 1515 under the prefix ``-fieldsplit_SPLITNAME_*``. 1516 1517 Parameters 1518 ---------- 1519 fields 1520 A sequence of tuples containing the split name and the `IS` that 1521 defines the elements in the split. 1522 1523 See Also 1524 -------- 1525 petsc_options, petsc.PCFieldSplitSetIS 1526 1527 """ 1528 cdef object name = None 1529 cdef IS field = None 1530 cdef const char *cname = NULL 1531 for name, field in fields: 1532 name = str2bytes(name, &cname) 1533 CHKERR(PCFieldSplitSetIS(self.pc, cname, field.iset)) 1534 1535 def setFieldSplitFields(self, bsize: int, *fields: tuple[str, Sequence[int]]) -> None: 1536 """Sets the elements for the field split. 1537 1538 Collective. 1539 1540 Parameters 1541 ---------- 1542 bsize 1543 The block size 1544 fields 1545 A sequence of tuples containing the split name and a sequence of 1546 integers that define the elements in the split. 1547 1548 See Also 1549 -------- 1550 petsc.PCFieldSplitSetBlockSize, petsc.PCFieldSplitSetFields 1551 1552 """ 1553 cdef PetscInt bs = asInt(bsize) 1554 CHKERR(PCFieldSplitSetBlockSize(self.pc, bs)) 1555 cdef object name = None 1556 cdef object field = None 1557 cdef const char *cname = NULL 1558 cdef PetscInt nfields = 0, *ifields = NULL 1559 for name, field in fields: 1560 name = str2bytes(name, &cname) 1561 field = iarray_i(field, &nfields, &ifields) 1562 CHKERR(PCFieldSplitSetFields(self.pc, cname, 1563 nfields, ifields, ifields)) 1564 1565 def getFieldSplitSubKSP(self) -> list[KSP]: 1566 """Return the `KSP` for all splits. 1567 1568 Not collective. 1569 1570 See Also 1571 -------- 1572 petsc.PCFieldSplitGetSubKSP 1573 1574 """ 1575 cdef PetscInt n = 0 1576 cdef PetscKSP *p = NULL 1577 cdef object subksp = None 1578 try: 1579 CHKERR(PCFieldSplitGetSubKSP(self.pc, &n, &p)) 1580 subksp = [ref_KSP(p[i]) for i from 0 <= i <n] 1581 finally: 1582 CHKERR(PetscFree(p)) 1583 return subksp 1584 1585 def getFieldSplitSchurGetSubKSP(self) -> list[KSP]: 1586 """Return the `KSP` for the Schur complement based splits. 1587 1588 Not collective. 1589 1590 See Also 1591 -------- 1592 petsc.PCFieldSplitSchurGetSubKSP, petsc.PCFieldSplitGetSubKSP 1593 1594 """ 1595 cdef PetscInt n = 0 1596 cdef PetscKSP *p = NULL 1597 cdef object subksp = None 1598 try: 1599 CHKERR(PCFieldSplitSchurGetSubKSP(self.pc, &n, &p)) 1600 subksp = [ref_KSP(p[i]) for i from 0 <= i <n] 1601 finally: 1602 CHKERR(PetscFree(p)) 1603 return subksp 1604 1605 def getFieldSplitSubIS(self, splitname: str) -> IS: 1606 """Return the `IS` associated with a given name. 1607 1608 Not collective. 1609 1610 See Also 1611 -------- 1612 petsc.PCFieldSplitGetIS 1613 1614 """ 1615 cdef PetscIS field = NULL 1616 cdef const char *cname = NULL 1617 splitname = str2bytes(splitname, &cname) 1618 CHKERR(PCFieldSplitGetIS(self.pc, cname, &field)) 1619 return ref_IS(field) 1620 1621 def setFieldSplitSchurFactType(self, ctype: FieldSplitSchurFactType) -> None: 1622 """Set the type of approximate block factorization. 1623 1624 Collective. 1625 1626 Parameters 1627 ---------- 1628 ctype 1629 The type indicating which blocks to retain. 1630 1631 See Also 1632 -------- 1633 petsc.PCFieldSplitSetSchurFactType 1634 1635 """ 1636 cdef PetscPCFieldSplitSchurFactType cval = ctype 1637 CHKERR(PCFieldSplitSetSchurFactType(self.pc, cval)) 1638 1639 def setFieldSplitSchurPreType( 1640 self, 1641 ptype: FieldSplitSchurPreType, 1642 Mat pre=None) -> None: 1643 """Set from what operator the `PC` is constructed. 1644 1645 Collective. 1646 1647 Parameters 1648 ---------- 1649 ptype 1650 The type of matrix to use for preconditioning the Schur complement. 1651 pre 1652 The optional matrix to use for preconditioning. 1653 1654 See Also 1655 -------- 1656 petsc.PCFieldSplitSetSchurPre 1657 1658 """ 1659 cdef PetscPCFieldSplitSchurPreType pval = ptype 1660 cdef PetscMat pmat = NULL 1661 if pre is not None: pmat = pre.mat 1662 CHKERR(PCFieldSplitSetSchurPre(self.pc, pval, pmat)) 1663 1664 # --- COMPOSITE --- 1665 1666 def setCompositeType(self, ctype: CompositeType) -> None: 1667 """Set the type of composite preconditioner. 1668 1669 Logically collective. 1670 1671 Parameters 1672 ---------- 1673 ctype 1674 The type of composition. 1675 1676 See Also 1677 -------- 1678 petsc.PCCompositeSetType 1679 1680 """ 1681 cdef PetscPCCompositeType cval = ctype 1682 CHKERR(PCCompositeSetType(self.pc, cval)) 1683 1684 def getCompositePC(self, n: int) -> None: 1685 """Return a component of the composite `PC`. 1686 1687 Not collective. 1688 1689 Parameters 1690 ---------- 1691 n 1692 The index of the `PC` in the composition. 1693 1694 See Also 1695 -------- 1696 petsc.PCCompositeGetPC 1697 1698 """ 1699 cdef PC pc = PC() 1700 cdef cn = asInt(n) 1701 CHKERR(PCCompositeGetPC(self.pc, cn, &pc.pc)) 1702 CHKERR(PetscINCREF(pc.obj)) 1703 return pc 1704 1705 def addCompositePCType(self, pc_type: Type | str) -> None: 1706 """Add a `PC` of the given type to the composite `PC`. 1707 1708 Collective. 1709 1710 Parameters 1711 ---------- 1712 pc_type 1713 The type of the preconditioner to add. 1714 1715 See Also 1716 -------- 1717 petsc.PCCompositeAddPCType 1718 1719 """ 1720 cdef PetscPCType cval = NULL 1721 pc_type = str2bytes(pc_type, &cval) 1722 CHKERR(PCCompositeAddPCType(self.pc, cval)) 1723 1724 # --- KSP --- 1725 1726 def getKSP(self) -> KSP: 1727 """Return the `KSP` if the `PC` is `Type.KSP`. 1728 1729 Not collective. 1730 1731 See Also 1732 -------- 1733 petsc.PCKSPGetKSP 1734 1735 """ 1736 cdef KSP ksp = KSP() 1737 CHKERR(PCKSPGetKSP(self.pc, &ksp.ksp)) 1738 CHKERR(PetscINCREF(ksp.obj)) 1739 return ksp 1740 1741 # --- MG --- 1742 1743 def getMGType(self) -> MGType: 1744 """Return the form of multigrid. 1745 1746 Logically collective. 1747 1748 See Also 1749 -------- 1750 petsc.PCMGGetType 1751 1752 """ 1753 cdef PetscPCMGType cval = PC_MG_ADDITIVE 1754 CHKERR(PCMGGetType(self.pc, &cval)) 1755 return cval 1756 1757 def setMGType(self, mgtype: MGType) -> None: 1758 """Set the form of multigrid. 1759 1760 Logically collective. 1761 1762 See Also 1763 -------- 1764 petsc.PCMGSetType 1765 1766 """ 1767 cdef PetscPCMGType cval = mgtype 1768 CHKERR(PCMGSetType(self.pc, cval)) 1769 1770 def getMGLevels(self) -> int: 1771 """Return the number of `MG` levels. 1772 1773 Not collective. 1774 1775 See Also 1776 -------- 1777 petsc.PCMGGetLevels 1778 1779 """ 1780 cdef PetscInt levels = 0 1781 CHKERR(PCMGGetLevels(self.pc, &levels)) 1782 return toInt(levels) 1783 1784 def setMGLevels(self, levels: int) -> None: 1785 """Set the number of `MG` levels. 1786 1787 Logically collective. 1788 1789 Parameters 1790 ---------- 1791 levels 1792 The number of levels 1793 1794 See Also 1795 -------- 1796 petsc.PCMGSetLevels 1797 1798 """ 1799 cdef PetscInt clevels = asInt(levels) 1800 CHKERR(PCMGSetLevels(self.pc, clevels, NULL)) 1801 1802 def getMGCoarseSolve(self) -> KSP: 1803 """Return the `KSP` used on the coarse grid. 1804 1805 Not collective. 1806 1807 See Also 1808 -------- 1809 petsc.PCMGGetCoarseSolve 1810 1811 """ 1812 cdef KSP ksp = KSP() 1813 CHKERR(PCMGGetCoarseSolve(self.pc, &ksp.ksp)) 1814 CHKERR(PetscINCREF(ksp.obj)) 1815 return ksp 1816 1817 def setMGInterpolation(self, level, Mat mat) -> None: 1818 """Set the interpolation operator for the given level. 1819 1820 Logically collective. 1821 1822 Parameters 1823 ---------- 1824 level 1825 The level where interpolation is defined from ``level-1`` to ``level``. 1826 mat 1827 The interpolation operator 1828 1829 See Also 1830 -------- 1831 petsc.PCMGSetInterpolation 1832 1833 """ 1834 cdef PetscInt clevel = asInt(level) 1835 CHKERR(PCMGSetInterpolation(self.pc, clevel, mat.mat)) 1836 1837 def getMGInterpolation(self, level: int) -> Mat: 1838 """Return the interpolation operator for the given level. 1839 1840 Logically collective. 1841 1842 Parameters 1843 ---------- 1844 level 1845 The level where interpolation is defined from ``level-1`` to ``level``. 1846 1847 See Also 1848 -------- 1849 petsc.PCMGGetInterpolation 1850 1851 """ 1852 cdef PetscInt clevel = asInt(level) 1853 cdef Mat interpolation = Mat() 1854 CHKERR(PCMGGetInterpolation(self.pc, clevel, &interpolation.mat)) 1855 CHKERR(PetscINCREF(interpolation.obj)) 1856 return interpolation 1857 1858 def setMGRestriction(self, level: int, Mat mat) -> None: 1859 """Set the restriction operator for the given level. 1860 1861 Logically collective. 1862 1863 Parameters 1864 ---------- 1865 level 1866 The level where restriction is defined from ``level`` to ``level-1``. 1867 mat 1868 The restriction operator 1869 1870 See Also 1871 -------- 1872 petsc.PCMGSetRestriction 1873 1874 """ 1875 cdef PetscInt clevel = asInt(level) 1876 CHKERR(PCMGSetRestriction(self.pc, clevel, mat.mat)) 1877 1878 def getMGRestriction(self, level: int) -> Mat: 1879 """Return the restriction operator for the given level. 1880 1881 Logically collective. 1882 1883 Parameters 1884 ---------- 1885 level 1886 The level where restriction is defined from ``level`` to ``level-1``. 1887 1888 See Also 1889 -------- 1890 petsc.PCMGGetRestriction 1891 1892 """ 1893 cdef PetscInt clevel = asInt(level) 1894 cdef Mat restriction = Mat() 1895 CHKERR(PCMGGetRestriction(self.pc, clevel, &restriction.mat)) 1896 CHKERR(PetscINCREF(restriction.obj)) 1897 return restriction 1898 1899 def setMGRScale(self, level: int, Vec rscale) -> None: 1900 """Set the pointwise scaling for the restriction operator on the given level. 1901 1902 Logically collective. 1903 1904 Parameters 1905 ---------- 1906 level 1907 The level where restriction is defined from ``level`` to ``level-1``. 1908 rscale 1909 The scaling vector. 1910 1911 See Also 1912 -------- 1913 petsc.PCMGSetRScale 1914 1915 """ 1916 cdef PetscInt clevel = asInt(level) 1917 CHKERR(PCMGSetRScale(self.pc, clevel, rscale.vec)) 1918 1919 def getMGRScale(self, level: int) -> Vec: 1920 """Return the pointwise scaling for the restriction operator on the given level. 1921 1922 Logically collective. 1923 1924 Parameters 1925 ---------- 1926 level 1927 The level where restriction is defined from ``level`` to ``level-1``. 1928 1929 See Also 1930 -------- 1931 petsc.PCMGGetRScale 1932 1933 """ 1934 cdef PetscInt clevel = asInt(level) 1935 cdef Vec rscale = Vec() 1936 CHKERR(PCMGGetRScale(self.pc, clevel, &rscale.vec)) 1937 CHKERR(PetscINCREF(rscale.obj)) 1938 return rscale 1939 1940 def getMGSmoother(self, level: int) -> KSP: 1941 """Return the `KSP` to be used as a smoother. 1942 1943 Not collective. 1944 1945 Parameters 1946 ---------- 1947 level 1948 The level of the smoother. 1949 1950 See Also 1951 -------- 1952 getMGSmootherDown, getMGSmootherUp, petsc.PCMGGetSmoother 1953 1954 """ 1955 cdef PetscInt clevel = asInt(level) 1956 cdef KSP ksp = KSP() 1957 CHKERR(PCMGGetSmoother(self.pc, clevel, &ksp.ksp)) 1958 CHKERR(PetscINCREF(ksp.obj)) 1959 return ksp 1960 1961 def getMGSmootherDown(self, level: int) -> KSP: 1962 """Return the `KSP` to be used as a smoother before coarse grid correction. 1963 1964 Not collective. 1965 1966 Parameters 1967 ---------- 1968 level 1969 The level of the smoother. 1970 1971 See Also 1972 -------- 1973 getMGSmoother, getMGSmootherUp, petsc.PCMGGetSmootherDown 1974 1975 """ 1976 cdef PetscInt clevel = asInt(level) 1977 cdef KSP ksp = KSP() 1978 CHKERR(PCMGGetSmootherDown(self.pc, clevel, &ksp.ksp)) 1979 CHKERR(PetscINCREF(ksp.obj)) 1980 return ksp 1981 1982 def getMGSmootherUp(self, level: int) -> KSP: 1983 """Return the `KSP` to be used as a smoother after coarse grid correction. 1984 1985 Not collective. 1986 1987 Parameters 1988 ---------- 1989 level 1990 The level of the smoother. 1991 1992 See Also 1993 -------- 1994 getMGSmootherDown, getMGSmoother, petsc.PCMGGetSmootherUp 1995 1996 """ 1997 cdef PetscInt clevel = asInt(level) 1998 cdef KSP ksp = KSP() 1999 CHKERR(PCMGGetSmootherUp(self.pc, clevel, &ksp.ksp)) 2000 CHKERR(PetscINCREF(ksp.obj)) 2001 return ksp 2002 2003 def setMGCycleType(self, cycle_type: MGCycleType) -> None: 2004 """Set the type of cycles. 2005 2006 Logically collective. 2007 2008 Parameters 2009 ---------- 2010 cycle_type 2011 The type of multigrid cycles to use. 2012 2013 See Also 2014 -------- 2015 setMGCycleTypeOnLevel, petsc.PCMGSetCycleType 2016 2017 """ 2018 cdef PetscPCMGCycleType ctype = cycle_type 2019 CHKERR(PCMGSetCycleType(self.pc, ctype)) 2020 2021 def setMGCycleTypeOnLevel(self, level: int, cycle_type: MGCycleType) -> None: 2022 """Set the type of cycle on the given level. 2023 2024 Logically collective. 2025 2026 Parameters 2027 ---------- 2028 level 2029 The level on which to set the cycle type. 2030 cycle_type 2031 The type of multigrid cycles to use. 2032 2033 See Also 2034 -------- 2035 setMGCycleType, petsc.PCMGSetCycleTypeOnLevel 2036 2037 """ 2038 cdef PetscInt clevel = asInt(level) 2039 cdef PetscPCMGCycleType ctype = cycle_type 2040 CHKERR(PCMGSetCycleTypeOnLevel(self.pc, clevel, ctype)) 2041 2042 def setMGRhs(self, level: int, Vec rhs) -> None: 2043 """Set the vector where the right-hand side is stored. 2044 2045 Logically collective. 2046 2047 If not provided, one will be set internally. Will 2048 be cleaned up in `destroy`. 2049 2050 Parameters 2051 ---------- 2052 level 2053 The level on which to set the right-hand side. 2054 rhs 2055 The vector where the right-hand side is stored. 2056 2057 See Also 2058 -------- 2059 petsc.PCMGSetRhs 2060 2061 """ 2062 cdef PetscInt clevel = asInt(level) 2063 CHKERR(PCMGSetRhs(self.pc, clevel, rhs.vec)) 2064 2065 def setMGX(self, level: int, Vec x) -> None: 2066 """Set the vector where the solution is stored. 2067 2068 Logically collective. 2069 2070 If not provided, one will be set internally. Will 2071 be cleaned up in `destroy`. 2072 2073 Parameters 2074 ---------- 2075 level 2076 The level on which to set the solution. 2077 x 2078 The vector where the solution is stored. 2079 2080 See Also 2081 -------- 2082 petsc.PCMGSetX 2083 2084 """ 2085 cdef PetscInt clevel = asInt(level) 2086 CHKERR(PCMGSetX(self.pc, clevel, x.vec)) 2087 2088 def setMGR(self, level: int, Vec r) -> None: 2089 """Set the vector where the residual is stored. 2090 2091 Logically collective. 2092 2093 If not provided, one will be set internally. Will 2094 be cleaned up in `destroy`. 2095 2096 Parameters 2097 ---------- 2098 level 2099 The level on which to set the residual. 2100 r 2101 The vector where the residual is stored. 2102 2103 See Also 2104 -------- 2105 petsc.PCMGSetR 2106 2107 """ 2108 cdef PetscInt clevel = asInt(level) 2109 CHKERR(PCMGSetR(self.pc, clevel, r.vec)) 2110 2111 # --- BDDC --- 2112 2113 def setBDDCLocalAdjacency(self, csr: CSRIndicesSpec) -> None: 2114 """Provide a custom connectivity graph for local dofs. 2115 2116 Not collective. 2117 2118 Parameters 2119 ---------- 2120 csr 2121 Compressed sparse row layout information. 2122 2123 See Also 2124 -------- 2125 petsc.PCBDDCSetLocalAdjacencyGraph 2126 2127 """ 2128 cdef object oi, oj 2129 cdef PetscInt ni=0, *i=NULL 2130 cdef PetscInt nj=0, *j=NULL 2131 oi, oj = csr 2132 oi = iarray_i(oi, &ni, &i) 2133 oj = iarray_i(oj, &nj, &j) 2134 if (i[0] != 0): 2135 raise ValueError("I[0] is %d, expected %d" % 2136 (toInt(i[0]), toInt(0))) 2137 if (i[ni-1] != nj): 2138 raise ValueError("size(J) is %d, expected %d" % 2139 (toInt(nj), toInt(i[ni-1]))) 2140 2141 CHKERR(PCBDDCSetLocalAdjacencyGraph(self.pc, ni - 1, i, j, PETSC_COPY_VALUES)) 2142 2143 def setBDDCDivergenceMat(self, Mat div, trans: bool = False, IS l2l=None) -> None: 2144 """Set the linear operator representing ∫ div(u)•p dx. 2145 2146 Collective. 2147 2148 Parameters 2149 ---------- 2150 div 2151 The matrix in `Mat.Type.IS` format. 2152 trans 2153 If `True`, the pressure/velocity is in the trial/test space 2154 respectively. If `False` the pressure/velocity is in the test/trial 2155 space. 2156 l2l 2157 Optional `IS` describing the local to local map for velocities. 2158 2159 See Also 2160 -------- 2161 petsc.PCBDDCSetDivergenceMat 2162 2163 """ 2164 cdef PetscBool ptrans = trans 2165 cdef PetscIS pl2l = NULL 2166 if l2l is not None: pl2l = l2l.iset 2167 CHKERR(PCBDDCSetDivergenceMat(self.pc, div.mat, ptrans, pl2l)) 2168 2169 def setBDDCDiscreteGradient( 2170 self, 2171 Mat G, 2172 order: int = 1, 2173 field: int = 1, 2174 gord: bool = True, 2175 conforming: bool = True) -> None: 2176 """Set the discrete gradient. 2177 2178 Collective. 2179 2180 Parameters 2181 ---------- 2182 G 2183 The discrete gradient matrix in `Mat.Type.AIJ` format. 2184 order 2185 The order of the Nedelec space. 2186 field 2187 The field number of the Nedelec degrees of freedom. This is not used 2188 if no fields have been specified. 2189 gord 2190 Enable to use global ordering in the rows of ``G``. 2191 conforming 2192 Enable if the mesh is conforming. 2193 2194 See Also 2195 -------- 2196 petsc.PCBDDCSetDiscreteGradient 2197 2198 """ 2199 cdef PetscInt porder = asInt(order) 2200 cdef PetscInt pfield = asInt(field) 2201 cdef PetscBool pgord = gord 2202 cdef PetscBool pconforming = conforming 2203 CHKERR(PCBDDCSetDiscreteGradient(self.pc, G.mat, porder, pfield, pgord, pconforming)) 2204 2205 def setBDDCChangeOfBasisMat(self, Mat T, interior: bool = False) -> None: 2206 """Set a user defined change of basis for degrees of freedom. 2207 2208 Collective. 2209 2210 Parameters 2211 ---------- 2212 T 2213 The matrix representing the change of basis. 2214 interior 2215 Enable to indicate the change of basis affects interior degrees of 2216 freedom. 2217 2218 See Also 2219 -------- 2220 petsc.PCBDDCSetChangeOfBasisMat 2221 2222 """ 2223 cdef PetscBool pinterior = interior 2224 CHKERR(PCBDDCSetChangeOfBasisMat(self.pc, T.mat, pinterior)) 2225 2226 def setBDDCPrimalVerticesIS(self, IS primv) -> None: 2227 """Set additional user defined primal vertices. 2228 2229 Collective. 2230 2231 Parameters 2232 ---------- 2233 primv 2234 The `IS` of primal vertices in global numbering. 2235 2236 See Also 2237 -------- 2238 petsc.PCBDDCSetPrimalVerticesIS 2239 2240 """ 2241 CHKERR(PCBDDCSetPrimalVerticesIS(self.pc, primv.iset)) 2242 2243 def setBDDCPrimalVerticesLocalIS(self, IS primv) -> None: 2244 """Set additional user defined primal vertices. 2245 2246 Collective. 2247 2248 Parameters 2249 ---------- 2250 primv 2251 The `IS` of primal vertices in local numbering. 2252 2253 See Also 2254 -------- 2255 petsc.PCBDDCSetPrimalVerticesLocalIS 2256 2257 """ 2258 CHKERR(PCBDDCSetPrimalVerticesLocalIS(self.pc, primv.iset)) 2259 2260 def setBDDCCoarseningRatio(self, cratio: int) -> None: 2261 """Set the coarsening ratio used in the multilevel version. 2262 2263 Logically collective. 2264 2265 Parameters 2266 ---------- 2267 cratio 2268 The coarsening ratio at the coarse level 2269 2270 See Also 2271 -------- 2272 petsc.PCBDDCSetCoarseningRatio 2273 2274 """ 2275 cdef PetscInt pcratio = asInt(cratio) 2276 CHKERR(PCBDDCSetCoarseningRatio(self.pc, pcratio)) 2277 2278 def setBDDCLevels(self, levels: int) -> None: 2279 """Set the maximum number of additional levels allowed. 2280 2281 Logically collective. 2282 2283 Parameters 2284 ---------- 2285 levels 2286 The maximum number of levels. 2287 2288 See Also 2289 -------- 2290 petsc.PCBDDCSetLevels 2291 2292 """ 2293 cdef PetscInt plevels = asInt(levels) 2294 CHKERR(PCBDDCSetLevels(self.pc, plevels)) 2295 2296 def setBDDCDirichletBoundaries(self, IS bndr) -> None: 2297 """Set the `IS` defining Dirichlet boundaries for the global problem. 2298 2299 Collective. 2300 2301 Parameters 2302 ---------- 2303 bndr 2304 The parallel `IS` defining Dirichlet boundaries. 2305 2306 See Also 2307 -------- 2308 petsc.PCBDDCSetDirichletBoundaries 2309 2310 """ 2311 CHKERR(PCBDDCSetDirichletBoundaries(self.pc, bndr.iset)) 2312 2313 def setBDDCDirichletBoundariesLocal(self, IS bndr) -> None: 2314 """Set the `IS` defining Dirichlet boundaries in local ordering. 2315 2316 Collective. 2317 2318 Parameters 2319 ---------- 2320 bndr 2321 The parallel `IS` defining Dirichlet boundaries in local ordering. 2322 2323 See Also 2324 -------- 2325 setBDDCDirichletBoundaries, petsc.PCBDDCSetDirichletBoundariesLocal 2326 2327 """ 2328 CHKERR(PCBDDCSetDirichletBoundariesLocal(self.pc, bndr.iset)) 2329 2330 def setBDDCNeumannBoundaries(self, IS bndr) -> None: 2331 """Set the `IS` defining Neumann boundaries for the global problem. 2332 2333 Collective. 2334 2335 Parameters 2336 ---------- 2337 bndr 2338 The parallel `IS` defining Neumann boundaries. 2339 2340 See Also 2341 -------- 2342 petsc.PCBDDCSetNeumannBoundaries 2343 2344 """ 2345 CHKERR(PCBDDCSetNeumannBoundaries(self.pc, bndr.iset)) 2346 2347 def setBDDCNeumannBoundariesLocal(self, IS bndr) -> None: 2348 """Set the `IS` defining Neumann boundaries in local ordering. 2349 2350 Collective. 2351 2352 Parameters 2353 ---------- 2354 bndr 2355 The parallel `IS` defining Neumann boundaries in local ordering. 2356 2357 See Also 2358 -------- 2359 setBDDCNeumannBoundaries, petsc.PCBDDCSetNeumannBoundariesLocal 2360 2361 """ 2362 CHKERR(PCBDDCSetNeumannBoundariesLocal(self.pc, bndr.iset)) 2363 2364 def setBDDCDofsSplitting(self, isfields: IS | Sequence[IS]) -> None: 2365 """Set the index set(s) defining fields of the global matrix. 2366 2367 Collective. 2368 2369 Parameters 2370 ---------- 2371 isfields 2372 The sequence of `IS` describing the fields in global ordering. 2373 2374 See Also 2375 -------- 2376 petsc.PCBDDCSetDofsSplitting 2377 2378 """ 2379 isfields = [isfields] if isinstance(isfields, IS) else list(isfields) 2380 cdef Py_ssize_t i, n = len(isfields) 2381 cdef PetscIS *cisfields = NULL 2382 cdef object unused = oarray_p(empty_p(<PetscInt>n), NULL, <void**>&cisfields) 2383 for i from 0 <= i < n: cisfields[i] = (<IS?>isfields[i]).iset 2384 CHKERR(PCBDDCSetDofsSplitting(self.pc, <PetscInt>n, cisfields)) 2385 2386 def setBDDCDofsSplittingLocal(self, isfields: IS | Sequence[IS]) -> None: 2387 """Set the index set(s) defining fields of the local subdomain matrix. 2388 2389 Collective. 2390 2391 Not all nodes need to be listed. Unlisted nodes will belong 2392 to the complement field. 2393 2394 Parameters 2395 ---------- 2396 isfields 2397 The sequence of `IS` describing the fields in local ordering. 2398 2399 See Also 2400 -------- 2401 petsc.PCBDDCSetDofsSplittingLocal 2402 2403 """ 2404 isfields = [isfields] if isinstance(isfields, IS) else list(isfields) 2405 cdef Py_ssize_t i, n = len(isfields) 2406 cdef PetscIS *cisfields = NULL 2407 cdef object unused = oarray_p(empty_p(<PetscInt>n), NULL, <void**>&cisfields) 2408 for i from 0 <= i < n: cisfields[i] = (<IS?>isfields[i]).iset 2409 CHKERR(PCBDDCSetDofsSplittingLocal(self.pc, <PetscInt>n, cisfields)) 2410 2411 # --- Patch --- 2412 def getPatchSubKSP(self) -> list[KSP]: 2413 """Return the local `KSP` object for all blocks on this process. 2414 2415 Not collective. 2416 2417 """ 2418 cdef PetscInt n = 0 2419 cdef PetscKSP *p = NULL 2420 CHKERR(PCPatchGetSubKSP(self.pc, &n, &p)) 2421 return [ref_KSP(p[i]) for i from 0 <= i <n] 2422 2423 def setPatchCellNumbering(self, Section sec) -> None: 2424 """Set the cell numbering.""" 2425 CHKERR(PCPatchSetCellNumbering(self.pc, sec.sec)) 2426 2427 def setPatchDiscretisationInfo(self, dms, bs, 2428 cellNodeMaps, 2429 subspaceOffsets, 2430 ghostBcNodes, 2431 globalBcNodes) -> None: 2432 """Set discretisation info.""" 2433 cdef PetscInt numSubSpaces = 0 2434 cdef PetscInt numGhostBcs = 0, numGlobalBcs = 0 2435 cdef PetscInt *nodesPerCell = NULL 2436 cdef const PetscInt **ccellNodeMaps = NULL 2437 cdef PetscDM *cdms = NULL 2438 cdef PetscInt *cbs = NULL 2439 cdef PetscInt *csubspaceOffsets = NULL 2440 cdef PetscInt *cghostBcNodes = NULL 2441 cdef PetscInt *cglobalBcNodes = NULL 2442 cdef PetscInt i = 0 2443 2444 bs = iarray_i(bs, &numSubSpaces, &cbs) 2445 ghostBcNodes = iarray_i(ghostBcNodes, &numGhostBcs, &cghostBcNodes) 2446 globalBcNodes = iarray_i(globalBcNodes, &numGlobalBcs, &cglobalBcNodes) 2447 subspaceOffsets = iarray_i(subspaceOffsets, NULL, &csubspaceOffsets) 2448 2449 CHKERR(PetscMalloc(<size_t>numSubSpaces*sizeof(PetscInt), &nodesPerCell)) 2450 CHKERR(PetscMalloc(<size_t>numSubSpaces*sizeof(PetscDM), &cdms)) 2451 CHKERR(PetscMalloc(<size_t>numSubSpaces*sizeof(PetscInt*), &ccellNodeMaps)) 2452 for i in range(numSubSpaces): 2453 cdms[i] = (<DM?>dms[i]).dm 2454 _, nodes = asarray(cellNodeMaps[i]).shape 2455 cellNodeMaps[i] = iarray_i(cellNodeMaps[i], NULL, <PetscInt**>&ccellNodeMaps[i]) 2456 nodesPerCell[i] = asInt(nodes) 2457 2458 # TODO: refactor on the PETSc side to take ISes? 2459 CHKERR(PCPatchSetDiscretisationInfo(self.pc, numSubSpaces, 2460 cdms, cbs, nodesPerCell, 2461 ccellNodeMaps, csubspaceOffsets, 2462 numGhostBcs, cghostBcNodes, 2463 numGlobalBcs, cglobalBcNodes)) 2464 CHKERR(PetscFree(nodesPerCell)) 2465 CHKERR(PetscFree(cdms)) 2466 CHKERR(PetscFree(ccellNodeMaps)) 2467 2468 def setPatchComputeOperator(self, operator, args=None, kargs=None) -> None: 2469 """Set compute operator callbacks.""" 2470 if args is None: args = () 2471 if kargs is None: kargs = {} 2472 context = (operator, args, kargs) 2473 self.set_attr("__patch_compute_operator__", context) 2474 CHKERR(PCPatchSetComputeOperator(self.pc, PCPatch_ComputeOperator, <void*>context)) 2475 2476 def setPatchComputeOperatorInteriorFacets(self, operator, args=None, kargs=None) -> None: 2477 """Set compute operator callbacks.""" 2478 if args is None: args = () 2479 if kargs is None: kargs = {} 2480 context = (operator, args, kargs) 2481 self.set_attr("__patch_compute_operator_interior_facets__", context) 2482 CHKERR(PCPatchSetComputeOperatorInteriorFacets(self.pc, PCPatch_ComputeOperatorInteriorFacets, <void*>context)) 2483 2484 def setPatchComputeFunction(self, function, args=None, kargs=None) -> None: 2485 """Set compute operator callbacks.""" 2486 if args is None: args = () 2487 if kargs is None: kargs = {} 2488 context = (function, args, kargs) 2489 self.set_attr("__patch_compute_function__", context) 2490 CHKERR(PCPatchSetComputeFunction(self.pc, PCPatch_ComputeFunction, <void*>context)) 2491 2492 def setPatchComputeFunctionInteriorFacets(self, function, args=None, kargs=None) -> None: 2493 """Set compute operator callbacks.""" 2494 if args is None: args = () 2495 if kargs is None: kargs = {} 2496 context = (function, args, kargs) 2497 self.set_attr("__patch_compute_function_interior_facets__", context) 2498 CHKERR(PCPatchSetComputeFunction(self.pc, PCPatch_ComputeFunctionInteriorFacets, <void*>context)) 2499 2500 def setPatchConstructType(self, typ, operator=None, args=None, kargs=None) -> None: 2501 """Set compute operator callbacks.""" 2502 if args is None: args = () 2503 if kargs is None: kargs = {} 2504 2505 if typ in {PC.PatchConstructType.PYTHON, PC.PatchConstructType.USER} and operator is None: 2506 raise ValueError("Must provide operator for USER or PYTHON type") 2507 if operator is not None: 2508 context = (operator, args, kargs) 2509 else: 2510 context = None 2511 self.set_attr("__patch_construction_operator__", context) 2512 CHKERR(PCPatchSetConstructType(self.pc, typ, PCPatch_UserConstructOperator, <void*>context)) 2513 2514 # --- HPDDM --- 2515 2516 def setHPDDMAuxiliaryMat(self, IS uis, Mat uaux) -> None: 2517 """Set the auxiliary matrix used by the preconditioner. 2518 2519 Logically collective. 2520 2521 Parameters 2522 ---------- 2523 uis 2524 The `IS` of the local auxiliary matrix 2525 uaux 2526 The auxiliary sequential matrix 2527 2528 See Also 2529 -------- 2530 petsc.PCHPDDMSetAuxiliaryMat 2531 2532 """ 2533 CHKERR(PCHPDDMSetAuxiliaryMat(self.pc, uis.iset, uaux.mat, NULL, <void*>NULL)) 2534 2535 def setHPDDMRHSMat(self, Mat B) -> None: 2536 """Set the right-hand side matrix of the preconditioner. 2537 2538 Logically collective. 2539 2540 Parameters 2541 ---------- 2542 B 2543 The right-hand side sequential matrix. 2544 2545 See Also 2546 -------- 2547 petsc.PCHPDDMSetRHSMat 2548 2549 """ 2550 CHKERR(PCHPDDMSetRHSMat(self.pc, B.mat)) 2551 2552 def getHPDDMComplexities(self) -> tuple[float, float]: 2553 """Compute the grid and operator complexities. 2554 2555 Collective. 2556 2557 See Also 2558 -------- 2559 petsc.PCHPDDMGetComplexities 2560 2561 """ 2562 cdef PetscReal gc = 0.0, oc = 0.0 2563 CHKERR(PCHPDDMGetComplexities(self.pc, &gc, &oc)) 2564 return (toReal(gc), toReal(oc)) 2565 2566 def setHPDDMHasNeumannMat(self, has: bool) -> None: 2567 """Set to indicate that the `Mat` passed to the `PC` is the local Neumann matrix. 2568 2569 Logically collective. 2570 2571 Parameters 2572 ---------- 2573 has 2574 Enable to indicate the matrix is the local Neumann matrix. 2575 2576 See Also 2577 -------- 2578 petsc.PCHPDDMHasNeumannMat 2579 2580 """ 2581 cdef PetscBool phas = has 2582 CHKERR(PCHPDDMHasNeumannMat(self.pc, phas)) 2583 2584 def setHPDDMCoarseCorrectionType(self, correction_type: HPDDMCoarseCorrectionType) -> None: 2585 """Set the coarse correction type. 2586 2587 Collective. 2588 2589 Parameters 2590 ---------- 2591 correction_type 2592 The type of coarse correction to apply. 2593 2594 See Also 2595 -------- 2596 petsc.PCHPDDMSetCoarseCorrectionType 2597 2598 """ 2599 cdef PetscPCHPDDMCoarseCorrectionType ctype = correction_type 2600 CHKERR(PCHPDDMSetCoarseCorrectionType(self.pc, ctype)) 2601 2602 def getHPDDMCoarseCorrectionType(self) -> HPDDMCoarseCorrectionType: 2603 """Return the coarse correction type. 2604 2605 Not collective. 2606 2607 See Also 2608 -------- 2609 petsc.PCHPDDMGetCoarseCorrectionType 2610 2611 """ 2612 cdef PetscPCHPDDMCoarseCorrectionType cval = PC_HPDDM_COARSE_CORRECTION_DEFLATED 2613 CHKERR(PCHPDDMGetCoarseCorrectionType(self.pc, &cval)) 2614 return cval 2615 2616 def getHPDDMSTShareSubKSP(self) -> bool: 2617 """Return true if the `KSP` in SLEPc ``ST`` and the subdomain solver is shared. 2618 2619 Not collective. 2620 2621 See Also 2622 -------- 2623 petsc.PCHPDDMGetSTShareSubKSP 2624 2625 """ 2626 cdef PetscBool cval = PETSC_FALSE 2627 CHKERR(PCHPDDMGetSTShareSubKSP(self.pc, &cval)) 2628 return toBool(cval) 2629 2630 def setHPDDMDeflationMat(self, IS uis, Mat U) -> None: 2631 """Set the deflation space used to assemble a coarse operator. 2632 2633 Logically collective. 2634 2635 Parameters 2636 ---------- 2637 uis 2638 The `IS` of the local deflation matrix. 2639 U 2640 The deflation sequential matrix of type `Mat.Type.DENSE`. 2641 2642 See Also 2643 -------- 2644 petsc.PCHPDDMSetDeflationMat 2645 2646 """ 2647 CHKERR(PCHPDDMSetDeflationMat(self.pc, uis.iset, U.mat)) 2648 2649 # --- SPAI --- 2650 2651 def setSPAIEpsilon(self, val: float) -> None: 2652 """Set the tolerance for the preconditioner. 2653 2654 Logically collective. 2655 2656 Parameters 2657 ---------- 2658 val 2659 The tolerance, defaults to ``0.4``. 2660 2661 See Also 2662 -------- 2663 petsc.PCSPAISetEpsilon 2664 2665 """ 2666 cdef PetscReal cval = asReal(val) 2667 CHKERR(PCSPAISetEpsilon(self.pc, cval)) 2668 2669 def setSPAINBSteps(self, nbsteps: int) -> None: 2670 """Set the maximum number of improvement steps per row. 2671 2672 Logically collective. 2673 2674 Parameters 2675 ---------- 2676 nbsteps 2677 The number of steps, defaults to ``5``. 2678 2679 See Also 2680 -------- 2681 petsc.PCSPAISetNBSteps 2682 2683 """ 2684 cdef PetscInt cval = asInt(nbsteps) 2685 CHKERR(PCSPAISetNBSteps(self.pc, cval)) 2686 2687 def setSPAIMax(self, maxval: int) -> None: 2688 """Set the size of working buffers in the preconditioner. 2689 2690 Logically collective. 2691 2692 Parameters 2693 ---------- 2694 maxval 2695 Number of entries in the work arrays to be allocated, defaults to 2696 ``5000``. 2697 2698 See Also 2699 -------- 2700 petsc.PCSPAISetMax 2701 2702 """ 2703 cdef PetscInt cval = asInt(maxval) 2704 CHKERR(PCSPAISetMax(self.pc, cval)) 2705 2706 def setSPAIMaxNew(self, maxval: int) -> None: 2707 """Set the maximum number of new non-zero candidates per step. 2708 2709 Logically collective. 2710 2711 Parameters 2712 ---------- 2713 maxval 2714 Number of entries allowed, defaults to ``5``. 2715 2716 See Also 2717 -------- 2718 petsc.PCSPAISetMaxNew 2719 2720 """ 2721 cdef PetscInt cval = asInt(maxval) 2722 CHKERR(PCSPAISetMaxNew(self.pc, cval)) 2723 2724 def setSPAIBlockSize(self, n: int) -> None: 2725 """Set the block size of the preconditioner. 2726 2727 Logically collective. 2728 2729 Parameters 2730 ---------- 2731 n 2732 The block size, defaults to ``1``. 2733 2734 See Also 2735 -------- 2736 petsc.PCSPAISetBlockSize 2737 2738 """ 2739 cdef PetscInt cval = asInt(n) 2740 CHKERR(PCSPAISetBlockSize(self.pc, cval)) 2741 2742 def setSPAICacheSize(self, size: int) -> None: 2743 """Set the cache size. 2744 2745 Logically collective. 2746 2747 Parameters 2748 ---------- 2749 size 2750 The size of the cache, defaults to ``5``. 2751 2752 See Also 2753 -------- 2754 petsc.PCSPAISetCacheSize 2755 2756 """ 2757 cdef PetscInt cval = asInt(size) 2758 CHKERR(PCSPAISetCacheSize(self.pc, cval)) 2759 2760 def setSPAIVerbose(self, level: int) -> None: 2761 """Set the verbosity level. 2762 2763 Logically collective. 2764 2765 Parameters 2766 ---------- 2767 level 2768 The level of verbosity, defaults to ``1``. 2769 2770 See Also 2771 -------- 2772 petsc.PCSPAISetVerbose 2773 2774 """ 2775 cdef PetscInt cval = asInt(level) 2776 CHKERR(PCSPAISetVerbose(self.pc, cval)) 2777 2778 def setSPAISp(self, sym: int) -> None: 2779 """Set to specify a symmetric sparsity pattern. 2780 2781 Logically collective. 2782 2783 Parameters 2784 ---------- 2785 sym 2786 Enable to indicate the matrix is symmetric. 2787 2788 See Also 2789 -------- 2790 petsc.PCSPAISetSp 2791 2792 """ 2793 cdef PetscInt cval = asInt(sym) 2794 CHKERR(PCSPAISetSp(self.pc, cval)) 2795 2796 # --- DEFLATION --- 2797 2798 def setDeflationInitOnly(self, flg: bool) -> None: 2799 """Set to only perform the initialization. 2800 2801 Logically collective. 2802 2803 Sets initial guess to the solution on the 2804 deflation space but does not apply the deflation preconditioner. The 2805 additional preconditioner is still applied. 2806 2807 Parameters 2808 ---------- 2809 flg 2810 Enable to only initialize the preconditioner. 2811 2812 See Also 2813 -------- 2814 petsc.PCDeflationSetInitOnly 2815 2816 """ 2817 cdef PetscBool cflg = asBool(flg) 2818 CHKERR(PCDeflationSetInitOnly(self.pc, cflg)) 2819 2820 def setDeflationLevels(self, levels: int) -> None: 2821 """Set the maximum level of deflation nesting. 2822 2823 Logically collective. 2824 2825 Parameters 2826 ---------- 2827 levels 2828 The maximum deflation level. 2829 2830 See Also 2831 -------- 2832 petsc.PCDeflationSetLevels 2833 2834 """ 2835 cdef PetscInt clevels = asInt(levels) 2836 CHKERR(PCDeflationSetLevels(self.pc, clevels)) 2837 2838 def setDeflationReductionFactor(self, red: int) -> None: 2839 """Set the reduction factor for the preconditioner. 2840 2841 Logically collective. 2842 2843 Parameters 2844 ---------- 2845 red 2846 The reduction factor or ``DEFAULT``. 2847 2848 See Also 2849 -------- 2850 petsc.PCDeflationSetReductionFactor 2851 2852 """ 2853 cdef PetscInt cred = asInt(red) 2854 CHKERR(PCDeflationSetReductionFactor(self.pc, cred)) 2855 2856 def setDeflationCorrectionFactor(self, fact: float) -> None: 2857 """Set the coarse problem correction factor. 2858 2859 Logically collective. 2860 2861 Parameters 2862 ---------- 2863 fact 2864 The correction factor. 2865 2866 See Also 2867 -------- 2868 petsc.PCDeflationSetCorrectionFactor 2869 2870 """ 2871 cdef PetscScalar cfact = asScalar(fact) 2872 CHKERR(PCDeflationSetCorrectionFactor(self.pc, cfact)) 2873 2874 def setDeflationSpaceToCompute(self, space_type: DeflationSpaceType, size: int) -> None: 2875 """Set the deflation space type. 2876 2877 Logically collective. 2878 2879 Parameters 2880 ---------- 2881 space_type 2882 The deflation space type. 2883 size 2884 The size of the space to compute 2885 2886 See Also 2887 -------- 2888 petsc.PCDeflationSetSpaceToCompute 2889 2890 """ 2891 cdef PetscInt csize = asInt(size) 2892 cdef PetscPCDeflationSpaceType ctype = space_type 2893 CHKERR(PCDeflationSetSpaceToCompute(self.pc, ctype, csize)) 2894 2895 def setDeflationSpace(self, Mat W, transpose: bool) -> None: 2896 """Set the deflation space matrix or its (Hermitian) transpose. 2897 2898 Logically collective. 2899 2900 Parameters 2901 ---------- 2902 W 2903 The deflation matrix. 2904 transpose 2905 Enable to indicate that ``W`` is an explicit transpose of the 2906 deflation matrix. 2907 2908 See Also 2909 -------- 2910 petsc.PCDeflationSetSpace 2911 2912 """ 2913 cdef PetscBool ctranspose = asBool(transpose) 2914 CHKERR(PCDeflationSetSpace(self.pc, W.mat, ctranspose)) 2915 2916 def setDeflationProjectionNullSpaceMat(self, Mat mat) -> None: 2917 """Set the projection null space matrix. 2918 2919 Collective. 2920 2921 Parameters 2922 ---------- 2923 mat 2924 The projection null space matrix. 2925 2926 See Also 2927 -------- 2928 petsc.PCDeflationSetProjectionNullSpaceMat 2929 2930 """ 2931 CHKERR(PCDeflationSetProjectionNullSpaceMat(self.pc, mat.mat)) 2932 2933 def setDeflationCoarseMat(self, Mat mat) -> None: 2934 """Set the coarse problem matrix. 2935 2936 Collective. 2937 2938 Parameters 2939 ---------- 2940 mat 2941 The coarse problem matrix. 2942 2943 See Also 2944 -------- 2945 petsc.PCDeflationSetCoarseMat 2946 2947 """ 2948 CHKERR(PCDeflationSetCoarseMat(self.pc, mat.mat)) 2949 2950 def getDeflationCoarseKSP(self) -> KSP: 2951 """Return the coarse problem `KSP`. 2952 2953 Not collective. 2954 2955 See Also 2956 -------- 2957 petsc.PCDeflationGetCoarseKSP 2958 2959 """ 2960 cdef KSP ksp = KSP() 2961 CHKERR(PCDeflationGetCoarseKSP(self.pc, &ksp.ksp)) 2962 CHKERR(PetscINCREF(ksp.obj)) 2963 return ksp 2964 2965 def getDeflationPC(self) -> PC: 2966 """Return the additional preconditioner. 2967 2968 Not collective. 2969 2970 See Also 2971 -------- 2972 petsc.PCDeflationGetPC 2973 2974 """ 2975 cdef PC apc = PC() 2976 CHKERR(PCDeflationGetPC(self.pc, &apc.pc)) 2977 CHKERR(PetscINCREF(apc.obj)) 2978 return apc 2979 2980# -------------------------------------------------------------------- 2981 2982del PCType 2983del PCSide 2984del PCASMType 2985del PCGASMType 2986del PCMGType 2987del PCMGCycleType 2988del PCGAMGType 2989del PCCompositeType 2990del PCFieldSplitSchurPreType 2991del PCFieldSplitSchurFactType 2992del PCPatchConstructType 2993del PCHPDDMCoarseCorrectionType 2994del PCDeflationSpaceType 2995del PCFailedReason 2996 2997# -------------------------------------------------------------------- 2998