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 */