Actual source code: petscpc.h
1: /*
2: Preconditioner module.
3: */
6: #include petscdm.h
7: PETSC_EXTERN_CXX_BEGIN
9: extern PetscErrorCode PCInitializePackage(const char[]);
11: /*
12: PCList contains the list of preconditioners currently registered
13: These are added with the PCRegisterDynamic() macro
14: */
15: extern PetscFList PCList;
17: /*S
18: PC - Abstract PETSc object that manages all preconditioners
20: Level: beginner
22: Concepts: preconditioners
24: .seealso: PCCreate(), PCSetType(), PCType (for list of available types)
25: S*/
26: typedef struct _p_PC* PC;
28: /*J
29: PCType - String with the name of a PETSc preconditioner method or the creation function
30: with an optional dynamic library name, for example
31: http://www.mcs.anl.gov/petsc/lib.a:mypccreate()
33: Level: beginner
35: Notes: Click on the links below to see details on a particular solver
37: .seealso: PCSetType(), PC, PCCreate()
38: J*/
39: #define PCType char*
40: #define PCNONE "none"
41: #define PCJACOBI "jacobi"
42: #define PCSOR "sor"
43: #define PCLU "lu"
44: #define PCSHELL "shell"
45: #define PCBJACOBI "bjacobi"
46: #define PCMG "mg"
47: #define PCEISENSTAT "eisenstat"
48: #define PCILU "ilu"
49: #define PCICC "icc"
50: #define PCASM "asm"
51: #define PCGASM "gasm"
52: #define PCKSP "ksp"
53: #define PCCOMPOSITE "composite"
54: #define PCREDUNDANT "redundant"
55: #define PCSPAI "spai"
56: #define PCNN "nn"
57: #define PCCHOLESKY "cholesky"
58: #define PCPBJACOBI "pbjacobi"
59: #define PCMAT "mat"
60: #define PCHYPRE "hypre"
61: #define PCPARMS "parms"
62: #define PCFIELDSPLIT "fieldsplit"
63: #define PCTFS "tfs"
64: #define PCML "ml"
65: #define PCPROMETHEUS "prometheus"
66: #define PCGALERKIN "galerkin"
67: #define PCEXOTIC "exotic"
68: #define PCHMPI "hmpi"
69: #define PCSUPPORTGRAPH "supportgraph"
70: #define PCASA "asa"
71: #define PCCP "cp"
72: #define PCBFBT "bfbt"
73: #define PCLSC "lsc"
74: #define PCPYTHON "python"
75: #define PCPFMG "pfmg"
76: #define PCSYSPFMG "syspfmg"
77: #define PCREDISTRIBUTE "redistribute"
78: #define PCSACUSP "sacusp"
79: #define PCSACUSPPOLY "sacusppoly"
80: #define PCBICGSTABCUSP "bicgstabcusp"
81: #define PCSVD "svd"
82: #define PCAINVCUSP "ainvcusp"
83: #define PCGAMG "gamg"
85: /* Logging support */
86: extern PetscClassId PC_CLASSID;
88: /*E
89: PCSide - If the preconditioner is to be applied to the left, right
90: or symmetrically around the operator.
92: Level: beginner
94: .seealso:
95: E*/
96: typedef enum { PC_SIDE_DEFAULT=-1,PC_LEFT,PC_RIGHT,PC_SYMMETRIC} PCSide;
97: #define PC_SIDE_MAX (PC_SYMMETRIC + 1)
98: extern const char *PCSides[];
100: extern PetscErrorCode PCCreate(MPI_Comm,PC*);
101: extern PetscErrorCode PCSetType(PC,const PCType);
102: extern PetscErrorCode PCSetUp(PC);
103: extern PetscErrorCode PCSetUpOnBlocks(PC);
104: extern PetscErrorCode PCApply(PC,Vec,Vec);
105: extern PetscErrorCode PCApplySymmetricLeft(PC,Vec,Vec);
106: extern PetscErrorCode PCApplySymmetricRight(PC,Vec,Vec);
107: extern PetscErrorCode PCApplyBAorAB(PC,PCSide,Vec,Vec,Vec);
108: extern PetscErrorCode PCApplyTranspose(PC,Vec,Vec);
109: extern PetscErrorCode PCApplyTransposeExists(PC,PetscBool *);
110: extern PetscErrorCode PCApplyBAorABTranspose(PC,PCSide,Vec,Vec,Vec);
112: /*E
113: PCRichardsonConvergedReason - reason a PCApplyRichardson method terminates
115: Level: advanced
117: Notes: this must match finclude/petscpc.h and the KSPConvergedReason values in petscksp.h
119: .seealso: PCApplyRichardson()
120: E*/
121: typedef enum {
122: PCRICHARDSON_CONVERGED_RTOL = 2,
123: PCRICHARDSON_CONVERGED_ATOL = 3,
124: PCRICHARDSON_CONVERGED_ITS = 4,
125: PCRICHARDSON_DIVERGED_DTOL = -4} PCRichardsonConvergedReason;
127: extern PetscErrorCode PCApplyRichardson(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool ,PetscInt*,PCRichardsonConvergedReason*);
128: extern PetscErrorCode PCApplyRichardsonExists(PC,PetscBool *);
129: extern PetscErrorCode PCSetInitialGuessNonzero(PC,PetscBool );
131: extern PetscErrorCode PCRegisterDestroy(void);
132: extern PetscErrorCode PCRegisterAll(const char[]);
133: extern PetscBool PCRegisterAllCalled;
135: extern PetscErrorCode PCRegister(const char[],const char[],const char[],PetscErrorCode(*)(PC));
137: /*MC
138: PCRegisterDynamic - Adds a method to the preconditioner package.
140: Synopsis:
141: PetscErrorCode PCRegisterDynamic(const char *name_solver,const char *path,const char *name_create,PetscErrorCode (*routine_create)(PC))
143: Not collective
145: Input Parameters:
146: + name_solver - name of a new user-defined solver
147: . path - path (either absolute or relative) the library containing this solver
148: . name_create - name of routine to create method context
149: - routine_create - routine to create method context
151: Notes:
152: PCRegisterDynamic() may be called multiple times to add several user-defined preconditioners.
154: If dynamic libraries are used, then the fourth input argument (routine_create)
155: is ignored.
157: Sample usage:
158: .vb
159: PCRegisterDynamic("my_solver","/home/username/my_lib/lib/libO/solaris/mylib",
160: "MySolverCreate",MySolverCreate);
161: .ve
163: Then, your solver can be chosen with the procedural interface via
164: $ PCSetType(pc,"my_solver")
165: or at runtime via the option
166: $ -pc_type my_solver
168: Level: advanced
170: Notes: ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR}, or ${any environmental variable}
171: occuring in pathname will be replaced with appropriate values.
172: If your function is not being put into a shared library then use PCRegister() instead
174: .keywords: PC, register
176: .seealso: PCRegisterAll(), PCRegisterDestroy()
177: M*/
178: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
179: #define PCRegisterDynamic(a,b,c,d) PCRegister(a,b,c,0)
180: #else
181: #define PCRegisterDynamic(a,b,c,d) PCRegister(a,b,c,d)
182: #endif
184: extern PetscErrorCode PCReset(PC);
185: extern PetscErrorCode PCDestroy(PC*);
186: extern PetscErrorCode PCSetFromOptions(PC);
187: extern PetscErrorCode PCGetType(PC,const PCType*);
189: extern PetscErrorCode PCFactorGetMatrix(PC,Mat*);
190: extern PetscErrorCode PCSetModifySubMatrices(PC,PetscErrorCode(*)(PC,PetscInt,const IS[],const IS[],Mat[],void*),void*);
191: extern PetscErrorCode PCModifySubMatrices(PC,PetscInt,const IS[],const IS[],Mat[],void*);
193: extern PetscErrorCode PCSetOperators(PC,Mat,Mat,MatStructure);
194: extern PetscErrorCode PCGetOperators(PC,Mat*,Mat*,MatStructure*);
195: extern PetscErrorCode PCGetOperatorsSet(PC,PetscBool *,PetscBool *);
197: extern PetscErrorCode PCView(PC,PetscViewer);
199: extern PetscErrorCode PCSetOptionsPrefix(PC,const char[]);
200: extern PetscErrorCode PCAppendOptionsPrefix(PC,const char[]);
201: extern PetscErrorCode PCGetOptionsPrefix(PC,const char*[]);
203: extern PetscErrorCode PCComputeExplicitOperator(PC,Mat*);
205: /*
206: These are used to provide extra scaling of preconditioned
207: operator for time-stepping schemes like in SUNDIALS
208: */
209: extern PetscErrorCode PCGetDiagonalScale(PC,PetscBool *);
210: extern PetscErrorCode PCDiagonalScaleLeft(PC,Vec,Vec);
211: extern PetscErrorCode PCDiagonalScaleRight(PC,Vec,Vec);
212: extern PetscErrorCode PCSetDiagonalScale(PC,Vec);
214: /* ------------- options specific to particular preconditioners --------- */
216: extern PetscErrorCode PCJacobiSetUseRowMax(PC);
217: extern PetscErrorCode PCJacobiSetUseRowSum(PC);
218: extern PetscErrorCode PCJacobiSetUseAbs(PC);
219: extern PetscErrorCode PCSORSetSymmetric(PC,MatSORType);
220: extern PetscErrorCode PCSORSetOmega(PC,PetscReal);
221: extern PetscErrorCode PCSORSetIterations(PC,PetscInt,PetscInt);
223: extern PetscErrorCode PCEisenstatSetOmega(PC,PetscReal);
224: extern PetscErrorCode PCEisenstatNoDiagonalScaling(PC);
226: #define USE_PRECONDITIONER_MATRIX 0
227: #define USE_TRUE_MATRIX 1
228: extern PetscErrorCode PCBJacobiSetUseTrueLocal(PC);
229: extern PetscErrorCode PCBJacobiSetTotalBlocks(PC,PetscInt,const PetscInt[]);
230: extern PetscErrorCode PCBJacobiSetLocalBlocks(PC,PetscInt,const PetscInt[]);
232: extern PetscErrorCode PCKSPSetUseTrue(PC);
234: extern PetscErrorCode PCShellSetApply(PC,PetscErrorCode (*)(PC,Vec,Vec));
235: extern PetscErrorCode PCShellSetApplyBA(PC,PetscErrorCode (*)(PC,PCSide,Vec,Vec,Vec));
236: extern PetscErrorCode PCShellSetApplyTranspose(PC,PetscErrorCode (*)(PC,Vec,Vec));
237: extern PetscErrorCode PCShellSetSetUp(PC,PetscErrorCode (*)(PC));
238: extern PetscErrorCode PCShellSetApplyRichardson(PC,PetscErrorCode (*)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool ,PetscInt*,PCRichardsonConvergedReason*));
239: extern PetscErrorCode PCShellSetView(PC,PetscErrorCode (*)(PC,PetscViewer));
240: extern PetscErrorCode PCShellSetDestroy(PC,PetscErrorCode (*)(PC));
241: extern PetscErrorCode PCShellGetContext(PC,void**);
242: extern PetscErrorCode PCShellSetContext(PC,void*);
243: extern PetscErrorCode PCShellSetName(PC,const char[]);
244: extern PetscErrorCode PCShellGetName(PC,char*[]);
246: extern PetscErrorCode PCFactorSetZeroPivot(PC,PetscReal);
248: extern PetscErrorCode PCFactorSetShiftType(PC,MatFactorShiftType);
249: extern PetscErrorCode PCFactorSetShiftAmount(PC,PetscReal);
251: extern PetscErrorCode PCFactorSetMatSolverPackage(PC,const MatSolverPackage);
252: extern PetscErrorCode PCFactorGetMatSolverPackage(PC,const MatSolverPackage*);
253: extern PetscErrorCode PCFactorSetUpMatSolverPackage(PC);
255: extern PetscErrorCode PCFactorSetFill(PC,PetscReal);
256: extern PetscErrorCode PCFactorSetColumnPivot(PC,PetscReal);
257: extern PetscErrorCode PCFactorReorderForNonzeroDiagonal(PC,PetscReal);
259: extern PetscErrorCode PCFactorSetMatOrderingType(PC,const MatOrderingType);
260: extern PetscErrorCode PCFactorSetReuseOrdering(PC,PetscBool );
261: extern PetscErrorCode PCFactorSetReuseFill(PC,PetscBool );
262: extern PetscErrorCode PCFactorSetUseInPlace(PC);
263: extern PetscErrorCode PCFactorSetAllowDiagonalFill(PC);
264: extern PetscErrorCode PCFactorSetPivotInBlocks(PC,PetscBool );
266: extern PetscErrorCode PCFactorSetLevels(PC,PetscInt);
267: extern PetscErrorCode PCFactorSetDropTolerance(PC,PetscReal,PetscReal,PetscInt);
269: extern PetscErrorCode PCASMSetLocalSubdomains(PC,PetscInt,IS[],IS[]);
270: extern PetscErrorCode PCASMSetTotalSubdomains(PC,PetscInt,IS[],IS[]);
271: extern PetscErrorCode PCASMSetOverlap(PC,PetscInt);
272: extern PetscErrorCode PCASMSetSortIndices(PC,PetscBool );
274: /*E
275: PCASMType - Type of additive Schwarz method to use
277: $ PC_ASM_BASIC - symmetric version where residuals from the ghost points are used
278: $ and computed values in ghost regions are added together. Classical
279: $ standard additive Schwarz
280: $ PC_ASM_RESTRICT - residuals from ghost points are used but computed values in ghost
281: $ region are discarded. Default
282: $ PC_ASM_INTERPOLATE - residuals from ghost points are not used, computed values in ghost
283: $ region are added back in
284: $ PC_ASM_NONE - ghost point residuals are not used, computed ghost values are discarded
285: $ not very good.
287: Level: beginner
289: .seealso: PCASMSetType()
290: E*/
291: typedef enum {PC_ASM_BASIC = 3,PC_ASM_RESTRICT = 1,PC_ASM_INTERPOLATE = 2,PC_ASM_NONE = 0} PCASMType;
292: extern const char *PCASMTypes[];
294: extern PetscErrorCode PCASMSetType(PC,PCASMType);
295: extern PetscErrorCode PCASMCreateSubdomains(Mat,PetscInt,IS*[]);
296: extern PetscErrorCode PCASMDestroySubdomains(PetscInt,IS[],IS[]);
297: extern PetscErrorCode PCASMCreateSubdomains2D(PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,IS**,IS**);
298: extern PetscErrorCode PCASMGetLocalSubdomains(PC,PetscInt*,IS*[],IS*[]);
299: extern PetscErrorCode PCASMGetLocalSubmatrices(PC,PetscInt*,Mat*[]);
301: /*E
302: PCGASMType - Type of generalized additive Schwarz method to use (differs from ASM in allowing multiple processors per domain)
304: $ PC_GASM_BASIC - symmetric version where residuals from the ghost points are used
305: $ and computed values in ghost regions are added together. Classical
306: $ standard additive Schwarz
307: $ PC_GASM_RESTRICT - residuals from ghost points are used but computed values in ghost
308: $ region are discarded. Default
309: $ PC_GASM_INTERPOLATE - residuals from ghost points are not used, computed values in ghost
310: $ region are added back in
311: $ PC_GASM_NONE - ghost point residuals are not used, computed ghost values are discarded
312: $ not very good.
314: Level: beginner
316: .seealso: PCGASMSetType()
317: E*/
318: typedef enum {PC_GASM_BASIC = 3,PC_GASM_RESTRICT = 1,PC_GASM_INTERPOLATE = 2,PC_GASM_NONE = 0} PCGASMType;
319: extern const char *PCGASMTypes[];
321: extern PetscErrorCode PCGASMSetLocalSubdomains(PC,PetscInt,IS[],IS[]);
322: extern PetscErrorCode PCGASMSetTotalSubdomains(PC,PetscInt);
323: extern PetscErrorCode PCGASMSetOverlap(PC,PetscInt);
324: extern PetscErrorCode PCGASMSetSortIndices(PC,PetscBool );
326: extern PetscErrorCode PCGASMSetType(PC,PCGASMType);
327: extern PetscErrorCode PCGASMCreateSubdomains(Mat,PetscInt,IS*[]);
328: extern PetscErrorCode PCGASMDestroySubdomains(PetscInt,IS[],IS[]);
329: extern PetscErrorCode PCGASMCreateSubdomains2D(PC,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,IS**,IS**);
330: extern PetscErrorCode PCGASMGetLocalSubdomains(PC,PetscInt*,IS*[],IS*[]);
331: extern PetscErrorCode PCGASMGetLocalSubmatrices(PC,PetscInt*,Mat*[]);
333: /*E
334: PCCompositeType - Determines how two or more preconditioner are composed
336: $ PC_COMPOSITE_ADDITIVE - results from application of all preconditioners are added together
337: $ PC_COMPOSITE_MULTIPLICATIVE - preconditioners are applied sequentially to the residual freshly
338: $ computed after the previous preconditioner application
339: $ PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE - preconditioners are applied sequentially to the residual freshly
340: $ computed from first preconditioner to last and then back (Use only for symmetric matrices and preconditions)
341: $ PC_COMPOSITE_SPECIAL - This is very special for a matrix of the form alpha I + R + S
342: $ where first preconditioner is built from alpha I + S and second from
343: $ alpha I + R
345: Level: beginner
347: .seealso: PCCompositeSetType()
348: E*/
349: typedef enum {PC_COMPOSITE_ADDITIVE,PC_COMPOSITE_MULTIPLICATIVE,PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE,PC_COMPOSITE_SPECIAL,PC_COMPOSITE_SCHUR} PCCompositeType;
350: extern const char *PCCompositeTypes[];
352: extern PetscErrorCode PCCompositeSetUseTrue(PC);
353: extern PetscErrorCode PCCompositeSetType(PC,PCCompositeType);
354: extern PetscErrorCode PCCompositeAddPC(PC,PCType);
355: extern PetscErrorCode PCCompositeGetPC(PC,PetscInt,PC *);
356: extern PetscErrorCode PCCompositeSpecialSetAlpha(PC,PetscScalar);
358: extern PetscErrorCode PCRedundantSetNumber(PC,PetscInt);
359: extern PetscErrorCode PCRedundantSetScatter(PC,VecScatter,VecScatter);
360: extern PetscErrorCode PCRedundantGetOperators(PC,Mat*,Mat*);
362: extern PetscErrorCode PCSPAISetEpsilon(PC,double);
363: extern PetscErrorCode PCSPAISetNBSteps(PC,PetscInt);
364: extern PetscErrorCode PCSPAISetMax(PC,PetscInt);
365: extern PetscErrorCode PCSPAISetMaxNew(PC,PetscInt);
366: extern PetscErrorCode PCSPAISetBlockSize(PC,PetscInt);
367: extern PetscErrorCode PCSPAISetCacheSize(PC,PetscInt);
368: extern PetscErrorCode PCSPAISetVerbose(PC,PetscInt);
369: extern PetscErrorCode PCSPAISetSp(PC,PetscInt);
371: extern PetscErrorCode PCHYPRESetType(PC,const char[]);
372: extern PetscErrorCode PCHYPREGetType(PC,const char*[]);
373: extern PetscErrorCode PCBJacobiGetLocalBlocks(PC,PetscInt*,const PetscInt*[]);
374: extern PetscErrorCode PCBJacobiGetTotalBlocks(PC,PetscInt*,const PetscInt*[]);
376: extern PetscErrorCode PCFieldSplitSetFields(PC,const char[],PetscInt,const PetscInt*);
377: extern PetscErrorCode PCFieldSplitSetType(PC,PCCompositeType);
378: extern PetscErrorCode PCFieldSplitSetBlockSize(PC,PetscInt);
379: extern PetscErrorCode PCFieldSplitSetIS(PC,const char[],IS);
380: extern PetscErrorCode PCFieldSplitGetIS(PC,const char[],IS*);
382: /*E
383: PCFieldSplitSchurPreType - Determines how to precondition Schur complement
385: Level: intermediate
387: .seealso: PCFieldSplitSchurPrecondition()
388: E*/
389: typedef enum {PC_FIELDSPLIT_SCHUR_PRE_SELF,PC_FIELDSPLIT_SCHUR_PRE_DIAG,PC_FIELDSPLIT_SCHUR_PRE_USER} PCFieldSplitSchurPreType;
390: extern const char *const PCFieldSplitSchurPreTypes[];
392: extern PetscErrorCode PCFieldSplitSchurPrecondition(PC,PCFieldSplitSchurPreType,Mat);
393: extern PetscErrorCode PCFieldSplitGetSchurBlocks(PC,Mat*,Mat*,Mat*,Mat*);
395: extern PetscErrorCode PCGalerkinSetRestriction(PC,Mat);
396: extern PetscErrorCode PCGalerkinSetInterpolation(PC,Mat);
398: extern PetscErrorCode PCSetCoordinates(PC,PetscInt,PetscReal*);
399: extern PetscErrorCode PCSASetVectors(PC,PetscInt,PetscReal *);
401: extern PetscErrorCode PCPythonSetType(PC,const char[]);
403: extern PetscErrorCode PCSetDM(PC,DM);
404: extern PetscErrorCode PCGetDM(PC,DM*);
406: extern PetscErrorCode PCSetApplicationContext(PC,void*);
407: extern PetscErrorCode PCGetApplicationContext(PC,void*);
409: extern PetscErrorCode PCBiCGStabCUSPSetTolerance(PC,PetscReal);
410: extern PetscErrorCode PCBiCGStabCUSPSetIterations(PC,PetscInt);
411: extern PetscErrorCode PCBiCGStabCUSPSetUseVerboseMonitor(PC,PetscBool);
413: extern PetscErrorCode PCAINVCUSPSetDropTolerance(PC,PetscReal);
414: extern PetscErrorCode PCAINVCUSPUseScaling(PC,PetscBool);
415: extern PetscErrorCode PCAINVCUSPSetNonzeros(PC,PetscInt);
416: extern PetscErrorCode PCAINVCUSPSetLinParameter(PC,PetscInt);
417: /*E
418: PCPARMSGlobalType - Determines the global preconditioner method in PARMS
420: Level: intermediate
422: .seealso: PCPARMSSetGlobal()
423: E*/
424: typedef enum {PC_PARMS_GLOBAL_RAS,PC_PARMS_GLOBAL_SCHUR,PC_PARMS_GLOBAL_BJ} PCPARMSGlobalType;
425: extern const char *PCPARMSGlobalTypes[];
426: /*E
427: PCPARMSLocalType - Determines the local preconditioner method in PARMS
429: Level: intermediate
431: .seealso: PCPARMSSetLocal()
432: E*/
433: typedef enum {PC_PARMS_LOCAL_ILU0,PC_PARMS_LOCAL_ILUK,PC_PARMS_LOCAL_ILUT,PC_PARMS_LOCAL_ARMS} PCPARMSLocalType;
434: extern const char *PCPARMSLocalTypes[];
436: extern PetscErrorCode PCPARMSSetGlobal(PC pc,PCPARMSGlobalType type);
437: extern PetscErrorCode PCPARMSSetLocal(PC pc,PCPARMSLocalType type);
438: extern PetscErrorCode PCPARMSSetSolveTolerances(PC pc,PetscReal tol,PetscInt maxits);
439: extern PetscErrorCode PCPARMSSetSolveRestart(PC pc,PetscInt restart);
440: extern PetscErrorCode PCPARMSSetNonsymPerm(PC pc,PetscBool nonsym);
441: extern PetscErrorCode PCPARMSSetFill(PC pc,PetscInt lfil0,PetscInt lfil1,PetscInt lfil2);
443: PETSC_EXTERN_CXX_END
445: #endif /* __PETSCPC_H */