Actual source code: petscpc.h
1: /* $Id: petscpc.h,v 1.122 2001/08/21 21:03:12 bsmith Exp $ */
3: /*
4: Preconditioner module.
5: */
8: #include petscmat.h
10: /*
11: PCList contains the list of preconditioners currently registered
12: These are added with the PCRegisterDynamic() macro
13: */
14: extern PetscFList PCList;
15: typedef char *PCType;
17: /*S
18: PC - Abstract PETSc object that manages all preconditioners
20: Level: beginner
22: Concepts: preconditioners
24: .seealso: PCCreate(), PCSetType(), PCType
25: S*/
26: typedef struct _p_PC* PC;
28: /*E
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: .seealso: PCSetType(), PC
36: E*/
37: #define PCNONE "none"
38: #define PCJACOBI "jacobi"
39: #define PCSOR "sor"
40: #define PCLU "lu"
41: #define PCSHELL "shell"
42: #define PCBJACOBI "bjacobi"
43: #define PCMG "mg"
44: #define PCEISENSTAT "eisenstat"
45: #define PCILU "ilu"
46: #define PCICC "icc"
47: #define PCASM "asm"
48: #define PCSLES "sles"
49: #define PCCOMPOSITE "composite"
50: #define PCREDUNDANT "redundant"
51: #define PCSPAI "spai"
52: #define PCMILU "milu"
53: #define PCNN "nn"
54: #define PCCHOLESKY "cholesky"
55: #define PCRAMG "ramg"
56: #define PCSAMG "samg"
57: #define PCPBJACOBI "pbjacobi"
58: #define PCMULTILEVEL "multilevel"
59: #define PCSCHUR "schur"
60: #define PCESI "esi"
61: #define PCPETSCESI "petscesi"
62: #define PCMAT "mat"
63: #define PCHYPRE "hypre"
65: /* Logging support */
66: extern int PC_COOKIE;
67: extern int PC_SetUp, PC_SetUpOnBlocks, PC_Apply, PC_ApplyCoarse, PC_ApplyMultiple, PC_ApplySymmetricLeft;
68: extern int PC_ApplySymmetricRight, PC_ModifySubMatrices;
70: /*E
71: PCSide - If the preconditioner is to be applied to the left, right
72: or symmetrically around the operator.
74: Level: beginner
76: .seealso:
77: E*/
78: typedef enum { PC_LEFT,PC_RIGHT,PC_SYMMETRIC } PCSide;
80: EXTERN int PCCreate(MPI_Comm,PC*);
81: EXTERN int PCSetType(PC,PCType);
82: EXTERN int PCSetUp(PC);
83: EXTERN int PCSetUpOnBlocks(PC);
84: EXTERN int PCApply(PC,Vec,Vec);
85: EXTERN int PCApplySymmetricLeft(PC,Vec,Vec);
86: EXTERN int PCApplySymmetricRight(PC,Vec,Vec);
87: EXTERN int PCApplyBAorAB(PC,PCSide,Vec,Vec,Vec);
88: EXTERN int PCApplyTranspose(PC,Vec,Vec);
89: EXTERN int PCApplyBAorABTranspose(PC,PCSide,Vec,Vec,Vec);
90: EXTERN int PCApplyRichardson(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,int);
91: EXTERN int PCApplyRichardsonExists(PC,PetscTruth*);
93: EXTERN int PCRegisterDestroy(void);
94: EXTERN int PCRegisterAll(char*);
95: extern PetscTruth PCRegisterAllCalled;
97: EXTERN int PCRegister(char*,char*,char*,int(*)(PC));
98: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
99: #define PCRegisterDynamic(a,b,c,d) PCRegister(a,b,c,0)
100: #else
101: #define PCRegisterDynamic(a,b,c,d) PCRegister(a,b,c,d)
102: #endif
104: EXTERN int PCDestroy(PC);
105: EXTERN int PCSetFromOptions(PC);
106: EXTERN int PCESISetFromOptions(PC);
107: EXTERN int PCGetType(PC,PCType*);
109: EXTERN int PCGetFactoredMatrix(PC,Mat*);
110: EXTERN int PCSetModifySubMatrices(PC,int(*)(PC,int,IS*,IS*,Mat*,void*),void*);
111: EXTERN int PCModifySubMatrices(PC,int,IS*,IS*,Mat*,void*);
113: EXTERN int PCSetOperators(PC,Mat,Mat,MatStructure);
114: EXTERN int PCGetOperators(PC,Mat*,Mat*,MatStructure*);
116: EXTERN int PCSetVector(PC,Vec);
117: EXTERN int PCGetVector(PC,Vec*);
118: EXTERN int PCView(PC,PetscViewer);
120: EXTERN int PCSetOptionsPrefix(PC,char*);
121: EXTERN int PCAppendOptionsPrefix(PC,char*);
122: EXTERN int PCGetOptionsPrefix(PC,char**);
124: EXTERN int PCNullSpaceAttach(PC,MatNullSpace);
126: EXTERN int PCComputeExplicitOperator(PC,Mat*);
128: /*
129: These are used to provide extra scaling of preconditioned
130: operator for time-stepping schemes like in PVODE
131: */
132: EXTERN int PCDiagonalScale(PC,PetscTruth*);
133: EXTERN int PCDiagonalScaleLeft(PC,Vec,Vec);
134: EXTERN int PCDiagonalScaleRight(PC,Vec,Vec);
135: EXTERN int PCDiagonalScaleSet(PC,Vec);
137: /* ------------- options specific to particular preconditioners --------- */
139: EXTERN int PCJacobiSetUseRowMax(PC);
140: EXTERN int PCSORSetSymmetric(PC,MatSORType);
141: EXTERN int PCSORSetOmega(PC,PetscReal);
142: EXTERN int PCSORSetIterations(PC,int,int);
144: EXTERN int PCEisenstatSetOmega(PC,PetscReal);
145: EXTERN int PCEisenstatNoDiagonalScaling(PC);
147: #define USE_PRECONDITIONER_MATRIX 0
148: #define USE_TRUE_MATRIX 1
149: EXTERN int PCBJacobiSetUseTrueLocal(PC);
150: EXTERN int PCBJacobiSetTotalBlocks(PC,int,int*);
151: EXTERN int PCBJacobiSetLocalBlocks(PC,int,int*);
153: EXTERN int PCSLESSetUseTrue(PC);
155: EXTERN int PCShellSetApply(PC,int (*)(void*,Vec,Vec),void*);
156: EXTERN int PCShellSetSetUp(PC,int (*)(void*));
157: EXTERN int PCShellSetApplyRichardson(PC,int (*)(void*,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,int),void*);
158: EXTERN int PCShellSetView(PC,int (*)(void*,PetscViewer));
159: EXTERN int PCShellSetName(PC,char*);
160: EXTERN int PCShellGetName(PC,char**);
162: EXTERN int PCLUSetMatOrdering(PC,MatOrderingType);
163: EXTERN int PCLUSetReuseOrdering(PC,PetscTruth);
164: EXTERN int PCLUSetReuseFill(PC,PetscTruth);
165: EXTERN int PCLUSetUseInPlace(PC);
166: EXTERN int PCLUSetFill(PC,PetscReal);
167: EXTERN int PCLUSetDamping(PC,PetscReal);
168: EXTERN int PCLUSetPivoting(PC,PetscReal);
169: EXTERN int PCLUSetPivotInBlocks(PC,PetscTruth);
171: EXTERN int PCCholeskySetMatOrdering(PC,MatOrderingType);
172: EXTERN int PCCholeskySetReuseOrdering(PC,PetscTruth);
173: EXTERN int PCCholeskySetReuseFill(PC,PetscTruth);
174: EXTERN int PCCholeskySetUseInPlace(PC);
175: EXTERN int PCCholeskySetFill(PC,PetscReal);
176: EXTERN int PCCholeskySetDamping(PC,PetscReal);
177: EXTERN int PCCholeskySetPivotInBlocks(PC,PetscTruth);
179: EXTERN int PCILUSetMatOrdering(PC,MatOrderingType);
180: EXTERN int PCILUSetUseInPlace(PC);
181: EXTERN int PCILUSetFill(PC,PetscReal);
182: EXTERN int PCILUSetLevels(PC,int);
183: EXTERN int PCILUSetReuseOrdering(PC,PetscTruth);
184: EXTERN int PCILUSetUseDropTolerance(PC,PetscReal,PetscReal,int);
185: EXTERN int PCILUDTSetReuseFill(PC,PetscTruth);
186: EXTERN int PCILUSetAllowDiagonalFill(PC);
187: EXTERN int PCILUSetDamping(PC,PetscReal);
188: EXTERN int PCILUSetSinglePrecisionSolves(PC,PetscTruth);
189: EXTERN int PCILUSetPivotInBlocks(PC,PetscTruth);
191: EXTERN int PCICCSetMatOrdering(PC,MatOrderingType);
192: EXTERN int PCICCSetFill(PC,PetscReal);
193: EXTERN int PCICCSetLevels(PC,int);
194: EXTERN int PCICCSetPivotInBlocks(PC,PetscTruth);
196: EXTERN int PCASMSetLocalSubdomains(PC,int,IS *);
197: EXTERN int PCASMSetTotalSubdomains(PC,int,IS *);
198: EXTERN int PCASMSetOverlap(PC,int);
199: /*E
200: PCASMType - Type of additive Schwarz method to use
202: $ PC_ASM_BASIC - symmetric version where residuals from the ghost points are used
203: $ and computed values in ghost regions are added together. Classical
204: $ standard additive Schwarz
205: $ PC_ASM_RESTRICT - residuals from ghost points are used but computed values in ghost
206: $ region are discarded. Default
207: $ PC_ASM_INTERPOLATE - residuals from ghost points are not used, computed values in ghost
208: $ region are added back in
209: $ PC_ASM_NONE - ghost point residuals are not used, computed ghost values are discarded
210: $ not very good.
212: Level: beginner
214: .seealso: PCASMSetType()
215: E*/
216: typedef enum {PC_ASM_BASIC = 3,PC_ASM_RESTRICT = 1,PC_ASM_INTERPOLATE = 2,PC_ASM_NONE = 0} PCASMType;
217: EXTERN int PCASMSetType(PC,PCASMType);
218: EXTERN int PCASMCreateSubdomains2D(int,int,int,int,int,int,int *,IS **);
219: EXTERN int PCASMSetUseInPlace(PC);
220: EXTERN int PCASMGetLocalSubdomains(PC,int*,IS**);
221: EXTERN int PCASMGetLocalSubmatrices(PC,int*,Mat**);
223: /*E
224: PCCompositeType - Determines how two or more preconditioner are composed
226: $ PC_COMPOSITE_ADDITIVE - results from application of all preconditioners are added together
227: $ PC_COMPOSITE_MULTIPLICATIVE - preconditioners are applied sequentially to the residual freshly
228: $ computed after the previous preconditioner application
229: $ PC_COMPOSITE_SPECIAL - This is very special for a matrix of the form alpha I + R + S
230: $ where first preconditioner is built from alpha I + S and second from
231: $ alpha I + R
233: Level: beginner
235: .seealso: PCCompositeSetType()
236: E*/
237: typedef enum {PC_COMPOSITE_ADDITIVE,PC_COMPOSITE_MULTIPLICATIVE,PC_COMPOSITE_SPECIAL} PCCompositeType;
238: EXTERN int PCCompositeSetUseTrue(PC);
239: EXTERN int PCCompositeSetType(PC,PCCompositeType);
240: EXTERN int PCCompositeAddPC(PC,PCType);
241: EXTERN int PCCompositeGetPC(PC pc,int n,PC *);
242: EXTERN int PCCompositeSpecialSetAlpha(PC,PetscScalar);
244: EXTERN int PCRedundantSetScatter(PC,VecScatter,VecScatter);
245: EXTERN int PCRedundantGetOperators(PC,Mat*,Mat*);
246: EXTERN int PCRedundantGetPC(PC,PC*);
247: EXTERN int MatGetOrderingList(PetscFList *list);
249: EXTERN int PCMultiLevelSetFields(PC, int, int);
250: EXTERN int PCMultiLevelSetNonlinearIterate(PC, Vec);
251: EXTERN int PCMultiLevelSetGradientOperator(PC, int, int, PetscScalar);
252: EXTERN int PCMultiLevelApplyGradient(PC, Vec, Vec);
253: EXTERN int PCMultiLevelApplyGradientTrans(PC, Vec, Vec);
254: EXTERN int PCMultiLevelBuildSolution(PC, Vec);
255: EXTERN int PCMultiLevelGetMultiplier(PC, Vec, Vec);
257: EXTERN int PCSchurSetGradientOperator(PC, int, int);
258: EXTERN int PCSchurGetIterationNumber(PC, int *, int *);
260: EXTERN int PCSPAISetEpsilon(PC,double);
261: EXTERN int PCSPAISetNBSteps(PC,int);
262: EXTERN int PCSPAISetMax(PC,int);
263: EXTERN int PCSPAISetMaxNew(PC,int);
264: EXTERN int PCSPAISetBlockSize(PC,int);
265: EXTERN int PCSPAISetCacheSize(PC,int);
266: EXTERN int PCSPAISetVerbose(PC,int);
267: EXTERN int PCSPAISetSp(PC,int);
273: #endif /* __PETSCPC_H */