1 #pragma once 2 3 /* SUBMANSEC = PC */ 4 5 /*S 6 PC - Abstract PETSc object that manages all preconditioners including direct solvers such as `PCLU` 7 8 Level: beginner 9 10 .seealso: [](doc_linsolve), [](sec_pc), `PCCreate()`, `PCSetType()`, `PCType` 11 S*/ 12 typedef struct _p_PC *PC; 13 14 /*J 15 PCType - String with the name of a PETSc preconditioner 16 17 Level: beginner 18 19 Note: 20 `PCRegister()` is used to register preconditioners that are then accessible via `PCSetType()` 21 22 .seealso: [](doc_linsolve), [](sec_pc), `PCSetType()`, `PC`, `PCCreate()`, `PCRegister()`, `PCSetFromOptions()`, `PCLU`, `PCJACOBI`, `PCBJACOBI` 23 J*/ 24 typedef const char *PCType; 25 #define PCNONE "none" 26 #define PCJACOBI "jacobi" 27 #define PCSOR "sor" 28 #define PCLU "lu" 29 #define PCQR "qr" 30 #define PCSHELL "shell" 31 #define PCAMGX "amgx" 32 #define PCBJACOBI "bjacobi" 33 #define PCMG "mg" 34 #define PCEISENSTAT "eisenstat" 35 #define PCILU "ilu" 36 #define PCICC "icc" 37 #define PCASM "asm" 38 #define PCGASM "gasm" 39 #define PCKSP "ksp" 40 #define PCBJKOKKOS "bjkokkos" 41 #define PCCOMPOSITE "composite" 42 #define PCREDUNDANT "redundant" 43 #define PCSPAI "spai" 44 #define PCNN "nn" 45 #define PCCHOLESKY "cholesky" 46 #define PCPBJACOBI "pbjacobi" 47 #define PCVPBJACOBI "vpbjacobi" 48 #define PCMAT "mat" 49 #define PCHYPRE "hypre" 50 #define PCPARMS "parms" 51 #define PCFIELDSPLIT "fieldsplit" 52 #define PCTFS "tfs" 53 #define PCML "ml" 54 #define PCGALERKIN "galerkin" 55 #define PCEXOTIC "exotic" 56 #define PCCP "cp" 57 #define PCBFBT "bfbt" 58 #define PCLSC "lsc" 59 #define PCPYTHON "python" 60 #define PCPFMG "pfmg" 61 #define PCSMG "smg" 62 #define PCSYSPFMG "syspfmg" 63 #define PCREDISTRIBUTE "redistribute" 64 #define PCSVD "svd" 65 #define PCGAMG "gamg" 66 #define PCCHOWILUVIENNACL "chowiluviennacl" 67 #define PCROWSCALINGVIENNACL "rowscalingviennacl" 68 #define PCSAVIENNACL "saviennacl" 69 #define PCBDDC "bddc" 70 #define PCKACZMARZ "kaczmarz" 71 #define PCTELESCOPE "telescope" 72 #define PCPATCH "patch" 73 #define PCLMVM "lmvm" 74 #define PCHMG "hmg" 75 #define PCDEFLATION "deflation" 76 #define PCHPDDM "hpddm" 77 #define PCH2OPUS "h2opus" 78 #define PCMPI "mpi" 79 80 /*E 81 PCSide - If the preconditioner is to be applied to the left, right 82 or symmetrically around the operator. 83 84 Values: 85 + `PC_LEFT` - applied after the operator is applied 86 . `PC_RIGHT` - applied before the operator is applied 87 - `PC_SYMMETRIC` - a portion of the preconditioner is applied before the operator and the transpose of this portion is applied after the operator is applied. 88 89 Level: beginner 90 91 Note: 92 Certain `KSPType` support only a subset of `PCSide` values 93 94 .seealso: [](sec_pc), `PC`, `KSPSetPCSide()`, `KSP`, `KSPType` 95 E*/ 96 typedef enum { 97 PC_SIDE_DEFAULT = -1, 98 PC_LEFT, 99 PC_RIGHT, 100 PC_SYMMETRIC 101 } PCSide; 102 #define PC_SIDE_MAX (PC_SYMMETRIC + 1) 103 104 /*E 105 PCRichardsonConvergedReason - reason a `PCRICHARDSON` `PCApplyRichardson()` method terminated 106 107 Level: advanced 108 109 .seealso: [](sec_pc), `PCRICHARDSON`, `PC`, `PCApplyRichardson()` 110 E*/ 111 typedef enum { 112 PCRICHARDSON_CONVERGED_RTOL = 2, 113 PCRICHARDSON_CONVERGED_ATOL = 3, 114 PCRICHARDSON_CONVERGED_ITS = 4, 115 PCRICHARDSON_DIVERGED_DTOL = -4 116 } PCRichardsonConvergedReason; 117 118 /*E 119 PCJacobiType - What elements of the matrix are used to form the Jacobi preconditioner 120 121 Values: 122 + `PC_JACOBI_DIAGONAL` - use the diagonal entry, if it is zero use one 123 . `PC_JACOBI_ROWMAX` - use the maximum absolute value in the row 124 - `PC_JACOBI_ROWSUM` - use the sum of the values in the row (not the absolute values) 125 126 Level: intermediate 127 128 .seealso: [](sec_pc), `PCJACOBI`, `PC` 129 E*/ 130 typedef enum { 131 PC_JACOBI_DIAGONAL, 132 PC_JACOBI_ROWMAX, 133 PC_JACOBI_ROWSUM 134 } PCJacobiType; 135 136 /*E 137 PCASMType - Type of additive Schwarz method to use 138 139 Values: 140 + `PC_ASM_BASIC` - Symmetric version where residuals from the ghost points are used 141 and computed values in ghost regions are added together. 142 Classical standard additive Schwarz as introduced in {cite}`dryja1987additive`. 143 . `PC_ASM_RESTRICT` - Residuals from ghost points are used but computed values in ghost 144 region are discarded {cite}`cs99`. Default. 145 . `PC_ASM_INTERPOLATE` - Residuals from ghost points are not used, computed values in ghost 146 region are added back in. 147 - `PC_ASM_NONE` - Residuals from ghost points are not used, computed ghost values are 148 discarded. Not very good. 149 150 Level: beginner 151 152 .seealso: [](sec_pc), `PC`, `PCASM`, `PCASMSetType()`, `PCGASMType` 153 E*/ 154 typedef enum { 155 PC_ASM_BASIC = 3, 156 PC_ASM_RESTRICT = 1, 157 PC_ASM_INTERPOLATE = 2, 158 PC_ASM_NONE = 0 159 } PCASMType; 160 161 /*E 162 PCGASMType - Type of generalized additive Schwarz method to use (differs from `PCASM` in allowing multiple processors per subdomain). 163 164 Values: 165 + `PC_GASM_BASIC` - Symmetric version where the full from the outer subdomain is used, and the resulting correction is applied 166 over the outer subdomains. As a result, points in the overlap will receive the sum of the corrections 167 from neighboring subdomains. Classical standard additive Schwarz {cite}`dryja1987additive`. 168 . `PC_GASM_RESTRICT` - Residual from the outer subdomain is used but the correction is restricted to the inner subdomain only 169 (i.e., zeroed out over the overlap portion of the outer subdomain before being applied). As a result, 170 each point will receive a correction only from the unique inner subdomain containing it (nonoverlapping covering 171 assumption) {cite}`cs99`. Default. 172 . `PC_GASM_INTERPOLATE` - Residual is zeroed out over the overlap portion of the outer subdomain, but the resulting correction is 173 applied over the outer subdomain. As a result, points in the overlap will receive the sum of the corrections 174 from neighboring subdomains. 175 - `PC_GASM_NONE` - Residuals and corrections are zeroed out outside the local subdomains. Not very good. 176 177 Level: beginner 178 179 Note: 180 Each subdomain has nested inner and outer parts. The inner subdomains are assumed to form a non-overlapping covering of the computational 181 domain, while the outer subdomains contain the inner subdomains and overlap with each other. This preconditioner will compute 182 a subdomain correction over each *outer* subdomain from a residual computed there, but its different variants will differ in 183 (a) how the outer subdomain residual is computed, and (b) how the outer subdomain correction is computed. 184 185 Developer Note: 186 Perhaps better to remove this since it matches `PCASMType` 187 188 .seealso: [](sec_pc), `PCGASM`, `PCASM`, `PC`, `PCGASMSetType()`, `PCASMType` 189 E*/ 190 typedef enum { 191 PC_GASM_BASIC = 3, 192 PC_GASM_RESTRICT = 1, 193 PC_GASM_INTERPOLATE = 2, 194 PC_GASM_NONE = 0 195 } PCGASMType; 196 197 /*E 198 PCCompositeType - Determines how two or more preconditioner are composed with the `PCType` of `PCCOMPOSITE` 199 200 Values: 201 + `PC_COMPOSITE_ADDITIVE` - results from application of all preconditioners are added together 202 . `PC_COMPOSITE_MULTIPLICATIVE` - preconditioners are applied sequentially to the residual freshly 203 computed after the previous preconditioner application 204 . `PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE` - preconditioners are applied sequentially to the residual freshly 205 computed from first preconditioner to last and then back (Use only for symmetric matrices and preconditioners) 206 . `PC_COMPOSITE_SPECIAL` - This is very special for a matrix of the form $ \alpha I + R + S$ 207 where the first preconditioner is built from $\alpha I + S$ and second from $\alpha I + R$ 208 . `PC_COMPOSITE_SCHUR` - composes the Schur complement of the matrix from two blocks, see `PCFIELDSPLIT` 209 - `PC_COMPOSITE_GKB` - the generalized Golub-Kahan bidiagonalization preconditioner, see `PCFIELDSPLIT` 210 211 Level: beginner 212 213 .seealso: [](sec_pc), `PCCOMPOSITE`, `PCFIELDSPLIT`, `PC`, `PCCompositeSetType()`, `SNESCompositeType`, `PCCompositeSpecialSetAlpha()` 214 E*/ 215 typedef enum { 216 PC_COMPOSITE_ADDITIVE, 217 PC_COMPOSITE_MULTIPLICATIVE, 218 PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE, 219 PC_COMPOSITE_SPECIAL, 220 PC_COMPOSITE_SCHUR, 221 PC_COMPOSITE_GKB 222 } PCCompositeType; 223 224 /*E 225 PCFieldSplitSchurPreType - Determines how to precondition a Schur complement 226 227 Values: 228 + `PC_FIELDSPLIT_SCHUR_PRE_SELF` - the preconditioner for the Schur complement is generated from the symbolic representation of the Schur complement matrix. 229 The only preconditioners that currently work with this symbolic representation matrix object are `PCLSC` and `PCHPDDM` 230 . `PC_FIELDSPLIT_SCHUR_PRE_SELFP` - the preconditioning for the Schur complement is generated from an explicitly-assembled approximation Sp = A11 - A10 inv(diag(A00)) A01. 231 This is only a good preconditioner when diag(A00) is a good preconditioner for A00. Optionally, A00 can be 232 lumped before extracting the diagonal using the additional option `-fieldsplit_1_mat_schur_complement_ainv_type lump` 233 . `PC_FIELDSPLIT_SCHUR_PRE_A11` - the preconditioner for the Schur complement is generated from the block diagonal part of the matrix used to define the preconditioner, 234 associated with the Schur complement (i.e. A11), not the Schur complement matrix 235 . `PC_FIELDSPLIT_SCHUR_PRE_USER` - the preconditioner for the Schur complement is generated from the user provided matrix (pre argument 236 to this function). 237 - `PC_FIELDSPLIT_SCHUR_PRE_FULL` - the preconditioner for the Schur complement is generated from the exact Schur complement matrix representation 238 computed internally by `PCFIELDSPLIT` (this is expensive) useful mostly as a test that the Schur complement approach can work for your problem 239 240 Level: intermediate 241 242 .seealso: [](sec_pc), `PCFIELDSPLIT`, `PCFieldSplitSetSchurPre()`, `PC` 243 E*/ 244 typedef enum { 245 PC_FIELDSPLIT_SCHUR_PRE_SELF, 246 PC_FIELDSPLIT_SCHUR_PRE_SELFP, 247 PC_FIELDSPLIT_SCHUR_PRE_A11, 248 PC_FIELDSPLIT_SCHUR_PRE_USER, 249 PC_FIELDSPLIT_SCHUR_PRE_FULL 250 } PCFieldSplitSchurPreType; 251 252 /*E 253 PCFieldSplitSchurFactType - determines which off-diagonal parts of the approximate block factorization to use 254 255 Values: 256 + `PC_FIELDSPLIT_SCHUR_FACT_DIAG` - the preconditioner is solving `D` 257 . `PC_FIELDSPLIT_SCHUR_FACT_LOWER` - the preconditioner is solving `L D` 258 . `PC_FIELDSPLIT_SCHUR_FACT_UPPER` - the preconditioner is solving `D U` 259 - `PC_FIELDSPLIT_SCHUR_FACT_FULL` - the preconditioner is solving `L(D U)` 260 261 where the matrix is factorized as 262 .vb 263 (A B) = (1 0) (A 0) (1 Ainv*B) = L D U 264 (C E) (C*Ainv 1) (0 S) (0 1) 265 .ve 266 267 Level: intermediate 268 269 .seealso: [](sec_pc), `PCFIELDSPLIT`, `PCFieldSplitSetSchurFactType()`, `PC` 270 E*/ 271 typedef enum { 272 PC_FIELDSPLIT_SCHUR_FACT_DIAG, 273 PC_FIELDSPLIT_SCHUR_FACT_LOWER, 274 PC_FIELDSPLIT_SCHUR_FACT_UPPER, 275 PC_FIELDSPLIT_SCHUR_FACT_FULL 276 } PCFieldSplitSchurFactType; 277 278 /*E 279 PCPARMSGlobalType - Determines the global preconditioner method in `PCPARMS` 280 281 Level: intermediate 282 283 .seealso: [](sec_pc), `PCPARMS`, `PCPARMSSetGlobal()`, `PC` 284 E*/ 285 typedef enum { 286 PC_PARMS_GLOBAL_RAS, 287 PC_PARMS_GLOBAL_SCHUR, 288 PC_PARMS_GLOBAL_BJ 289 } PCPARMSGlobalType; 290 291 /*E 292 PCPARMSLocalType - Determines the local preconditioner method in `PCPARMS` 293 294 Level: intermediate 295 296 .seealso: [](sec_pc), `PCPARMS`, `PCPARMSSetLocal()`, `PC` 297 E*/ 298 typedef enum { 299 PC_PARMS_LOCAL_ILU0, 300 PC_PARMS_LOCAL_ILUK, 301 PC_PARMS_LOCAL_ILUT, 302 PC_PARMS_LOCAL_ARMS 303 } PCPARMSLocalType; 304 305 /*J 306 PCGAMGType - type of generalized algebraic multigrid `PCGAMG` method 307 308 Values: 309 + `PCGAMGAGG` - (the default) smoothed aggregation algorithm, robust, very well tested 310 . `PCGAMGGEO` - geometric coarsening, uses mesh generator to produce coarser meshes, limited to triangles, not supported, reference implementation (2D) 311 - `PCGAMGCLASSICAL` - classical algebraic multigrid preconditioner, incomplete, not supported, reference implementation 312 313 Level: intermediate 314 315 .seealso: [](sec_pc), `PCGAMG`, `PCMG`, `PC`, `PCSetType()`, `PCGAMGSetThreshold()`, `PCGAMGSetThreshold()`, `PCGAMGSetReuseInterpolation()` 316 J*/ 317 typedef const char *PCGAMGType; 318 #define PCGAMGAGG "agg" 319 #define PCGAMGGEO "geo" 320 #define PCGAMGCLASSICAL "classical" 321 322 typedef const char *PCGAMGClassicalType; 323 #define PCGAMGCLASSICALDIRECT "direct" 324 #define PCGAMGCLASSICALSTANDARD "standard" 325 326 /*E 327 PCMGType - Determines the type of multigrid method that is run. 328 329 Values: 330 + `PC_MG_MULTIPLICATIVE` (default) - traditional V or W cycle as determined by `PCMGSetCycleType()` 331 . `PC_MG_ADDITIVE` - the additive multigrid preconditioner where all levels are 332 smoothed before updating the residual. This only uses the 333 down smoother, in the preconditioner the upper smoother is ignored 334 . `PC_MG_FULL` - same as multiplicative except one also performs grid sequencing, 335 that is starts on the coarsest grid, performs a cycle, interpolates 336 to the next, performs a cycle etc. This is much like the F-cycle presented in "Multigrid" by Trottenberg, Oosterlee, Schuller page 49, but that 337 algorithm supports smoothing on before the restriction on each level in the initial restriction to the coarsest stage. In addition that algorithm 338 calls the V-cycle only on the coarser level and has a post-smoother instead. 339 - `PC_MG_KASKADE` - like full multigrid except one never goes back to a coarser level from a finer 340 341 Level: beginner 342 343 .seealso: [](sec_pc), `PCMG`, `PC`, `PCMGSetType()`, `PCMGSetCycleType()`, `PCMGSetCycleTypeOnLevel()` 344 E*/ 345 typedef enum { 346 PC_MG_MULTIPLICATIVE, 347 PC_MG_ADDITIVE, 348 PC_MG_FULL, 349 PC_MG_KASKADE 350 } PCMGType; 351 #define PC_MG_CASCADE PC_MG_KASKADE; 352 353 /*E 354 PCMGCycleType - Use V-cycle or W-cycle 355 356 Values: 357 + `PC_MG_V_CYCLE` - use the V cycle 358 - `PC_MG_W_CYCLE` - use the W cycle 359 360 Level: beginner 361 362 .seealso: [](sec_pc), `PCMG`, `PC`, `PCMGSetCycleType()` 363 E*/ 364 typedef enum { 365 PC_MG_CYCLE_V = 1, 366 PC_MG_CYCLE_W = 2 367 } PCMGCycleType; 368 369 /*E 370 PCMGalerkinType - Determines if the coarse grid operators are computed via the Galerkin process 371 372 Values: 373 + `PC_MG_GALERKIN_PMAT` - computes the pmat (matrix from which the preconditioner is built) via the Galerkin process from the finest grid 374 . `PC_MG_GALERKIN_MAT` - computes the mat (matrix used to apply the operator) via the Galerkin process from the finest grid 375 . `PC_MG_GALERKIN_BOTH` - computes both the mat and pmat via the Galerkin process (if pmat == mat the construction is only done once 376 - `PC_MG_GALERKIN_NONE` - neither operator is computed via the Galerkin process, the user must provide the operator 377 378 Level: beginner 379 380 Note: 381 Users should never set `PC_MG_GALERKIN_EXTERNAL`, it is used by `PCHYPRE` and `PCML` 382 383 .seealso: [](sec_pc), `PCMG`, `PC`, `PCMGSetCycleType()` 384 E*/ 385 typedef enum { 386 PC_MG_GALERKIN_BOTH, 387 PC_MG_GALERKIN_PMAT, 388 PC_MG_GALERKIN_MAT, 389 PC_MG_GALERKIN_NONE, 390 PC_MG_GALERKIN_EXTERNAL 391 } PCMGGalerkinType; 392 393 /*E 394 PCExoticType - Face based or wirebasket based coarse grid space 395 396 Level: beginner 397 398 .seealso: [](sec_pc), `PCExoticSetType()`, `PCEXOTIC` 399 E*/ 400 typedef enum { 401 PC_EXOTIC_FACE, 402 PC_EXOTIC_WIREBASKET 403 } PCExoticType; 404 405 /*E 406 PCBDDCInterfaceExtType - Defines how interface balancing is extended into the interior of subdomains. 407 408 Values: 409 + `PC_BDDC_INTERFACE_EXT_DIRICHLET` - solves Dirichlet interior problem; this is the standard BDDC algorithm 410 - `PC_BDDC_INTERFACE_EXT_LUMP` - skips interior solve; sometimes called M_1 and associated with "lumped FETI-DP" 411 412 Level: intermediate 413 414 .seealso: [](sec_pc), `PCBDDC`, `PC` 415 E*/ 416 typedef enum { 417 PC_BDDC_INTERFACE_EXT_DIRICHLET, 418 PC_BDDC_INTERFACE_EXT_LUMP 419 } PCBDDCInterfaceExtType; 420 421 /*E 422 PCMGCoarseSpaceType - Function space for coarse space for adaptive interpolation 423 424 Level: beginner 425 426 .seealso: [](sec_pc), `PCMGSetAdaptCoarseSpaceType()`, `PCMG`, `PC` 427 E*/ 428 typedef enum { 429 PCMG_ADAPT_NONE, 430 PCMG_ADAPT_POLYNOMIAL, 431 PCMG_ADAPT_HARMONIC, 432 PCMG_ADAPT_EIGENVECTOR, 433 PCMG_ADAPT_GENERALIZED_EIGENVECTOR, 434 PCMG_ADAPT_GDSW 435 } PCMGCoarseSpaceType; 436 437 /*E 438 PCPatchConstructType - The algorithm used to construct patches for the `PCPATCH` preconditioner 439 440 Level: beginner 441 442 .seealso: [](sec_pc), `PCPatchSetConstructType()`, `PCPATCH`, `PC` 443 E*/ 444 typedef enum { 445 PC_PATCH_STAR, 446 PC_PATCH_VANKA, 447 PC_PATCH_PARDECOMP, 448 PC_PATCH_USER, 449 PC_PATCH_PYTHON 450 } PCPatchConstructType; 451 452 /*E 453 PCDeflationSpaceType - Type of deflation 454 455 Values: 456 + `PC_DEFLATION_SPACE_HAAR` - directly assembled based on Haar (db2) wavelet with overflowed filter cuted-off 457 . `PC_DEFLATION_SPACE_DB2` - `MATCOMPOSITE` of 1-lvl matices based on db2 (2 coefficient Daubechies / Haar wavelet) 458 . `PC_DEFLATION_SPACE_DB4` - same as above, but with db4 (4 coefficient Daubechies) 459 . `PC_DEFLATION_SPACE_DB8` - same as above, but with db8 (8 coefficient Daubechies) 460 . `PC_DEFLATION_SPACE_DB16` - same as above, but with db16 (16 coefficient Daubechies) 461 . `PC_DEFLATION_SPACE_BIORTH22` - same as above, but with biorthogonal 2.2 (6 coefficients) 462 . `PC_DEFLATION_SPACE_MEYER` - same as above, but with Meyer/FIR (62 coefficients) 463 . `PC_DEFLATION_SPACE_AGGREGATION` - aggregates local indices (given by operator matrix distribution) into a subdomain 464 - `PC_DEFLATION_SPACE_USER` - indicates space set by user 465 466 Level: intermediate 467 468 Note: 469 Wavelet-based space (except Haar) can be used in multilevel deflation. 470 471 .seealso: [](sec_pc), `PCDeflationSetSpaceToCompute()`, `PCDEFLATION`, `PC` 472 E*/ 473 typedef enum { 474 PC_DEFLATION_SPACE_HAAR, 475 PC_DEFLATION_SPACE_DB2, 476 PC_DEFLATION_SPACE_DB4, 477 PC_DEFLATION_SPACE_DB8, 478 PC_DEFLATION_SPACE_DB16, 479 PC_DEFLATION_SPACE_BIORTH22, 480 PC_DEFLATION_SPACE_MEYER, 481 PC_DEFLATION_SPACE_AGGREGATION, 482 PC_DEFLATION_SPACE_USER 483 } PCDeflationSpaceType; 484 485 /*E 486 PCHPDDMCoarseCorrectionType - Type of coarse correction used by `PCHPDDM` 487 488 Values: 489 + `PC_HPDDM_COARSE_CORRECTION_DEFLATED` (default) - eq. (1) in `PCHPDDMShellApply()` 490 . `PC_HPDDM_COARSE_CORRECTION_ADDITIVE` - eq. (2) 491 - `PC_HPDDM_COARSE_CORRECTION_BALANCED` - eq. (3) 492 493 Level: intermediate 494 495 .seealso: [](sec_pc), `PCHPDDM`, `PC`, `PCSetType()`, `PCHPDDMShellApply()` 496 E*/ 497 typedef enum { 498 PC_HPDDM_COARSE_CORRECTION_DEFLATED, 499 PC_HPDDM_COARSE_CORRECTION_ADDITIVE, 500 PC_HPDDM_COARSE_CORRECTION_BALANCED 501 } PCHPDDMCoarseCorrectionType; 502 503 /*E 504 PCHPDDMSchurPreType - Type of `PCHPDDM` preconditioner for a `MATSCHURCOMPLEMENT` generated by `PCFIELDSPLIT` with `PCFieldSplitSchurPreType` set to `PC_FIELDSPLIT_SCHUR_PRE_SELF` 505 506 Values: 507 + `PC_HPDDM_SCHUR_PRE_LEAST_SQUARES` (default) - only with a near-zero A11 block and A10 = A01^T; a preconditioner for solving A01^T A00^-1 A01 x = b is built by approximating the Schur complement with (inv(sqrt(diag(A00))) A01)^T (inv(sqrt(diag(A00))) A01) and by considering the associated linear least squares problem 508 - `PC_HPDDM_SCHUR_PRE_GENEO` - only with A10 = A01^T, `PCHPDDMSetAuxiliaryMat()` called on the `PC` of the A00 block, and if A11 is nonzero, then `PCHPDDMSetAuxiliaryMat()` must be called on the associated `PC` as well (it is built automatically for the user otherwise); the Schur complement `PC` is set internally to `PCKSP`, with the prefix `-fieldsplit_1_pc_hpddm_`; the operator associated to the `PC` is spectrally equivalent to the original Schur complement 509 510 Level: advanced 511 512 .seealso: [](sec_pc), `PCHPDDM`, `PC`, `PCFIELDSPLIT`, `PC_FIELDSPLIT_SCHUR_PRE_SELF`, `PCFieldSplitSetSchurPre()`, `PCHPDDMSetAuxiliaryMat()` 513 E*/ 514 typedef enum { 515 PC_HPDDM_SCHUR_PRE_LEAST_SQUARES, 516 PC_HPDDM_SCHUR_PRE_GENEO, 517 } PCHPDDMSchurPreType; 518 519 /*E 520 PCFailedReason - indicates type of `PC` failure 521 522 Level: beginner 523 524 .seealso: [](sec_pc), `PC` 525 E*/ 526 typedef enum { 527 PC_SETUP_ERROR = -1, 528 PC_NOERROR, 529 PC_FACTOR_STRUCT_ZEROPIVOT, 530 PC_FACTOR_NUMERIC_ZEROPIVOT, 531 PC_FACTOR_OUTMEMORY, 532 PC_FACTOR_OTHER, 533 PC_INCONSISTENT_RHS, 534 PC_SUBPC_ERROR 535 } PCFailedReason; 536 537 /*E 538 PCGAMGLayoutType - Layout for reduced grids 539 540 Level: intermediate 541 542 .seealso: [](sec_pc), `PCGAMG`, `PC`, `PCGAMGSetCoarseGridLayoutType()` 543 E*/ 544 typedef enum { 545 PCGAMG_LAYOUT_COMPACT, 546 PCGAMG_LAYOUT_SPREAD 547 } PCGAMGLayoutType; 548