Actual source code: petscpctypes.h

  1: #pragma once

  3: /* SUBMANSEC = PC */

  5: /*S
  6:    PC - Abstract PETSc object that manages all preconditioners including direct solvers such as `PCLU`

  8:    Level: beginner

 10: .seealso: [](doc_linsolve), [](sec_pc), `PCCreate()`, `PCSetType()`, `PCType`
 11: S*/
 12: typedef struct _p_PC *PC;

 14: /*J
 15:    PCType - String with the name of a PETSc preconditioner

 17:    Level: beginner

 19:    Note:
 20:    `PCRegister()` is used to register preconditioners that are then accessible via `PCSetType()`

 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"

 80: /*E
 81:     PCSide - If the preconditioner is to be applied to the left, right
 82:              or symmetrically around the operator.

 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.

 89:    Level: beginner

 91:    Note:
 92:    Certain `KSPType` support only a subset of `PCSide` values

 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)

104: /*E
105:     PCRichardsonConvergedReason - reason a `PCRICHARDSON` `PCApplyRichardson()` method terminated

107:    Level: advanced

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;

118: /*E
119:     PCJacobiType - What elements of the matrix are used to form the Jacobi preconditioner

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)

126:    Level: intermediate

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;

136: /*E
137:     PCASMType - Type of additive Schwarz method to use

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.

150:    Level: beginner

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;

161: /*E
162:     PCGASMType - Type of generalized additive Schwarz method to use (differs from `PCASM` in allowing multiple processors per subdomain).

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.

177:    Level: beginner

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.

185:    Developer Note:
186:    Perhaps better to remove this since it matches `PCASMType`

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;

197: /*E
198:     PCCompositeType - Determines how two or more preconditioner are composed with the `PCType` of `PCCOMPOSITE`

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`

211:    Level: beginner

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;

224: /*E
225:     PCFieldSplitSchurPreType - Determines how to precondition a Schur complement

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

240:     Level: intermediate

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;

252: /*E
253:     PCFieldSplitSchurFactType - determines which off-diagonal parts of the approximate block factorization to use

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)`

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

267:     Level: intermediate

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;

278: /*E
279:     PCPARMSGlobalType - Determines the global preconditioner method in `PCPARMS`

281:     Level: intermediate

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;

291: /*E
292:     PCPARMSLocalType - Determines the local preconditioner method in `PCPARMS`

294:     Level: intermediate

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;

305: /*J
306:     PCGAMGType - type of generalized algebraic multigrid `PCGAMG` method

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

313:      Level: intermediate

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"

322: typedef const char *PCGAMGClassicalType;
323: #define PCGAMGCLASSICALDIRECT   "direct"
324: #define PCGAMGCLASSICALSTANDARD "standard"

326: /*E
327:    PCMGType - Determines the type of multigrid method that is run.

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

341:    Level: beginner

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;

353: /*E
354:    PCMGCycleType - Use V-cycle or W-cycle

356:    Values:
357: +  `PC_MG_V_CYCLE` - use the V cycle
358: -  `PC_MG_W_CYCLE` - use the W cycle

360:    Level: beginner

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;

369: /*E
370:     PCMGalerkinType - Determines if the coarse grid operators are computed via the Galerkin process

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

378:    Level: beginner

380:    Note:
381:    Users should never set `PC_MG_GALERKIN_EXTERNAL`, it is used by `PCHYPRE` and `PCML`

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;

393: /*E
394:     PCExoticType - Face based or wirebasket based coarse grid space

396:    Level: beginner

398: .seealso: [](sec_pc), `PCExoticSetType()`, `PCEXOTIC`
399: E*/
400: typedef enum {
401:   PC_EXOTIC_FACE,
402:   PC_EXOTIC_WIREBASKET
403: } PCExoticType;

405: /*E
406:    PCBDDCInterfaceExtType - Defines how interface balancing is extended into the interior of subdomains.

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"

412:    Level: intermediate

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;

421: /*E
422:   PCMGCoarseSpaceType - Function space for coarse space for adaptive interpolation

424:   Level: beginner

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;

437: /*E
438:     PCPatchConstructType - The algorithm used to construct patches for the `PCPATCH` preconditioner

440:    Level: beginner

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;

452: /*E
453:     PCDeflationSpaceType - Type of deflation

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

466:     Level: intermediate

468:     Note:
469:     Wavelet-based space (except Haar) can be used in multilevel deflation.

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;

485: /*E
486:     PCHPDDMCoarseCorrectionType - Type of coarse correction used by `PCHPDDM`

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)

493:     Level: intermediate

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;

503: /*E
504:     PCHPDDMSchurPreType - Type of `PCHPDDM` preconditioner for a `MATSCHURCOMPLEMENT` generated by `PCFIELDSPLIT` with `PCFieldSplitSchurPreType` set to `PC_FIELDSPLIT_SCHUR_PRE_SELF`

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

510:     Level: advanced

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;

519: /*E
520:     PCFailedReason - indicates type of `PC` failure

522:     Level: beginner

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;

537: /*E
538:     PCGAMGLayoutType - Layout for reduced grids

540:     Level: intermediate

542: .seealso: [](sec_pc), `PCGAMG`, `PC`, `PCGAMGSetCoarseGridLayoutType()`
543: E*/
544: typedef enum {
545:   PCGAMG_LAYOUT_COMPACT,
546:   PCGAMG_LAYOUT_SPREAD
547: } PCGAMGLayoutType;