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 (for list of available types)
 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:    Notes: Click on the links below to see details on a particular solver

 37: .seealso: PCSetType(), PC, PCCreate()
 38: E*/
 39: #define PCNONE      "none"
 40: #define PCJACOBI    "jacobi"
 41: #define PCSOR       "sor"
 42: #define PCLU        "lu"
 43: #define PCSHELL     "shell"
 44: #define PCBJACOBI   "bjacobi"
 45: #define PCMG        "mg"
 46: #define PCEISENSTAT "eisenstat"
 47: #define PCILU       "ilu"
 48: #define PCICC       "icc"
 49: #define PCASM       "asm"
 50: #define PCSLES      "sles"
 51: #define PCCOMPOSITE "composite"
 52: #define PCREDUNDANT "redundant"
 53: #define PCSPAI      "spai"
 54: #define PCMILU      "milu"
 55: #define PCNN        "nn"
 56: #define PCCHOLESKY  "cholesky"
 57: #define PCRAMG      "ramg"
 58: #define PCSAMG      "samg"
 59: #define PCPBJACOBI  "pbjacobi"
 60: #define PCMULTILEVEL "multilevel"
 61: #define PCSCHUR      "schur"
 62: #define PCESI        "esi"
 63: #define PCPETSCESI   "petscesi"
 64: #define PCMAT        "mat"
 65: #define PCHYPRE      "hypre"

 67: /* Logging support */
 68: extern int PC_COOKIE;
 69: extern int PC_SetUp, PC_SetUpOnBlocks, PC_Apply, PC_ApplyCoarse, PC_ApplyMultiple, PC_ApplySymmetricLeft;
 70: extern int PC_ApplySymmetricRight, PC_ModifySubMatrices;

 72: /*E
 73:     PCSide - If the preconditioner is to be applied to the left, right
 74:      or symmetrically around the operator.

 76:    Level: beginner

 78: .seealso: 
 79: E*/
 80: typedef enum { PC_LEFT,PC_RIGHT,PC_SYMMETRIC } PCSide;

 82: EXTERN int PCCreate(MPI_Comm,PC*);
 83: EXTERN int PCSetType(PC,PCType);
 84: EXTERN int PCSetUp(PC);
 85: EXTERN int PCSetUpOnBlocks(PC);
 86: EXTERN int PCApply(PC,Vec,Vec,PCSide);
 87: EXTERN int PCApplySymmetricLeft(PC,Vec,Vec);
 88: EXTERN int PCApplySymmetricRight(PC,Vec,Vec);
 89: EXTERN int PCApplyBAorAB(PC,PCSide,Vec,Vec,Vec);
 90: EXTERN int PCApplyTranspose(PC,Vec,Vec);
 91: EXTERN int PCApplyBAorABTranspose(PC,PCSide,Vec,Vec,Vec);
 92: EXTERN int PCApplyRichardson(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,int);
 93: EXTERN int PCApplyRichardsonExists(PC,PetscTruth*);

 95: EXTERN int        PCRegisterDestroy(void);
 96: EXTERN int        PCRegisterAll(const char[]);
 97: extern PetscTruth PCRegisterAllCalled;

 99: EXTERN int PCRegister(const char[],const char[],const char[],int(*)(PC));

101: /*MC
102:    PCRegisterDynamic - Adds a method to the preconditioner package.

104:    Synopsis:
105:    int PCRegisterDynamic(char *name_solver,char *path,char *name_create,int (*routine_create)(PC))

107:    Not collective

109:    Input Parameters:
110: +  name_solver - name of a new user-defined solver
111: .  path - path (either absolute or relative) the library containing this solver
112: .  name_create - name of routine to create method context
113: -  routine_create - routine to create method context

115:    Notes:
116:    PCRegisterDynamic() may be called multiple times to add several user-defined preconditioners.

118:    If dynamic libraries are used, then the fourth input argument (routine_create)
119:    is ignored.

121:    Sample usage:
122: .vb
123:    PCRegisterDynamic("my_solver","/home/username/my_lib/lib/libO/solaris/mylib",
124:               "MySolverCreate",MySolverCreate);
125: .ve

127:    Then, your solver can be chosen with the procedural interface via
128: $     PCSetType(pc,"my_solver")
129:    or at runtime via the option
130: $     -pc_type my_solver

132:    Level: advanced

134:    Notes: ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR}, ${BOPT}, or ${any environmental variable}
135:            occuring in pathname will be replaced with appropriate values.
136:          If your function is not being put into a shared library then use PCRegister() instead

138: .keywords: PC, register

140: .seealso: PCRegisterAll(), PCRegisterDestroy()
141: M*/
142: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
143: #define PCRegisterDynamic(a,b,c,d) PCRegister(a,b,c,0)
144: #else
145: #define PCRegisterDynamic(a,b,c,d) PCRegister(a,b,c,d)
146: #endif

148: EXTERN int PCDestroy(PC);
149: EXTERN int PCSetFromOptions(PC);
150: EXTERN int PCESISetFromOptions(PC);
151: EXTERN int PCGetType(PC,PCType*);

153: EXTERN int PCGetFactoredMatrix(PC,Mat*);
154: EXTERN int PCSetModifySubMatrices(PC,int(*)(PC,int,const IS[],const IS[],Mat[],void*),void*);
155: EXTERN int PCModifySubMatrices(PC,int,const IS[],const IS[],Mat[],void*);

157: EXTERN int PCSetOperators(PC,Mat,Mat,MatStructure);
158: EXTERN int PCGetOperators(PC,Mat*,Mat*,MatStructure*);

160: EXTERN int PCSetVector(PC,Vec);
161: EXTERN int PCGetVector(PC,Vec*);
162: EXTERN int PCView(PC,PetscViewer);

164: EXTERN int PCSetOptionsPrefix(PC,const char[]);
165: EXTERN int PCAppendOptionsPrefix(PC,const char[]);
166: EXTERN int PCGetOptionsPrefix(PC,char*[]);

168: EXTERN int PCNullSpaceAttach(PC,MatNullSpace);

170: EXTERN int PCComputeExplicitOperator(PC,Mat*);

172: /*
173:       These are used to provide extra scaling of preconditioned 
174:    operator for time-stepping schemes like in PVODE 
175: */
176: EXTERN int PCDiagonalScale(PC,PetscTruth*);
177: EXTERN int PCDiagonalScaleLeft(PC,Vec,Vec);
178: EXTERN int PCDiagonalScaleRight(PC,Vec,Vec);
179: EXTERN int PCDiagonalScaleSet(PC,Vec);

181: /* ------------- options specific to particular preconditioners --------- */

183: EXTERN int PCJacobiSetUseRowMax(PC);
184: EXTERN int PCSORSetSymmetric(PC,MatSORType);
185: EXTERN int PCSORSetOmega(PC,PetscReal);
186: EXTERN int PCSORSetIterations(PC,int,int);

188: EXTERN int PCEisenstatSetOmega(PC,PetscReal);
189: EXTERN int PCEisenstatNoDiagonalScaling(PC);

191: #define USE_PRECONDITIONER_MATRIX 0
192: #define USE_TRUE_MATRIX           1
193: EXTERN int PCBJacobiSetUseTrueLocal(PC);
194: EXTERN int PCBJacobiSetTotalBlocks(PC,int,const int[]);
195: EXTERN int PCBJacobiSetLocalBlocks(PC,int,const int[]);

197: EXTERN int PCSLESSetUseTrue(PC);

199: EXTERN int PCShellSetApply(PC,int (*)(void*,Vec,Vec),void*);
200: EXTERN int PCShellSetApplyTranspose(PC,int (*)(void*,Vec,Vec));
201: EXTERN int PCShellSetSetUp(PC,int (*)(void*));
202: EXTERN int PCShellSetApplyRichardson(PC,int (*)(void*,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,int),void*);
203: EXTERN int PCShellSetView(PC,int (*)(void*,PetscViewer));
204: EXTERN int PCShellSetName(PC,const char[]);
205: EXTERN int PCShellGetName(PC,char*[]);

207: EXTERN int PCLUSetMatOrdering(PC,MatOrderingType);
208: EXTERN int PCLUSetReuseOrdering(PC,PetscTruth);
209: EXTERN int PCLUSetReuseFill(PC,PetscTruth);
210: EXTERN int PCLUSetUseInPlace(PC);
211: EXTERN int PCLUSetFill(PC,PetscReal);
212: EXTERN int PCLUSetDamping(PC,PetscReal);
213: EXTERN int PCLUSetShift(PC,PetscTruth);
214: EXTERN int PCLUSetPivoting(PC,PetscReal);
215: EXTERN int PCLUSetPivotInBlocks(PC,PetscTruth);
216: EXTERN int PCLUSetZeroPivot(PC,PetscReal);

218: EXTERN int PCCholeskySetMatOrdering(PC,MatOrderingType);
219: EXTERN int PCCholeskySetReuseOrdering(PC,PetscTruth);
220: EXTERN int PCCholeskySetReuseFill(PC,PetscTruth);
221: EXTERN int PCCholeskySetUseInPlace(PC);
222: EXTERN int PCCholeskySetFill(PC,PetscReal);
223: EXTERN int PCCholeskySetDamping(PC,PetscReal);
224: EXTERN int PCCholeskySetShift(PC,PetscTruth);
225: EXTERN int PCCholeskySetPivotInBlocks(PC,PetscTruth);

227: EXTERN int PCILUSetMatOrdering(PC,MatOrderingType);
228: EXTERN int PCILUSetUseInPlace(PC);
229: EXTERN int PCILUSetFill(PC,PetscReal);
230: EXTERN int PCILUSetLevels(PC,int);
231: EXTERN int PCILUSetReuseOrdering(PC,PetscTruth);
232: EXTERN int PCILUSetUseDropTolerance(PC,PetscReal,PetscReal,int);
233: EXTERN int PCILUDTSetReuseFill(PC,PetscTruth);
234: EXTERN int PCILUSetAllowDiagonalFill(PC);
235: EXTERN int PCILUSetDamping(PC,PetscReal);
236: EXTERN int PCILUSetShift(PC,PetscTruth);
237: EXTERN int PCILUSetPivotInBlocks(PC,PetscTruth);
238: EXTERN int PCILUSetZeroPivot(PC,PetscReal);

240: EXTERN int PCICCSetMatOrdering(PC,MatOrderingType);
241: EXTERN int PCICCSetFill(PC,PetscReal);
242: EXTERN int PCICCSetLevels(PC,int);
243: EXTERN int PCICCSetDamping(PC,PetscReal);
244: EXTERN int PCICCSetShift(PC,PetscTruth);
245: EXTERN int PCICCSetPivotInBlocks(PC,PetscTruth);
246: EXTERN int PCICCSetZeroPivot(PC,PetscReal);

248: EXTERN int PCASMSetLocalSubdomains(PC,int,IS[]);
249: EXTERN int PCASMSetTotalSubdomains(PC,int,IS[]);
250: EXTERN int PCASMSetOverlap(PC,int);
251: /*E
252:     PCASMType - Type of additive Schwarz method to use

254: $  PC_ASM_BASIC - symmetric version where residuals from the ghost points are used
255: $                 and computed values in ghost regions are added together. Classical
256: $                 standard additive Schwarz
257: $  PC_ASM_RESTRICT - residuals from ghost points are used but computed values in ghost
258: $                    region are discarded. Default
259: $  PC_ASM_INTERPOLATE - residuals from ghost points are not used, computed values in ghost
260: $                       region are added back in
261: $  PC_ASM_NONE - ghost point residuals are not used, computed ghost values are discarded
262: $                not very good.                

264:    Level: beginner

266: .seealso: PCASMSetType()
267: E*/
268: typedef enum {PC_ASM_BASIC = 3,PC_ASM_RESTRICT = 1,PC_ASM_INTERPOLATE = 2,PC_ASM_NONE = 0} PCASMType;
269: EXTERN int PCASMSetType(PC,PCASMType);
270: EXTERN int PCASMCreateSubdomains2D(int,int,int,int,int,int,int *,IS **);
271: EXTERN int PCASMSetUseInPlace(PC);
272: EXTERN int PCASMGetLocalSubdomains(PC,int*,IS*[]);
273: EXTERN int PCASMGetLocalSubmatrices(PC,int*,Mat*[]);

275: /*E
276:     PCCompositeType - Determines how two or more preconditioner are composed

278: $  PC_COMPOSITE_ADDITIVE - results from application of all preconditioners are added together
279: $  PC_COMPOSITE_MULTIPLICATIVE - preconditioners are applied sequentially to the residual freshly
280: $                                computed after the previous preconditioner application
281: $  PC_COMPOSITE_SPECIAL - This is very special for a matrix of the form alpha I + R + S
282: $                         where first preconditioner is built from alpha I + S and second from
283: $                         alpha I + R

285:    Level: beginner

287: .seealso: PCCompositeSetType()
288: E*/
289: typedef enum {PC_COMPOSITE_ADDITIVE,PC_COMPOSITE_MULTIPLICATIVE,PC_COMPOSITE_SPECIAL} PCCompositeType;
290: EXTERN int PCCompositeSetUseTrue(PC);
291: EXTERN int PCCompositeSetType(PC,PCCompositeType);
292: EXTERN int PCCompositeAddPC(PC,PCType);
293: EXTERN int PCCompositeGetPC(PC pc,int n,PC *);
294: EXTERN int PCCompositeSpecialSetAlpha(PC,PetscScalar);

296: EXTERN int PCRedundantSetScatter(PC,VecScatter,VecScatter);
297: EXTERN int PCRedundantGetOperators(PC,Mat*,Mat*);
298: EXTERN int PCRedundantGetPC(PC,PC*);
299: EXTERN int MatGetOrderingList(PetscFList *list);

301: EXTERN int PCMultiLevelSetFields(PC, int, int);
302: EXTERN int PCMultiLevelSetNonlinearIterate(PC, Vec);
303: EXTERN int PCMultiLevelSetGradientOperator(PC, int, int, PetscScalar);
304: EXTERN int PCMultiLevelApplyGradient(PC, Vec, Vec);
305: EXTERN int PCMultiLevelApplyGradientTrans(PC, Vec, Vec);
306: EXTERN int PCMultiLevelBuildSolution(PC, Vec);
307: EXTERN int PCMultiLevelGetMultiplier(PC, Vec, Vec);

309: EXTERN int PCSchurSetGradientOperator(PC, int, int);
310: EXTERN int PCSchurGetIterationNumber(PC, int *, int *);

312: EXTERN int PCSPAISetEpsilon(PC,double);
313: EXTERN int PCSPAISetNBSteps(PC,int);
314: EXTERN int PCSPAISetMax(PC,int);
315: EXTERN int PCSPAISetMaxNew(PC,int);
316: EXTERN int PCSPAISetBlockSize(PC,int);
317: EXTERN int PCSPAISetCacheSize(PC,int);
318: EXTERN int PCSPAISetVerbose(PC,int);
319: EXTERN int PCSPAISetSp(PC,int);

321: EXTERN int PCHYPRESetType(PC,const char[]);
322: EXTERN int PCBJacobiGetLocalBlocks(PC,int*,const int*[]);
323: EXTERN int PCBJacobiGetTotalBlocks(PC,int*,const int*[]);

325: #endif /* __PETSCPC_H */