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()` 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 Developer Note: 110 This must match `include/petsc/finclude/petscpc.h` and the `KSPConvergedReason` values in `include/petscksp.h 111 112 .seealso: [](sec_pc), `PCRICHARDSON`, `PC`, `PCApplyRichardson()` 113 E*/ 114 typedef enum { 115 PCRICHARDSON_CONVERGED_RTOL = 2, 116 PCRICHARDSON_CONVERGED_ATOL = 3, 117 PCRICHARDSON_CONVERGED_ITS = 4, 118 PCRICHARDSON_DIVERGED_DTOL = -4 119 } PCRichardsonConvergedReason; 120 121 /*E 122 PCJacobiType - What elements of the matrix are used to form the Jacobi preconditioner 123 124 Values: 125 + `PC_JACOBI_DIAGONAL` - use the diagonal entry, if it is zero use one 126 . `PC_JACOBI_ROWMAX` - use the maximum absolute value in the row 127 - `PC_JACOBI_ROWSUM` - use the sum of the values in the row (not the absolute values) 128 129 Level: intermediate 130 131 .seealso: [](sec_pc), `PCJACOBI`, `PC` 132 E*/ 133 typedef enum { 134 PC_JACOBI_DIAGONAL, 135 PC_JACOBI_ROWMAX, 136 PC_JACOBI_ROWSUM 137 } PCJacobiType; 138 139 /*E 140 PCASMType - Type of additive Schwarz method to use 141 142 Values: 143 + `PC_ASM_BASIC` - Symmetric version where residuals from the ghost points are used 144 and computed values in ghost regions are added together. 145 Classical standard additive Schwarz. 146 . `PC_ASM_RESTRICT` - Residuals from ghost points are used but computed values in ghost 147 region are discarded. 148 Default. 149 . `PC_ASM_INTERPOLATE` - Residuals from ghost points are not used, computed values in ghost 150 region are added back in. 151 - `PC_ASM_NONE` - Residuals from ghost points are not used, computed ghost values are 152 discarded. 153 Not very good. 154 155 Level: beginner 156 157 .seealso: [](sec_pc), `PC`, `PCASM`, `PCASMSetType()`, `PCGASMType` 158 E*/ 159 typedef enum { 160 PC_ASM_BASIC = 3, 161 PC_ASM_RESTRICT = 1, 162 PC_ASM_INTERPOLATE = 2, 163 PC_ASM_NONE = 0 164 } PCASMType; 165 166 /*E 167 PCGASMType - Type of generalized additive Schwarz method to use (differs from `PCASM` in allowing multiple processors per subdomain). 168 169 Values: 170 + `PC_GASM_BASIC` - Symmetric version where the full from the outer subdomain is used, and the resulting correction is applied 171 over the outer subdomains. As a result, points in the overlap will receive the sum of the corrections 172 from neighboring subdomains. 173 Classical standard additive Schwarz. 174 . `PC_GASM_RESTRICT` - Residual from the outer subdomain is used but the correction is restricted to the inner subdomain only 175 (i.e., zeroed out over the overlap portion of the outer subdomain before being applied). As a result, 176 each point will receive a correction only from the unique inner subdomain containing it (nonoverlapping covering 177 assumption). 178 Default. 179 . `PC_GASM_INTERPOLATE` - Residual is zeroed out over the overlap portion of the outer subdomain, but the resulting correction is 180 applied over the outer subdomain. As a result, points in the overlap will receive the sum of the corrections 181 from neighboring subdomains. 182 - `PC_GASM_NONE` - Residuals and corrections are zeroed out outside the local subdomains. 183 Not very good. 184 185 Level: beginner 186 187 Note: 188 Each subdomain has nested inner and outer parts. The inner subdomains are assumed to form a non-overlapping covering of the computational 189 domain, while the outer subdomains contain the inner subdomains and overlap with each other. This preconditioner will compute 190 a subdomain correction over each *outer* subdomain from a residual computed there, but its different variants will differ in 191 (a) how the outer subdomain residual is computed, and (b) how the outer subdomain correction is computed. 192 193 .seealso: [](sec_pc), `PCGASM`, `PCASM`, `PC`, `PCGASMSetType()`, `PCASMType` 194 E*/ 195 typedef enum { 196 PC_GASM_BASIC = 3, 197 PC_GASM_RESTRICT = 1, 198 PC_GASM_INTERPOLATE = 2, 199 PC_GASM_NONE = 0 200 } PCGASMType; 201 202 /*E 203 PCCompositeType - Determines how two or more preconditioner are composed with the `PCType` of `PCCOMPOSITE` 204 205 Values: 206 + `PC_COMPOSITE_ADDITIVE` - results from application of all preconditioners are added together 207 . `PC_COMPOSITE_MULTIPLICATIVE` - preconditioners are applied sequentially to the residual freshly 208 computed after the previous preconditioner application 209 . `PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE` - preconditioners are applied sequentially to the residual freshly 210 computed from first preconditioner to last and then back (Use only for symmetric matrices and preconditioners) 211 . `PC_COMPOSITE_SPECIAL` - This is very special for a matrix of the form alpha I + R + S 212 where first preconditioner is built from alpha I + S and second from 213 alpha I + R 214 . `PC_COMPOSITE_SCHUR` - composes the Schur complement of the matrix from two blocks, see `PCFIELDSPLIT` 215 - `PC_COMPOSITE_GKB` - the generalized Golub-Kahan bidiagonalization preconditioner, see `PCFIELDSPLIT` 216 217 Level: beginner 218 219 .seealso: [](sec_pc), `PCCOMPOSITE`, `PCFIELDSPLIT`, `PC`, `PCCompositeSetType()` 220 E*/ 221 typedef enum { 222 PC_COMPOSITE_ADDITIVE, 223 PC_COMPOSITE_MULTIPLICATIVE, 224 PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE, 225 PC_COMPOSITE_SPECIAL, 226 PC_COMPOSITE_SCHUR, 227 PC_COMPOSITE_GKB 228 } PCCompositeType; 229 230 /*E 231 PCFieldSplitSchurPreType - Determines how to precondition a Schur complement 232 233 Values: 234 + `PC_FIELDSPLIT_SCHUR_PRE_SELF` - the preconditioner for the Schur complement is generated from the symbolic representation of the Schur complement matrix. 235 The only preconditioners that currently work with this symbolic representation matrix object are `PCLSC` and `PCHPDDM` 236 . `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. 237 This is only a good preconditioner when diag(A00) is a good preconditioner for A00. Optionally, A00 can be 238 lumped before extracting the diagonal using the additional option `-fieldsplit_1_mat_schur_complement_ainv_type lump` 239 . `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, 240 associated with the Schur complement (i.e. A11), not the Schur complement matrix 241 . `PC_FIELDSPLIT_SCHUR_PRE_USER` - the preconditioner for the Schur complement is generated from the user provided matrix (pre argument 242 to this function). 243 - `PC_FIELDSPLIT_SCHUR_PRE_FULL` - the preconditioner for the Schur complement is generated from the exact Schur complement matrix representation 244 computed internally by `PCFIELDSPLIT` (this is expensive) useful mostly as a test that the Schur complement approach can work for your problem 245 246 Level: intermediate 247 248 .seealso: [](sec_pc), `PCFIELDSPLIT`, `PCFieldSplitSetSchurPre()`, `PC` 249 E*/ 250 typedef enum { 251 PC_FIELDSPLIT_SCHUR_PRE_SELF, 252 PC_FIELDSPLIT_SCHUR_PRE_SELFP, 253 PC_FIELDSPLIT_SCHUR_PRE_A11, 254 PC_FIELDSPLIT_SCHUR_PRE_USER, 255 PC_FIELDSPLIT_SCHUR_PRE_FULL 256 } PCFieldSplitSchurPreType; 257 258 /*E 259 PCFieldSplitSchurFactType - determines which off-diagonal parts of the approximate block factorization to use 260 261 Values: 262 + `PC_FIELDSPLIT_SCHUR_FACT_DIAG` - the preconditioner is solving `D` 263 . `PC_FIELDSPLIT_SCHUR_FACT_LOWER` - the preconditioner is solving `L D` 264 . `PC_FIELDSPLIT_SCHUR_FACT_UPPER` - the preconditioner is solving `D U` 265 - `PC_FIELDSPLIT_SCHUR_FACT_FULL` - the preconditioner is solving `L(D U)` 266 267 where the matrix is factorized as 268 .vb 269 (A B) = (1 0) (A 0) (1 Ainv*B) = L D U 270 (C E) (C*Ainv 1) (0 S) (0 1) 271 .ve 272 273 Level: intermediate 274 275 .seealso: [](sec_pc), `PCFIELDSPLIT`, `PCFieldSplitSetSchurFactType()`, `PC` 276 E*/ 277 typedef enum { 278 PC_FIELDSPLIT_SCHUR_FACT_DIAG, 279 PC_FIELDSPLIT_SCHUR_FACT_LOWER, 280 PC_FIELDSPLIT_SCHUR_FACT_UPPER, 281 PC_FIELDSPLIT_SCHUR_FACT_FULL 282 } PCFieldSplitSchurFactType; 283 284 /*E 285 PCPARMSGlobalType - Determines the global preconditioner method in `PCPARMS` 286 287 Level: intermediate 288 289 .seealso: [](sec_pc), `PCPARMS`, `PCPARMSSetGlobal()`, `PC` 290 E*/ 291 typedef enum { 292 PC_PARMS_GLOBAL_RAS, 293 PC_PARMS_GLOBAL_SCHUR, 294 PC_PARMS_GLOBAL_BJ 295 } PCPARMSGlobalType; 296 297 /*E 298 PCPARMSLocalType - Determines the local preconditioner method in `PCPARMS` 299 300 Level: intermediate 301 302 .seealso: [](sec_pc), `PCPARMS`, `PCPARMSSetLocal()`, `PC` 303 E*/ 304 typedef enum { 305 PC_PARMS_LOCAL_ILU0, 306 PC_PARMS_LOCAL_ILUK, 307 PC_PARMS_LOCAL_ILUT, 308 PC_PARMS_LOCAL_ARMS 309 } PCPARMSLocalType; 310 311 /*J 312 PCGAMGType - type of generalized algebraic multigrid `PCGAMG` method 313 314 Values: 315 + `PCGAMGAGG` - (the default) smoothed aggregation algorithm, robust, very well tested 316 . `PCGAMGGEO` - geometric coarsening, uses mesh generator to produce coarser meshes, limited to triangles, not supported, reference implementation (2D) 317 - `PCGAMGCLASSICAL` - classical algebraic multigrid preconditioner, incomplete, not supported, reference implementation 318 319 Level: intermediate 320 321 .seealso: [](sec_pc), `PCGAMG`, `PCMG`, `PC`, `PCSetType()`, `PCGAMGSetThreshold()`, `PCGAMGSetThreshold()`, `PCGAMGSetReuseInterpolation()` 322 J*/ 323 typedef const char *PCGAMGType; 324 #define PCGAMGAGG "agg" 325 #define PCGAMGGEO "geo" 326 #define PCGAMGCLASSICAL "classical" 327 328 typedef const char *PCGAMGClassicalType; 329 #define PCGAMGCLASSICALDIRECT "direct" 330 #define PCGAMGCLASSICALSTANDARD "standard" 331 332 /*E 333 PCMGType - Determines the type of multigrid method that is run. 334 335 Values: 336 + `PC_MG_MULTIPLICATIVE` (default) - traditional V or W cycle as determined by `PCMGSetCycleType()` 337 . `PC_MG_ADDITIVE` - the additive multigrid preconditioner where all levels are 338 smoothed before updating the residual. This only uses the 339 down smoother, in the preconditioner the upper smoother is ignored 340 . `PC_MG_FULL` - same as multiplicative except one also performs grid sequencing, 341 that is starts on the coarsest grid, performs a cycle, interpolates 342 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 343 algorithm supports smoothing on before the restriction on each level in the initial restriction to the coarsest stage. In addition that algorithm 344 calls the V-cycle only on the coarser level and has a post-smoother instead. 345 - `PC_MG_KASKADE` - like full multigrid except one never goes back to a coarser level from a finer 346 347 Level: beginner 348 349 .seealso: [](sec_pc), `PCMG`, `PC`, `PCMGSetType()`, `PCMGSetCycleType()`, `PCMGSetCycleTypeOnLevel()` 350 E*/ 351 typedef enum { 352 PC_MG_MULTIPLICATIVE, 353 PC_MG_ADDITIVE, 354 PC_MG_FULL, 355 PC_MG_KASKADE 356 } PCMGType; 357 #define PC_MG_CASCADE PC_MG_KASKADE; 358 359 /*E 360 PCMGCycleType - Use V-cycle or W-cycle 361 362 Values: 363 + `PC_MG_V_CYCLE` - use the v cycle 364 - `PC_MG_W_CYCLE` - use the w cycle 365 366 Level: beginner 367 368 .seealso: [](sec_pc), `PCMG`, `PC`, `PCMGSetCycleType()` 369 E*/ 370 typedef enum { 371 PC_MG_CYCLE_V = 1, 372 PC_MG_CYCLE_W = 2 373 } PCMGCycleType; 374 375 /*E 376 PCMGalerkinType - Determines if the coarse grid operators are computed via the Galerkin process 377 378 Values: 379 + `PC_MG_GALERKIN_PMAT` - computes the pmat (matrix from which the preconditioner is built) via the Galerkin process from the finest grid 380 . `PC_MG_GALERKIN_MAT` - computes the mat (matrix used to apply the operator) via the Galerkin process from the finest grid 381 . `PC_MG_GALERKIN_BOTH` - computes both the mat and pmat via the Galerkin process (if pmat == mat the construction is only done once 382 - `PC_MG_GALERKIN_NONE` - neither operator is computed via the Galerkin process, the user must provide the operator 383 384 Level: beginner 385 386 Note: 387 Users should never set `PC_MG_GALERKIN_EXTERNAL`, it is used by `PCHYPRE` and `PCML` 388 389 .seealso: [](sec_pc), `PCMG`, `PC`, `PCMGSetCycleType()` 390 E*/ 391 typedef enum { 392 PC_MG_GALERKIN_BOTH, 393 PC_MG_GALERKIN_PMAT, 394 PC_MG_GALERKIN_MAT, 395 PC_MG_GALERKIN_NONE, 396 PC_MG_GALERKIN_EXTERNAL 397 } PCMGGalerkinType; 398 399 /*E 400 PCExoticType - Face based or wirebasket based coarse grid space 401 402 Level: beginner 403 404 .seealso: [](sec_pc), `PCExoticSetType()`, `PCEXOTIC` 405 E*/ 406 typedef enum { 407 PC_EXOTIC_FACE, 408 PC_EXOTIC_WIREBASKET 409 } PCExoticType; 410 411 /*E 412 PCBDDCInterfaceExtType - Defines how interface balancing is extended into the interior of subdomains. 413 414 Values: 415 + `PC_BDDC_INTERFACE_EXT_DIRICHLET` - solves Dirichlet interior problem; this is the standard BDDC algorithm 416 - `PC_BDDC_INTERFACE_EXT_LUMP` - skips interior solve; sometimes called M_1 and associated with "lumped FETI-DP" 417 418 Level: intermediate 419 420 .seealso: [](sec_pc), `PCBDDC`, `PC` 421 E*/ 422 typedef enum { 423 PC_BDDC_INTERFACE_EXT_DIRICHLET, 424 PC_BDDC_INTERFACE_EXT_LUMP 425 } PCBDDCInterfaceExtType; 426 427 /*E 428 PCMGCoarseSpaceType - Function space for coarse space for adaptive interpolation 429 430 Level: beginner 431 432 .seealso: [](sec_pc), `PCMGSetAdaptCoarseSpaceType()`, `PCMG`, `PC` 433 E*/ 434 typedef enum { 435 PCMG_ADAPT_NONE, 436 PCMG_ADAPT_POLYNOMIAL, 437 PCMG_ADAPT_HARMONIC, 438 PCMG_ADAPT_EIGENVECTOR, 439 PCMG_ADAPT_GENERALIZED_EIGENVECTOR, 440 PCMG_ADAPT_GDSW 441 } PCMGCoarseSpaceType; 442 443 /*E 444 PCPatchConstructType - The algorithm used to construct patches for the `PCPATCH` preconditioner 445 446 Level: beginner 447 448 .seealso: [](sec_pc), `PCPatchSetConstructType()`, `PCPATCH`, `PC` 449 E*/ 450 typedef enum { 451 PC_PATCH_STAR, 452 PC_PATCH_VANKA, 453 PC_PATCH_PARDECOMP, 454 PC_PATCH_USER, 455 PC_PATCH_PYTHON 456 } PCPatchConstructType; 457 458 /*E 459 PCDeflationSpaceType - Type of deflation 460 461 Values: 462 + `PC_DEFLATION_SPACE_HAAR` - directly assembled based on Haar (db2) wavelet with overflowed filter cuted-off 463 . `PC_DEFLATION_SPACE_DB2` - `MATCOMPOSITE` of 1-lvl matices based on db2 (2 coefficient Daubechies / Haar wavelet) 464 . `PC_DEFLATION_SPACE_DB4` - same as above, but with db4 (4 coefficient Daubechies) 465 . `PC_DEFLATION_SPACE_DB8` - same as above, but with db8 (8 coefficient Daubechies) 466 . `PC_DEFLATION_SPACE_DB16` - same as above, but with db16 (16 coefficient Daubechies) 467 . `PC_DEFLATION_SPACE_BIORTH22` - same as above, but with biorthogonal 2.2 (6 coefficients) 468 . `PC_DEFLATION_SPACE_MEYER` - same as above, but with Meyer/FIR (62 coefficients) 469 . `PC_DEFLATION_SPACE_AGGREGATION` - aggregates local indices (given by operator matrix distribution) into a subdomain 470 - `PC_DEFLATION_SPACE_USER` - indicates space set by user 471 472 Level: intermediate 473 474 Note: 475 Wavelet-based space (except Haar) can be used in multilevel deflation. 476 477 .seealso: [](sec_pc), `PCDeflationSetSpaceToCompute()`, `PCDEFLATION`, `PC` 478 E*/ 479 typedef enum { 480 PC_DEFLATION_SPACE_HAAR, 481 PC_DEFLATION_SPACE_DB2, 482 PC_DEFLATION_SPACE_DB4, 483 PC_DEFLATION_SPACE_DB8, 484 PC_DEFLATION_SPACE_DB16, 485 PC_DEFLATION_SPACE_BIORTH22, 486 PC_DEFLATION_SPACE_MEYER, 487 PC_DEFLATION_SPACE_AGGREGATION, 488 PC_DEFLATION_SPACE_USER 489 } PCDeflationSpaceType; 490 491 /*E 492 PCHPDDMCoarseCorrectionType - Type of coarse correction used by `PCHPDDM` 493 494 Values: 495 + `PC_HPDDM_COARSE_CORRECTION_DEFLATED` (default) - eq. (1) in `PCHPDDMShellApply()` 496 . `PC_HPDDM_COARSE_CORRECTION_ADDITIVE` - eq. (2) 497 - `PC_HPDDM_COARSE_CORRECTION_BALANCED` - eq. (3) 498 499 Level: intermediate 500 501 .seealso: [](sec_pc), `PCHPDDM`, `PC`, `PCSetType()`, `PCHPDDMShellApply()` 502 E*/ 503 typedef enum { 504 PC_HPDDM_COARSE_CORRECTION_DEFLATED, 505 PC_HPDDM_COARSE_CORRECTION_ADDITIVE, 506 PC_HPDDM_COARSE_CORRECTION_BALANCED 507 } PCHPDDMCoarseCorrectionType; 508 509 /*E 510 PCHPDDMSchurPreType - Type of `PCHPDDM` preconditioner for a `MATSCHURCOMPLEMENT` generated by `PCFIELDSPLIT` with `PCFieldSplitSchurPreType` set to `PC_FIELDSPLIT_SCHUR_PRE_SELF` 511 512 Values: 513 + `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 514 - `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 515 516 Level: advanced 517 518 .seealso: [](sec_pc), `PCHPDDM`, `PC`, `PCFIELDSPLIT`, `PC_FIELDSPLIT_SCHUR_PRE_SELF`, `PCFieldSplitSetSchurPre()`, `PCHPDDMSetAuxiliaryMat()` 519 E*/ 520 typedef enum { 521 PC_HPDDM_SCHUR_PRE_LEAST_SQUARES, 522 PC_HPDDM_SCHUR_PRE_GENEO, 523 } PCHPDDMSchurPreType; 524 525 /*E 526 PCFailedReason - indicates type of `PC` failure 527 528 Level: beginner 529 530 Developer Note: 531 Any additions/changes here MUST also be made in `include/petsc/finclude/petscpc.h` 532 533 .seealso: [](sec_pc), `PC` 534 E*/ 535 typedef enum { 536 PC_SETUP_ERROR = -1, 537 PC_NOERROR, 538 PC_FACTOR_STRUCT_ZEROPIVOT, 539 PC_FACTOR_NUMERIC_ZEROPIVOT, 540 PC_FACTOR_OUTMEMORY, 541 PC_FACTOR_OTHER, 542 PC_INCONSISTENT_RHS, 543 PC_SUBPC_ERROR 544 } PCFailedReason; 545 546 /*E 547 PCGAMGLayoutType - Layout for reduced grids 548 549 Level: intermediate 550 551 Developer Note: 552 Any additions/changes here MUST also be made in `include/petsc/finclude/petscpc.h` 553 554 .seealso: [](sec_pc), `PCGAMG`, `PC`, `PCGAMGSetCoarseGridLayoutType()` 555 E*/ 556 typedef enum { 557 PCGAMG_LAYOUT_COMPACT, 558 PCGAMG_LAYOUT_SPREAD 559 } PCGAMGLayoutType; 560