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