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 PCSHELL           "shell"
 31: #define PCBJACOBI         "bjacobi"
 32: #define PCMG              "mg"
 33: #define PCEISENSTAT       "eisenstat"
 34: #define PCILU             "ilu"
 35: #define PCICC             "icc"
 36: #define PCASM             "asm"
 37: #define PCGASM            "gasm"
 38: #define PCKSP             "ksp"
 39: #define PCCOMPOSITE       "composite"
 40: #define PCREDUNDANT       "redundant"
 41: #define PCSPAI            "spai"
 42: #define PCNN              "nn"
 43: #define PCCHOLESKY        "cholesky"
 44: #define PCPBJACOBI        "pbjacobi"
 45: #define PCVPBJACOBI       "vpbjacobi"
 46: #define PCMAT             "mat"
 47: #define PCHYPRE           "hypre"
 48: #define PCPARMS           "parms"
 49: #define PCFIELDSPLIT      "fieldsplit"
 50: #define PCTFS             "tfs"
 51: #define PCML              "ml"
 52: #define PCGALERKIN        "galerkin"
 53: #define PCEXOTIC          "exotic"
 54: #define PCCP              "cp"
 55: #define PCBFBT            "bfbt"
 56: #define PCLSC             "lsc"
 57: #define PCPYTHON          "python"
 58: #define PCPFMG            "pfmg"
 59: #define PCSYSPFMG         "syspfmg"
 60: #define PCREDISTRIBUTE    "redistribute"
 61: #define PCSVD             "svd"
 62: #define PCGAMG            "gamg"
 63: #define PCCHOWILUVIENNACL "chowiluviennacl"
 64: #define PCROWSCALINGVIENNACL "rowscalingviennacl"
 65: #define PCSAVIENNACL      "saviennacl"
 66: #define PCBDDC            "bddc"
 67: #define PCKACZMARZ        "kaczmarz"
 68: #define PCTELESCOPE       "telescope"
 69: #define PCPATCH           "patch"
 70: #define PCLMVM            "lmvm"
 71: #define PCHMG             "hmg"
 72: #define PCDEFLATION       "deflation"
 73: #define PCHPDDM           "hpddm"
 74: #define PCHARA            "hara"

 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)

 87: /*E
 88:     PCRichardsonConvergedReason - reason a PCApplyRichardson method terminates

 90:    Level: advanced

 92:    Notes:
 93:     this must match petsc/finclude/petscpc.h and the KSPConvergedReason values in petscksp.h

 95: .seealso: PCApplyRichardson()
 96: E*/
 97: typedef enum {
 98:               PCRICHARDSON_CONVERGED_RTOL               =  2,
 99:               PCRICHARDSON_CONVERGED_ATOL               =  3,
100:               PCRICHARDSON_CONVERGED_ITS                =  4,
101:               PCRICHARDSON_DIVERGED_DTOL                = -4} PCRichardsonConvergedReason;

103: /*E
104:     PCJacobiType - What elements are used to form the Jacobi preconditioner

106:    Level: intermediate

108: .seealso:
109: E*/
110: typedef enum { PC_JACOBI_DIAGONAL,PC_JACOBI_ROWMAX,PC_JACOBI_ROWSUM} PCJacobiType;

112: /*E
113:     PCASMType - Type of additive Schwarz method to use

115: $  PC_ASM_BASIC        - Symmetric version where residuals from the ghost points are used
116: $                        and computed values in ghost regions are added together.
117: $                        Classical standard additive Schwarz.
118: $  PC_ASM_RESTRICT     - Residuals from ghost points are used but computed values in ghost
119: $                        region are discarded.
120: $                        Default.
121: $  PC_ASM_INTERPOLATE  - Residuals from ghost points are not used, computed values in ghost
122: $                        region are added back in.
123: $  PC_ASM_NONE         - Residuals from ghost points are not used, computed ghost values are
124: $                        discarded.
125: $                        Not very good.

127:    Level: beginner

129: .seealso: PCASMSetType()
130: E*/
131: typedef enum {PC_ASM_BASIC = 3,PC_ASM_RESTRICT = 1,PC_ASM_INTERPOLATE = 2,PC_ASM_NONE = 0} PCASMType;

133: /*E
134:     PCGASMType - Type of generalized additive Schwarz method to use (differs from ASM in allowing multiple processors per subdomain).

136:    Each subdomain has nested inner and outer parts.  The inner subdomains are assumed to form a non-overlapping covering of the computational
137:    domain, while the outer subdomains contain the inner subdomains and overlap with each other.  This preconditioner will compute
138:    a subdomain correction over each *outer* subdomain from a residual computed there, but its different variants will differ in
139:    (a) how the outer subdomain residual is computed, and (b) how the outer subdomain correction is computed.

141: $  PC_GASM_BASIC       - Symmetric version where the full from the outer subdomain is used, and the resulting correction is applied
142: $                        over the outer subdomains.  As a result, points in the overlap will receive the sum of the corrections
143: $                        from neighboring subdomains.
144: $                        Classical standard additive Schwarz.
145: $  PC_GASM_RESTRICT    - Residual from the outer subdomain is used but the correction is restricted to the inner subdomain only
146: $                        (i.e., zeroed out over the overlap portion of the outer subdomain before being applied).  As a result,
147: $                        each point will receive a correction only from the unique inner subdomain containing it (nonoverlapping covering
148: $                        assumption).
149: $                        Default.
150: $  PC_GASM_INTERPOLATE - Residual is zeroed out over the overlap portion of the outer subdomain, but the resulting correction is
151: $                        applied over the outer subdomain. As a result, points in the overlap will receive the sum of the corrections
152: $                        from neighboring subdomains.
153: $
154: $  PC_GASM_NONE        - Residuals and corrections are zeroed out outside the local subdomains.
155: $                        Not very good.

157:    Level: beginner

159: .seealso: PCGASMSetType()
160: E*/
161: typedef enum {PC_GASM_BASIC = 3,PC_GASM_RESTRICT = 1,PC_GASM_INTERPOLATE = 2,PC_GASM_NONE = 0} PCGASMType;

163: /*E
164:     PCCompositeType - Determines how two or more preconditioner are composed

166: $  PC_COMPOSITE_ADDITIVE - results from application of all preconditioners are added together
167: $  PC_COMPOSITE_MULTIPLICATIVE - preconditioners are applied sequentially to the residual freshly
168: $                                computed after the previous preconditioner application
169: $  PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE - preconditioners are applied sequentially to the residual freshly
170: $                                computed from first preconditioner to last and then back (Use only for symmetric matrices and preconditioners)
171: $  PC_COMPOSITE_SPECIAL - This is very special for a matrix of the form alpha I + R + S
172: $                         where first preconditioner is built from alpha I + S and second from
173: $                         alpha I + R

175:    Level: beginner

177: .seealso: PCCompositeSetType()
178: E*/
179: typedef enum {PC_COMPOSITE_ADDITIVE,PC_COMPOSITE_MULTIPLICATIVE,PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE,PC_COMPOSITE_SPECIAL,PC_COMPOSITE_SCHUR,PC_COMPOSITE_GKB} PCCompositeType;

181: /*E
182:     PCFieldSplitSchurPreType - Determines how to precondition Schur complement

184:     Level: intermediate

186: .seealso: PCFieldSplitSetSchurPre()
187: E*/
188: 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;

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

193:     Level: intermediate

195: .seealso: PCFieldSplitSetSchurFactType()
196: E*/
197: typedef enum {
198:   PC_FIELDSPLIT_SCHUR_FACT_DIAG,
199:   PC_FIELDSPLIT_SCHUR_FACT_LOWER,
200:   PC_FIELDSPLIT_SCHUR_FACT_UPPER,
201:   PC_FIELDSPLIT_SCHUR_FACT_FULL
202: } PCFieldSplitSchurFactType;

204: /*E
205:     PCPARMSGlobalType - Determines the global preconditioner method in PARMS

207:     Level: intermediate

209: .seealso: PCPARMSSetGlobal()
210: E*/
211: typedef enum {PC_PARMS_GLOBAL_RAS,PC_PARMS_GLOBAL_SCHUR,PC_PARMS_GLOBAL_BJ} PCPARMSGlobalType;

213: /*E
214:     PCPARMSLocalType - Determines the local preconditioner method in PARMS

216:     Level: intermediate

218: .seealso: PCPARMSSetLocal()
219: E*/
220: typedef enum {PC_PARMS_LOCAL_ILU0,PC_PARMS_LOCAL_ILUK,PC_PARMS_LOCAL_ILUT,PC_PARMS_LOCAL_ARMS} PCPARMSLocalType;

222: /*J
223:     PCGAMGType - type of generalized algebraic multigrid (PCGAMG) method

225:     Level: intermediate

227: $   PCGAMGAGG - (the default) smoothed aggregation algorithm, robust, very well tested
228: $   PCGAMGGEO - geometric coarsening, uses mesh generator to produce coarser meshes, limited to triangles, not well tested
229: $   PCGAMGCLASSICAL - classical algebraic multigrid preconditioner, incomplete, poorly tested

231: .seealso: PCMG, PCSetType(), PCGAMGSetThreshold(), PCGAMGSetThreshold(), PCGAMGSetReuseInterpolation()
232: J*/
233: typedef const char *PCGAMGType;
234: #define PCGAMGAGG         "agg"
235: #define PCGAMGGEO         "geo"
236: #define PCGAMGCLASSICAL   "classical"

238: typedef const char *PCGAMGClassicalType;
239: #define PCGAMGCLASSICALDIRECT   "direct"
240: #define PCGAMGCLASSICALSTANDARD "standard"

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

245:    Level: beginner

247:    Values:
248: +  PC_MG_MULTIPLICATIVE (default) - traditional V or W cycle as determined by PCMGSetCycleType()
249: .  PC_MG_ADDITIVE - the additive multigrid preconditioner where all levels are
250:                 smoothed before updating the residual. This only uses the
251:                 down smoother, in the preconditioner the upper smoother is ignored
252: .  PC_MG_FULL - same as multiplicative except one also performs grid sequencing,
253:             that is starts on the coarsest grid, performs a cycle, interpolates
254:             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
255:             algorithm supports smoothing on before the restriction on each level in the initial restriction to the coarsest stage. In addition that algorithm
256:             calls the V-cycle only on the coarser level and has a post-smoother instead.
257: -  PC_MG_KASKADE - like full multigrid except one never goes back to a coarser level
258:                from a finer

260: .seealso: PCMGSetType(), PCMGSetCycleType(), PCMGSetCycleTypeOnLevel()

262: E*/
263: typedef enum { PC_MG_MULTIPLICATIVE,PC_MG_ADDITIVE,PC_MG_FULL,PC_MG_KASKADE } PCMGType;
264: #define PC_MG_CASCADE PC_MG_KASKADE;

266: /*E
267:     PCMGCycleType - Use V-cycle or W-cycle

269:    Level: beginner

271:    Values:
272: +  PC_MG_V_CYCLE
273: -  PC_MG_W_CYCLE

275: .seealso: PCMGSetCycleType()

277: E*/
278: typedef enum { PC_MG_CYCLE_V = 1,PC_MG_CYCLE_W = 2 } PCMGCycleType;

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

283:    Level: beginner

285:    Values:
286: +  PC_MG_GALERKIN_PMAT - computes the pmat (matrix from which the preconditioner is built) via the Galerkin process from the finest grid
287: .  PC_MG_GALERKIN_MAT -  computes the mat (matrix used to apply the operator) via the Galerkin process from the finest grid
288: .  PC_MG_GALERKIN_BOTH - computes both the mat and pmat via the Galerkin process (if pmat == mat the construction is only done once
289: -  PC_MG_GALERKIN_NONE - neither operator is computed via the Galerkin process, the user must provide the operator

291:    Users should never set PC_MG_GALERKIN_EXTERNAL, it is used by GAMG and ML

293: .seealso: PCMGSetCycleType()

295: E*/
296: typedef enum { PC_MG_GALERKIN_BOTH,PC_MG_GALERKIN_PMAT,PC_MG_GALERKIN_MAT, PC_MG_GALERKIN_NONE, PC_MG_GALERKIN_EXTERNAL} PCMGGalerkinType;

298: /*E
299:     PCExoticType - Face based or wirebasket based coarse grid space

301:    Level: beginner

303: .seealso: PCExoticSetType(), PCEXOTIC
304: E*/
305: typedef enum { PC_EXOTIC_FACE,PC_EXOTIC_WIREBASKET } PCExoticType;

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

310:    Level: intermediate

312:    Values:
313: +  PC_BDDC_INTERFACE_EXT_DIRICHLET - solves Dirichlet interior problem; this is the standard BDDC algorithm
314: -  PC_BDDC_INTERFACE_EXT_LUMP - skips interior solve; sometimes called M_1 and associated with "lumped FETI-DP"

316: E*/
317: typedef enum {
318:   PC_BDDC_INTERFACE_EXT_DIRICHLET,
319:   PC_BDDC_INTERFACE_EXT_LUMP
320: } PCBDDCInterfaceExtType;

322: /*E
323:   PCMGCoarseSpaceType - Function space for coarse space for adaptive interpolation

325:   Level: beginner

327: .seealso: PCMGSetAdaptCoarseSpaceType(), PCMG
328: E*/
329: typedef enum { PCMG_POLYNOMIAL, PCMG_HARMONIC, PCMG_EIGENVECTOR, PCMG_GENERALIZED_EIGENVECTOR } PCMGCoarseSpaceType;

331: /*E
332:     PCPatchConstructType - The algorithm used to construct patches for the preconditioner

334:    Level: beginner

336: .seealso: PCPatchSetConstructType(), PCEXOTIC
337: E*/
338: typedef enum {PC_PATCH_STAR, PC_PATCH_VANKA, PC_PATCH_PARDECOMP, PC_PATCH_USER, PC_PATCH_PYTHON} PCPatchConstructType;

340: /*E
341:     PCDeflationSpaceType - Type of deflation

343:     Values:
344: +   PC_DEFLATION_SPACE_HAAR        - directly assembled based on Haar (db2) wavelet with overflowed filter cuted-off
345: .   PC_DEFLATION_SPACE_DB2         - MATCOMPOSITE of 1-lvl matices based on db2 (2 coefficient Daubechies / Haar wavelet)
346: .   PC_DEFLATION_SPACE_DB4         - same as above, but with db4 (4 coefficient Daubechies)
347: .   PC_DEFLATION_SPACE_DB8         - same as above, but with db8 (8 coefficient Daubechies)
348: .   PC_DEFLATION_SPACE_DB16        - same as above, but with db16 (16 coefficient Daubechies)
349: .   PC_DEFLATION_SPACE_BIORTH22    - same as above, but with biorthogonal 2.2 (6 coefficients)
350: .   PC_DEFLATION_SPACE_MEYER       - same as above, but with Meyer/FIR (62 coefficients)
351: .   PC_DEFLATION_SPACE_AGGREGATION - aggregates local indices (given by operator matix distribution) into a subdomain
352: -   PC_DEFLATION_SPACE_USER        - indicates space set by user

354:     Notes:
355:       Wavelet-based space (except Haar) can be used in multilevel deflation.

357:     Level: intermediate

359: .seealso: PCDeflationSetSpaceToCompute(), PCDEFLATION
360: E*/
361: typedef enum {
362:   PC_DEFLATION_SPACE_HAAR,
363:   PC_DEFLATION_SPACE_DB2,
364:   PC_DEFLATION_SPACE_DB4,
365:   PC_DEFLATION_SPACE_DB8,
366:   PC_DEFLATION_SPACE_DB16,
367:   PC_DEFLATION_SPACE_BIORTH22,
368:   PC_DEFLATION_SPACE_MEYER,
369:   PC_DEFLATION_SPACE_AGGREGATION,
370:   PC_DEFLATION_SPACE_USER
371: } PCDeflationSpaceType;

373: /*E
374:     PCHPDDMCoarseCorrectionType - Type of coarse correction used by PCHPDDM

376:     Level: intermediate

378:     Values:
379: +   PC_HPDDM_COARSE_CORRECTION_DEFLATED (default) - eq. (1) in PCHPDDMShellApply()
380: .   PC_HPDDM_COARSE_CORRECTION_ADDITIVE - eq. (2)
381: -   PC_HPDDM_COARSE_CORRECTION_BALANCED - eq. (3)

383: .seealso: PCHPDDM, PCSetType(), PCHPDDMShellApply()
384: E*/
385: typedef enum { PC_HPDDM_COARSE_CORRECTION_DEFLATED, PC_HPDDM_COARSE_CORRECTION_ADDITIVE, PC_HPDDM_COARSE_CORRECTION_BALANCED } PCHPDDMCoarseCorrectionType;

387: /*E
388:     PCFailedReason - indicates type of PC failure

390:     Level: beginner

392:     Any additions/changes here MUST also be made in include/petsc/finclude/petscpc.h
393: E*/
394: 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;

396: /*E
397:     PCGAMGLayoutType - Layout for reduced grids

399:     Level: intermediate

401: .seealso: PCGAMGSetCoarseGridLayoutType()
402:     Any additions/changes here MUST also be made in include/petsc/finclude/petscpc.h
403: E*/
404: typedef enum {PCGAMG_LAYOUT_COMPACT,PCGAMG_LAYOUT_SPREAD} PCGAMGLayoutType;

406: #endif