Actual source code: petscpctypes.h
petsc-3.10.5 2019-03-28
1: #if !defined(_PETSCPCTYPES_H)
2: #define _PETSCPCTYPES_H
4: #include <petscdmtypes.h>
6: /*S
7: PC - Abstract PETSc object that manages all preconditioners including direct solvers such as PCLU
9: Level: beginner
11: Concepts: preconditioners
13: .seealso: PCCreate(), PCSetType(), PCType (for list of available types)
14: S*/
15: typedef struct _p_PC* PC;
17: /*J
18: PCType - String with the name of a PETSc preconditioner method.
20: Level: beginner
22: Notes:
23: Click on the links above to see details on a particular solver
25: PCRegister() is used to register preconditioners that are then accessible via PCSetType()
27: .seealso: PCSetType(), PC, PCCreate(), PCRegister(), PCSetFromOptions()
28: J*/
29: typedef const char* PCType;
30: #define PCNONE "none"
31: #define PCJACOBI "jacobi"
32: #define PCSOR "sor"
33: #define PCLU "lu"
34: #define PCSHELL "shell"
35: #define PCBJACOBI "bjacobi"
36: #define PCMG "mg"
37: #define PCEISENSTAT "eisenstat"
38: #define PCILU "ilu"
39: #define PCICC "icc"
40: #define PCASM "asm"
41: #define PCGASM "gasm"
42: #define PCKSP "ksp"
43: #define PCCOMPOSITE "composite"
44: #define PCREDUNDANT "redundant"
45: #define PCSPAI "spai"
46: #define PCNN "nn"
47: #define PCCHOLESKY "cholesky"
48: #define PCPBJACOBI "pbjacobi"
49: #define PCVPBJACOBI "vpbjacobi"
50: #define PCMAT "mat"
51: #define PCHYPRE "hypre"
52: #define PCPARMS "parms"
53: #define PCFIELDSPLIT "fieldsplit"
54: #define PCTFS "tfs"
55: #define PCML "ml"
56: #define PCGALERKIN "galerkin"
57: #define PCEXOTIC "exotic"
58: #define PCCP "cp"
59: #define PCBFBT "bfbt"
60: #define PCLSC "lsc"
61: #define PCPYTHON "python"
62: #define PCPFMG "pfmg"
63: #define PCSYSPFMG "syspfmg"
64: #define PCREDISTRIBUTE "redistribute"
65: #define PCSVD "svd"
66: #define PCGAMG "gamg"
67: #define PCCHOWILUVIENNACL "chowiluviennacl"
68: #define PCROWSCALINGVIENNACL "rowscalingviennacl"
69: #define PCSAVIENNACL "saviennacl"
70: #define PCBDDC "bddc"
71: #define PCKACZMARZ "kaczmarz"
72: #define PCTELESCOPE "telescope"
73: #define PCPATCH "patch"
74: #define PCLMVM "lmvm"
76: /*E
77: PCSide - If the preconditioner is to be applied to the left, right
78: or symmetrically around the operator.
80: Level: beginner
82: .seealso:
83: E*/
84: typedef enum { PC_SIDE_DEFAULT=-1,PC_LEFT,PC_RIGHT,PC_SYMMETRIC} PCSide;
85: #define PC_SIDE_MAX (PC_SYMMETRIC + 1)
86: PETSC_EXTERN const char *const *const PCSides;
88: /*E
89: PCRichardsonConvergedReason - reason a PCApplyRichardson method terminates
91: Level: advanced
93: Notes:
94: this must match petsc/finclude/petscpc.h and the KSPConvergedReason values in petscksp.h
96: .seealso: PCApplyRichardson()
97: E*/
98: typedef enum {
99: PCRICHARDSON_CONVERGED_RTOL = 2,
100: PCRICHARDSON_CONVERGED_ATOL = 3,
101: PCRICHARDSON_CONVERGED_ITS = 4,
102: PCRICHARDSON_DIVERGED_DTOL = -4} PCRichardsonConvergedReason;
104: /*E
105: PCJacobiType - What elements are used to form the Jacobi preconditioner
107: Level: intermediate
109: .seealso:
110: E*/
111: typedef enum { PC_JACOBI_DIAGONAL,PC_JACOBI_ROWMAX,PC_JACOBI_ROWSUM} PCJacobiType;
112: PETSC_EXTERN const char *const PCJacobiTypes[];
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;
134: PETSC_EXTERN const char *const PCASMTypes[];
136: /*E
137: PCGASMType - Type of generalized additive Schwarz method to use (differs from ASM in allowing multiple processors per subdomain).
139: Each subdomain has nested inner and outer parts. The inner subdomains are assumed to form a non-overlapping covering of the computational
140: domain, while the outer subdomains contain the inner subdomains and overlap with each other. This preconditioner will compute
141: a subdomain correction over each *outer* subdomain from a residual computed there, but its different variants will differ in
142: (a) how the outer subdomain residual is computed, and (b) how the outer subdomain correction is computed.
144: $ PC_GASM_BASIC - Symmetric version where the full from the outer subdomain is used, and the resulting correction is applied
145: $ over the outer subdomains. As a result, points in the overlap will receive the sum of the corrections
146: $ from neighboring subdomains.
147: $ Classical standard additive Schwarz.
148: $ PC_GASM_RESTRICT - Residual from the outer subdomain is used but the correction is restricted to the inner subdomain only
149: $ (i.e., zeroed out over the overlap portion of the outer subdomain before being applied). As a result,
150: $ each point will receive a correction only from the unique inner subdomain containing it (nonoverlapping covering
151: $ assumption).
152: $ Default.
153: $ PC_GASM_INTERPOLATE - Residual is zeroed out over the overlap portion of the outer subdomain, but the resulting correction is
154: $ applied over the outer subdomain. As a result, points in the overlap will receive the sum of the corrections
155: $ from neighboring subdomains.
156: $
157: $ PC_GASM_NONE - Residuals and corrections are zeroed out outside the local subdomains.
158: $ Not very good.
160: Level: beginner
162: .seealso: PCGASMSetType()
163: E*/
164: typedef enum {PC_GASM_BASIC = 3,PC_GASM_RESTRICT = 1,PC_GASM_INTERPOLATE = 2,PC_GASM_NONE = 0} PCGASMType;
165: PETSC_EXTERN const char *const PCGASMTypes[];
167: /*E
168: PCCompositeType - Determines how two or more preconditioner are composed
170: $ PC_COMPOSITE_ADDITIVE - results from application of all preconditioners are added together
171: $ PC_COMPOSITE_MULTIPLICATIVE - preconditioners are applied sequentially to the residual freshly
172: $ computed after the previous preconditioner application
173: $ PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE - preconditioners are applied sequentially to the residual freshly
174: $ computed from first preconditioner to last and then back (Use only for symmetric matrices and preconditioners)
175: $ PC_COMPOSITE_SPECIAL - This is very special for a matrix of the form alpha I + R + S
176: $ where first preconditioner is built from alpha I + S and second from
177: $ alpha I + R
179: Level: beginner
181: .seealso: PCCompositeSetType()
182: E*/
183: typedef enum {PC_COMPOSITE_ADDITIVE,PC_COMPOSITE_MULTIPLICATIVE,PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE,PC_COMPOSITE_SPECIAL,PC_COMPOSITE_SCHUR} PCCompositeType;
184: PETSC_EXTERN const char *const PCCompositeTypes[];
186: /*E
187: PCFieldSplitSchurPreType - Determines how to precondition Schur complement
189: Level: intermediate
191: .seealso: PCFieldSplitSetSchurPre()
192: E*/
193: 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;
194: PETSC_EXTERN const char *const PCFieldSplitSchurPreTypes[];
196: /*E
197: PCFieldSplitSchurFactType - determines which off-diagonal parts of the approximate block factorization to use
199: Level: intermediate
201: .seealso: PCFieldSplitSetSchurFactType()
202: E*/
203: typedef enum {
204: PC_FIELDSPLIT_SCHUR_FACT_DIAG,
205: PC_FIELDSPLIT_SCHUR_FACT_LOWER,
206: PC_FIELDSPLIT_SCHUR_FACT_UPPER,
207: PC_FIELDSPLIT_SCHUR_FACT_FULL
208: } PCFieldSplitSchurFactType;
209: PETSC_EXTERN const char *const PCFieldSplitSchurFactTypes[];
211: /*E
212: PCPARMSGlobalType - Determines the global preconditioner method in PARMS
214: Level: intermediate
216: .seealso: PCPARMSSetGlobal()
217: E*/
218: typedef enum {PC_PARMS_GLOBAL_RAS,PC_PARMS_GLOBAL_SCHUR,PC_PARMS_GLOBAL_BJ} PCPARMSGlobalType;
219: PETSC_EXTERN const char *const PCPARMSGlobalTypes[];
220: /*E
221: PCPARMSLocalType - Determines the local preconditioner method in PARMS
223: Level: intermediate
225: .seealso: PCPARMSSetLocal()
226: E*/
227: typedef enum {PC_PARMS_LOCAL_ILU0,PC_PARMS_LOCAL_ILUK,PC_PARMS_LOCAL_ILUT,PC_PARMS_LOCAL_ARMS} PCPARMSLocalType;
228: PETSC_EXTERN const char *const PCPARMSLocalTypes[];
230: /*E
231: PCGAMGType - type of generalized algebraic multigrid (PCGAMG) method
233: Level: intermediate
235: .seealso: PCMG, PCSetType(), PCGAMGSetThreshold(), PCGAMGSetThreshold(), PCGAMGSetReuseInterpolation()
236: E*/
237: typedef const char *PCGAMGType;
238: #define PCGAMGAGG "agg"
239: #define PCGAMGGEO "geo"
240: #define PCGAMGCLASSICAL "classical"
242: typedef const char *PCGAMGClassicalType;
243: #define PCGAMGCLASSICALDIRECT "direct"
244: #define PCGAMGCLASSICALSTANDARD "standard"
246: /*E
247: PCMGType - Determines the type of multigrid method that is run.
249: Level: beginner
251: Values:
252: + PC_MG_MULTIPLICATIVE (default) - traditional V or W cycle as determined by PCMGSetCycleType()
253: . PC_MG_ADDITIVE - the additive multigrid preconditioner where all levels are
254: smoothed before updating the residual. This only uses the
255: down smoother, in the preconditioner the upper smoother is ignored
256: . PC_MG_FULL - same as multiplicative except one also performs grid sequencing,
257: that is starts on the coarsest grid, performs a cycle, interpolates
258: 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
259: algorithm supports smoothing on before the restriction on each level in the initial restriction to the coarsest stage. In addition that algorithm
260: calls the V-cycle only on the coarser level and has a post-smoother instead.
261: - PC_MG_KASKADE - like full multigrid except one never goes back to a coarser level
262: from a finer
264: .seealso: PCMGSetType(), PCMGSetCycleType(), PCMGSetCycleTypeOnLevel()
266: E*/
267: typedef enum { PC_MG_MULTIPLICATIVE,PC_MG_ADDITIVE,PC_MG_FULL,PC_MG_KASKADE } PCMGType;
268: PETSC_EXTERN const char *const PCMGTypes[];
269: #define PC_MG_CASCADE PC_MG_KASKADE;
271: /*E
272: PCMGCycleType - Use V-cycle or W-cycle
274: Level: beginner
276: Values:
277: + PC_MG_V_CYCLE
278: - PC_MG_W_CYCLE
280: .seealso: PCMGSetCycleType()
282: E*/
283: typedef enum { PC_MG_CYCLE_V = 1,PC_MG_CYCLE_W = 2 } PCMGCycleType;
284: PETSC_EXTERN const char *const PCMGCycleTypes[];
286: /*E
287: PCMGalerkinType - Determines if the coarse grid operators are computed via the Galerkin process
289: Level: beginner
291: Values:
292: + PC_MG_GALERKIN_PMAT - computes the pmat (matrix from which the preconditioner is built) via the Galerkin process from the finest grid
293: . PC_MG_GALERKIN_MAT - computes the mat (matrix used to apply the operator) via the Galerkin process from the finest grid
294: . PC_MG_GALERKIN_BOTH - computes both the mat and pmat via the Galerkin process (if pmat == mat the construction is only done once
295: - PC_MG_GALERKIN_NONE - neither operator is computed via the Galerkin process, the user must provide the operator
297: Users should never set PC_MG_GALERKIN_EXTERNAL, it is used by GAMG and ML
299: .seealso: PCMGSetCycleType()
301: E*/
302: typedef enum { PC_MG_GALERKIN_BOTH,PC_MG_GALERKIN_PMAT,PC_MG_GALERKIN_MAT, PC_MG_GALERKIN_NONE, PC_MG_GALERKIN_EXTERNAL} PCMGGalerkinType;
303: PETSC_EXTERN const char *const PCMGGalerkinTypes[];
305: /*E
306: PCExoticType - Face based or wirebasket based coarse grid space
308: Level: beginner
310: .seealso: PCExoticSetType(), PCEXOTIC
311: E*/
312: typedef enum { PC_EXOTIC_FACE,PC_EXOTIC_WIREBASKET } PCExoticType;
313: PETSC_EXTERN const char *const PCExoticTypes[];
314: PETSC_EXTERN PetscErrorCode PCExoticSetType(PC,PCExoticType);
316: /*E
317: PCPatchConstructType - The algorithm used to construct patches for the preconditioner
319: Level: beginner
321: .seealso: PCPatchSetConstructType(), PCEXOTIC
322: E*/
323: typedef enum {PC_PATCH_STAR, PC_PATCH_VANKA, PC_PATCH_USER, PC_PATCH_PYTHON} PCPatchConstructType;
324: PETSC_EXTERN const char *const PCPatchConstructTypes[];
326: /*E
327: PCFailedReason - indicates type of PC failure
329: Level: beginner
331: Any additions/changes here MUST also be made in include/petsc/finclude/petscpc.h
332: E*/
333: typedef enum {PC_NOERROR,PC_FACTOR_STRUCT_ZEROPIVOT,PC_FACTOR_NUMERIC_ZEROPIVOT,PC_FACTOR_OUTMEMORY,PC_FACTOR_OTHER,PC_SUBPC_ERROR} PCFailedReason;
334: PETSC_EXTERN const char *const PCFailedReasons[];
335: #endif