Actual source code: petscpctypes.h
1: #if !defined(PETSCPCTYPES_H)
2: #define PETSCPCTYPES_H
4: /*S
5: PC - Abstract PETSc object that manages all preconditioners including direct solvers such as PCLU
7: Level: beginner
9: .seealso: PCCreate(), PCSetType(), PCType (for list of available types)
10: S*/
11: typedef struct _p_PC* PC;
13: /*J
14: PCType - String with the name of a PETSc preconditioner method.
16: Level: beginner
18: Notes:
19: Click on the links above to see details on a particular solver
21: PCRegister() is used to register preconditioners that are then accessible via PCSetType()
23: .seealso: PCSetType(), PC, PCCreate(), PCRegister(), PCSetFromOptions()
24: J*/
25: typedef const char* PCType;
26: #define PCNONE "none"
27: #define PCJACOBI "jacobi"
28: #define PCSOR "sor"
29: #define PCLU "lu"
30: #define PCQR "qr"
31: #define PCSHELL "shell"
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 PCSYSPFMG "syspfmg"
62: #define PCREDISTRIBUTE "redistribute"
63: #define PCSVD "svd"
64: #define PCGAMG "gamg"
65: #define PCCHOWILUVIENNACL "chowiluviennacl"
66: #define PCROWSCALINGVIENNACL "rowscalingviennacl"
67: #define PCSAVIENNACL "saviennacl"
68: #define PCBDDC "bddc"
69: #define PCKACZMARZ "kaczmarz"
70: #define PCTELESCOPE "telescope"
71: #define PCPATCH "patch"
72: #define PCLMVM "lmvm"
73: #define PCHMG "hmg"
74: #define PCDEFLATION "deflation"
75: #define PCHPDDM "hpddm"
76: #define PCH2OPUS "h2opus"
78: /*E
79: PCSide - If the preconditioner is to be applied to the left, right
80: or symmetrically around the operator.
82: Level: beginner
84: .seealso:
85: E*/
86: typedef enum { PC_SIDE_DEFAULT=-1,PC_LEFT,PC_RIGHT,PC_SYMMETRIC} PCSide;
87: #define PC_SIDE_MAX (PC_SYMMETRIC + 1)
89: /*E
90: PCRichardsonConvergedReason - reason a PCApplyRichardson method terminates
92: Level: advanced
94: Notes:
95: this must match petsc/finclude/petscpc.h and the KSPConvergedReason values in petscksp.h
97: .seealso: PCApplyRichardson()
98: E*/
99: typedef enum {
100: PCRICHARDSON_CONVERGED_RTOL = 2,
101: PCRICHARDSON_CONVERGED_ATOL = 3,
102: PCRICHARDSON_CONVERGED_ITS = 4,
103: PCRICHARDSON_DIVERGED_DTOL = -4} PCRichardsonConvergedReason;
105: /*E
106: PCJacobiType - What elements are used to form the Jacobi preconditioner
108: Level: intermediate
110: .seealso:
111: E*/
112: typedef enum { PC_JACOBI_DIAGONAL,PC_JACOBI_ROWMAX,PC_JACOBI_ROWSUM} PCJacobiType;
114: /*E
115: PCASMType - Type of additive Schwarz method to use
117: $ PC_ASM_BASIC - Symmetric version where residuals from the ghost points are used
118: $ and computed values in ghost regions are added together.
119: $ Classical standard additive Schwarz.
120: $ PC_ASM_RESTRICT - Residuals from ghost points are used but computed values in ghost
121: $ region are discarded.
122: $ Default.
123: $ PC_ASM_INTERPOLATE - Residuals from ghost points are not used, computed values in ghost
124: $ region are added back in.
125: $ PC_ASM_NONE - Residuals from ghost points are not used, computed ghost values are
126: $ discarded.
127: $ Not very good.
129: Level: beginner
131: .seealso: PCASMSetType()
132: E*/
133: typedef enum {PC_ASM_BASIC = 3,PC_ASM_RESTRICT = 1,PC_ASM_INTERPOLATE = 2,PC_ASM_NONE = 0} PCASMType;
135: /*E
136: PCGASMType - Type of generalized additive Schwarz method to use (differs from ASM in allowing multiple processors per subdomain).
138: Each subdomain has nested inner and outer parts. The inner subdomains are assumed to form a non-overlapping covering of the computational
139: domain, while the outer subdomains contain the inner subdomains and overlap with each other. This preconditioner will compute
140: a subdomain correction over each *outer* subdomain from a residual computed there, but its different variants will differ in
141: (a) how the outer subdomain residual is computed, and (b) how the outer subdomain correction is computed.
143: $ PC_GASM_BASIC - Symmetric version where the full from the outer subdomain is used, and the resulting correction is applied
144: $ over the outer subdomains. As a result, points in the overlap will receive the sum of the corrections
145: $ from neighboring subdomains.
146: $ Classical standard additive Schwarz.
147: $ PC_GASM_RESTRICT - Residual from the outer subdomain is used but the correction is restricted to the inner subdomain only
148: $ (i.e., zeroed out over the overlap portion of the outer subdomain before being applied). As a result,
149: $ each point will receive a correction only from the unique inner subdomain containing it (nonoverlapping covering
150: $ assumption).
151: $ Default.
152: $ PC_GASM_INTERPOLATE - Residual is zeroed out over the overlap portion of the outer subdomain, but the resulting correction is
153: $ applied over the outer subdomain. As a result, points in the overlap will receive the sum of the corrections
154: $ from neighboring subdomains.
155: $
156: $ PC_GASM_NONE - Residuals and corrections are zeroed out outside the local subdomains.
157: $ Not very good.
159: Level: beginner
161: .seealso: PCGASMSetType()
162: E*/
163: typedef enum {PC_GASM_BASIC = 3,PC_GASM_RESTRICT = 1,PC_GASM_INTERPOLATE = 2,PC_GASM_NONE = 0} PCGASMType;
165: /*E
166: PCCompositeType - Determines how two or more preconditioner are composed
168: $ PC_COMPOSITE_ADDITIVE - results from application of all preconditioners are added together
169: $ PC_COMPOSITE_MULTIPLICATIVE - preconditioners are applied sequentially to the residual freshly
170: $ computed after the previous preconditioner application
171: $ PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE - preconditioners are applied sequentially to the residual freshly
172: $ computed from first preconditioner to last and then back (Use only for symmetric matrices and preconditioners)
173: $ PC_COMPOSITE_SPECIAL - This is very special for a matrix of the form alpha I + R + S
174: $ where first preconditioner is built from alpha I + S and second from
175: $ alpha I + R
177: Level: beginner
179: .seealso: PCCompositeSetType()
180: E*/
181: typedef enum {PC_COMPOSITE_ADDITIVE,PC_COMPOSITE_MULTIPLICATIVE,PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE,PC_COMPOSITE_SPECIAL,PC_COMPOSITE_SCHUR,PC_COMPOSITE_GKB} PCCompositeType;
183: /*E
184: PCFieldSplitSchurPreType - Determines how to precondition Schur complement
186: Level: intermediate
188: .seealso: PCFieldSplitSetSchurPre()
189: E*/
190: typedef enum {PC_FIELDSPLIT_SCHUR_PRE_SELF,PC_FIELDSPLIT_SCHUR_PRE_SELFP,PC_FIELDSPLIT_SCHUR_PRE_A11,PC_FIELDSPLIT_SCHUR_PRE_USER,PC_FIELDSPLIT_SCHUR_PRE_FULL} PCFieldSplitSchurPreType;
192: /*E
193: PCFieldSplitSchurFactType - determines which off-diagonal parts of the approximate block factorization to use
195: Level: intermediate
197: .seealso: PCFieldSplitSetSchurFactType()
198: E*/
199: typedef enum {
200: PC_FIELDSPLIT_SCHUR_FACT_DIAG,
201: PC_FIELDSPLIT_SCHUR_FACT_LOWER,
202: PC_FIELDSPLIT_SCHUR_FACT_UPPER,
203: PC_FIELDSPLIT_SCHUR_FACT_FULL
204: } PCFieldSplitSchurFactType;
206: /*E
207: PCPARMSGlobalType - Determines the global preconditioner method in PARMS
209: Level: intermediate
211: .seealso: PCPARMSSetGlobal()
212: E*/
213: typedef enum {PC_PARMS_GLOBAL_RAS,PC_PARMS_GLOBAL_SCHUR,PC_PARMS_GLOBAL_BJ} PCPARMSGlobalType;
215: /*E
216: PCPARMSLocalType - Determines the local preconditioner method in PARMS
218: Level: intermediate
220: .seealso: PCPARMSSetLocal()
221: E*/
222: typedef enum {PC_PARMS_LOCAL_ILU0,PC_PARMS_LOCAL_ILUK,PC_PARMS_LOCAL_ILUT,PC_PARMS_LOCAL_ARMS} PCPARMSLocalType;
224: /*J
225: PCGAMGType - type of generalized algebraic multigrid (PCGAMG) method
227: Level: intermediate
229: $ PCGAMGAGG - (the default) smoothed aggregation algorithm, robust, very well tested
230: $ PCGAMGGEO - geometric coarsening, uses mesh generator to produce coarser meshes, limited to triangles, not well tested
231: $ PCGAMGCLASSICAL - classical algebraic multigrid preconditioner, incomplete, poorly tested
233: .seealso: PCMG, PCSetType(), PCGAMGSetThreshold(), PCGAMGSetThreshold(), PCGAMGSetReuseInterpolation()
234: J*/
235: typedef const char *PCGAMGType;
236: #define PCGAMGAGG "agg"
237: #define PCGAMGGEO "geo"
238: #define PCGAMGCLASSICAL "classical"
240: typedef const char *PCGAMGClassicalType;
241: #define PCGAMGCLASSICALDIRECT "direct"
242: #define PCGAMGCLASSICALSTANDARD "standard"
244: /*E
245: PCMGType - Determines the type of multigrid method that is run.
247: Level: beginner
249: Values:
250: + PC_MG_MULTIPLICATIVE (default) - traditional V or W cycle as determined by PCMGSetCycleType()
251: . PC_MG_ADDITIVE - the additive multigrid preconditioner where all levels are
252: smoothed before updating the residual. This only uses the
253: down smoother, in the preconditioner the upper smoother is ignored
254: . PC_MG_FULL - same as multiplicative except one also performs grid sequencing,
255: that is starts on the coarsest grid, performs a cycle, interpolates
256: 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
257: algorithm supports smoothing on before the restriction on each level in the initial restriction to the coarsest stage. In addition that algorithm
258: calls the V-cycle only on the coarser level and has a post-smoother instead.
259: - PC_MG_KASKADE - like full multigrid except one never goes back to a coarser level
260: from a finer
262: .seealso: PCMGSetType(), PCMGSetCycleType(), PCMGSetCycleTypeOnLevel()
264: E*/
265: typedef enum { PC_MG_MULTIPLICATIVE,PC_MG_ADDITIVE,PC_MG_FULL,PC_MG_KASKADE } PCMGType;
266: #define PC_MG_CASCADE PC_MG_KASKADE;
268: /*E
269: PCMGCycleType - Use V-cycle or W-cycle
271: Level: beginner
273: Values:
274: + PC_MG_V_CYCLE - use the v cycle
275: - PC_MG_W_CYCLE - use the w cycle
277: .seealso: PCMGSetCycleType()
279: E*/
280: typedef enum { PC_MG_CYCLE_V = 1,PC_MG_CYCLE_W = 2 } PCMGCycleType;
282: /*E
283: PCMGalerkinType - Determines if the coarse grid operators are computed via the Galerkin process
285: Level: beginner
287: Values:
288: + PC_MG_GALERKIN_PMAT - computes the pmat (matrix from which the preconditioner is built) via the Galerkin process from the finest grid
289: . PC_MG_GALERKIN_MAT - computes the mat (matrix used to apply the operator) via the Galerkin process from the finest grid
290: . PC_MG_GALERKIN_BOTH - computes both the mat and pmat via the Galerkin process (if pmat == mat the construction is only done once
291: - PC_MG_GALERKIN_NONE - neither operator is computed via the Galerkin process, the user must provide the operator
293: Users should never set PC_MG_GALERKIN_EXTERNAL, it is used by GAMG and ML
295: .seealso: PCMGSetCycleType()
297: E*/
298: typedef enum { PC_MG_GALERKIN_BOTH,PC_MG_GALERKIN_PMAT,PC_MG_GALERKIN_MAT, PC_MG_GALERKIN_NONE, PC_MG_GALERKIN_EXTERNAL} PCMGGalerkinType;
300: /*E
301: PCExoticType - Face based or wirebasket based coarse grid space
303: Level: beginner
305: .seealso: PCExoticSetType(), PCEXOTIC
306: E*/
307: typedef enum { PC_EXOTIC_FACE,PC_EXOTIC_WIREBASKET } PCExoticType;
309: /*E
310: PCBDDCInterfaceExtType - Defines how interface balancing is extended into the interior of subdomains.
312: Level: intermediate
314: Values:
315: + PC_BDDC_INTERFACE_EXT_DIRICHLET - solves Dirichlet interior problem; this is the standard BDDC algorithm
316: - PC_BDDC_INTERFACE_EXT_LUMP - skips interior solve; sometimes called M_1 and associated with "lumped FETI-DP"
318: E*/
319: typedef enum {
320: PC_BDDC_INTERFACE_EXT_DIRICHLET,
321: PC_BDDC_INTERFACE_EXT_LUMP
322: } PCBDDCInterfaceExtType;
324: /*E
325: PCMGCoarseSpaceType - Function space for coarse space for adaptive interpolation
327: Level: beginner
329: .seealso: PCMGSetAdaptCoarseSpaceType(), PCMG
330: E*/
331: typedef enum { PCMG_POLYNOMIAL, PCMG_HARMONIC, PCMG_EIGENVECTOR, PCMG_GENERALIZED_EIGENVECTOR } PCMGCoarseSpaceType;
333: /*E
334: PCPatchConstructType - The algorithm used to construct patches for the preconditioner
336: Level: beginner
338: .seealso: PCPatchSetConstructType(), PCEXOTIC
339: E*/
340: typedef enum {PC_PATCH_STAR, PC_PATCH_VANKA, PC_PATCH_PARDECOMP, PC_PATCH_USER, PC_PATCH_PYTHON} PCPatchConstructType;
342: /*E
343: PCDeflationSpaceType - Type of deflation
345: Values:
346: + PC_DEFLATION_SPACE_HAAR - directly assembled based on Haar (db2) wavelet with overflowed filter cuted-off
347: . PC_DEFLATION_SPACE_DB2 - MATCOMPOSITE of 1-lvl matices based on db2 (2 coefficient Daubechies / Haar wavelet)
348: . PC_DEFLATION_SPACE_DB4 - same as above, but with db4 (4 coefficient Daubechies)
349: . PC_DEFLATION_SPACE_DB8 - same as above, but with db8 (8 coefficient Daubechies)
350: . PC_DEFLATION_SPACE_DB16 - same as above, but with db16 (16 coefficient Daubechies)
351: . PC_DEFLATION_SPACE_BIORTH22 - same as above, but with biorthogonal 2.2 (6 coefficients)
352: . PC_DEFLATION_SPACE_MEYER - same as above, but with Meyer/FIR (62 coefficients)
353: . PC_DEFLATION_SPACE_AGGREGATION - aggregates local indices (given by operator matix distribution) into a subdomain
354: - PC_DEFLATION_SPACE_USER - indicates space set by user
356: Notes:
357: Wavelet-based space (except Haar) can be used in multilevel deflation.
359: Level: intermediate
361: .seealso: PCDeflationSetSpaceToCompute(), PCDEFLATION
362: E*/
363: typedef enum {
364: PC_DEFLATION_SPACE_HAAR,
365: PC_DEFLATION_SPACE_DB2,
366: PC_DEFLATION_SPACE_DB4,
367: PC_DEFLATION_SPACE_DB8,
368: PC_DEFLATION_SPACE_DB16,
369: PC_DEFLATION_SPACE_BIORTH22,
370: PC_DEFLATION_SPACE_MEYER,
371: PC_DEFLATION_SPACE_AGGREGATION,
372: PC_DEFLATION_SPACE_USER
373: } PCDeflationSpaceType;
375: /*E
376: PCHPDDMCoarseCorrectionType - Type of coarse correction used by PCHPDDM
378: Level: intermediate
380: Values:
381: + PC_HPDDM_COARSE_CORRECTION_DEFLATED (default) - eq. (1) in PCHPDDMShellApply()
382: . PC_HPDDM_COARSE_CORRECTION_ADDITIVE - eq. (2)
383: - PC_HPDDM_COARSE_CORRECTION_BALANCED - eq. (3)
385: .seealso: PCHPDDM, PCSetType(), PCHPDDMShellApply()
386: E*/
387: typedef enum { PC_HPDDM_COARSE_CORRECTION_DEFLATED, PC_HPDDM_COARSE_CORRECTION_ADDITIVE, PC_HPDDM_COARSE_CORRECTION_BALANCED } PCHPDDMCoarseCorrectionType;
389: /*E
390: PCFailedReason - indicates type of PC failure
392: Level: beginner
394: Any additions/changes here MUST also be made in include/petsc/finclude/petscpc.h
395: E*/
396: typedef enum {PC_SETUP_ERROR = -1,PC_NOERROR,PC_FACTOR_STRUCT_ZEROPIVOT,PC_FACTOR_NUMERIC_ZEROPIVOT,PC_FACTOR_OUTMEMORY,PC_FACTOR_OTHER,PC_SUBPC_ERROR} PCFailedReason;
398: /*E
399: PCGAMGLayoutType - Layout for reduced grids
401: Level: intermediate
403: .seealso: PCGAMGSetCoarseGridLayoutType()
404: Any additions/changes here MUST also be made in include/petsc/finclude/petscpc.h
405: E*/
406: typedef enum {PCGAMG_LAYOUT_COMPACT,PCGAMG_LAYOUT_SPREAD} PCGAMGLayoutType;
408: #endif