1 /* 2 Preconditioner module. 3 */ 4 #if !defined(__PETSCPC_H) 5 #define __PETSCPC_H 6 #include "petscdm.h" 7 PETSC_EXTERN_CXX_BEGIN 8 9 extern PetscErrorCode PCInitializePackage(const char[]); 10 11 /* 12 PCList contains the list of preconditioners currently registered 13 These are added with the PCRegisterDynamic() macro 14 */ 15 extern PetscFList PCList; 16 17 /*S 18 PC - Abstract PETSc object that manages all preconditioners 19 20 Level: beginner 21 22 Concepts: preconditioners 23 24 .seealso: PCCreate(), PCSetType(), PCType (for list of available types) 25 S*/ 26 typedef struct _p_PC* PC; 27 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() 32 33 Level: beginner 34 35 Notes: Click on the links below to see details on a particular solver 36 37 .seealso: PCSetType(), PC, PCCreate() 38 E*/ 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 PCFIELDSPLIT "fieldsplit" 62 #define PCTFS "tfs" 63 #define PCML "ml" 64 #define PCPROMETHEUS "prometheus" 65 #define PCGALERKIN "galerkin" 66 #define PCEXOTIC "exotic" 67 #define PCOPENMP "openmp" 68 #define PCSUPPORTGRAPH "supportgraph" 69 #define PCASA "asa" 70 #define PCCP "cp" 71 #define PCBFBT "bfbt" 72 #define PCLSC "lsc" 73 #define PCPYTHON "python" 74 #define PCPFMG "pfmg" 75 #define PCSYSPFMG "syspfmg" 76 #define PCREDISTRIBUTE "redistribute" 77 #define PCSACUSP "sacusp" 78 #define PCSACUSPPOLY "sacusppoly" 79 #define PCBICGSTABCUSP "bicgstabcusp" 80 #define PCSVD "svd" 81 #define PCAINVCUSP "ainvcusp" 82 83 /* Logging support */ 84 extern PetscClassId PC_CLASSID; 85 86 /*E 87 PCSide - If the preconditioner is to be applied to the left, right 88 or symmetrically around the operator. 89 90 Level: beginner 91 92 .seealso: 93 E*/ 94 typedef enum { PC_LEFT,PC_RIGHT,PC_SYMMETRIC } PCSide; 95 extern const char *PCSides[]; 96 97 extern PetscErrorCode PCCreate(MPI_Comm,PC*); 98 extern PetscErrorCode PCSetType(PC,const PCType); 99 extern PetscErrorCode PCSetUp(PC); 100 extern PetscErrorCode PCSetUpOnBlocks(PC); 101 extern PetscErrorCode PCApply(PC,Vec,Vec); 102 extern PetscErrorCode PCApplySymmetricLeft(PC,Vec,Vec); 103 extern PetscErrorCode PCApplySymmetricRight(PC,Vec,Vec); 104 extern PetscErrorCode PCApplyBAorAB(PC,PCSide,Vec,Vec,Vec); 105 extern PetscErrorCode PCApplyTranspose(PC,Vec,Vec); 106 extern PetscErrorCode PCApplyTransposeExists(PC,PetscBool *); 107 extern PetscErrorCode PCApplyBAorABTranspose(PC,PCSide,Vec,Vec,Vec); 108 109 /*E 110 PCRichardsonConvergedReason - reason a PCApplyRichardson method terminates 111 112 Level: advanced 113 114 Notes: this must match finclude/petscpc.h and the KSPConvergedReason values in petscksp.h 115 116 .seealso: PCApplyRichardson() 117 E*/ 118 typedef enum { 119 PCRICHARDSON_CONVERGED_RTOL = 2, 120 PCRICHARDSON_CONVERGED_ATOL = 3, 121 PCRICHARDSON_CONVERGED_ITS = 4, 122 PCRICHARDSON_DIVERGED_DTOL = -4} PCRichardsonConvergedReason; 123 124 extern PetscErrorCode PCApplyRichardson(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool ,PetscInt*,PCRichardsonConvergedReason*); 125 extern PetscErrorCode PCApplyRichardsonExists(PC,PetscBool *); 126 extern PetscErrorCode PCSetInitialGuessNonzero(PC,PetscBool ); 127 128 extern PetscErrorCode PCRegisterDestroy(void); 129 extern PetscErrorCode PCRegisterAll(const char[]); 130 extern PetscBool PCRegisterAllCalled; 131 132 extern PetscErrorCode PCRegister(const char[],const char[],const char[],PetscErrorCode(*)(PC)); 133 134 /*MC 135 PCRegisterDynamic - Adds a method to the preconditioner package. 136 137 Synopsis: 138 PetscErrorCode PCRegisterDynamic(const char *name_solver,const char *path,const char *name_create,PetscErrorCode (*routine_create)(PC)) 139 140 Not collective 141 142 Input Parameters: 143 + name_solver - name of a new user-defined solver 144 . path - path (either absolute or relative) the library containing this solver 145 . name_create - name of routine to create method context 146 - routine_create - routine to create method context 147 148 Notes: 149 PCRegisterDynamic() may be called multiple times to add several user-defined preconditioners. 150 151 If dynamic libraries are used, then the fourth input argument (routine_create) 152 is ignored. 153 154 Sample usage: 155 .vb 156 PCRegisterDynamic("my_solver","/home/username/my_lib/lib/libO/solaris/mylib", 157 "MySolverCreate",MySolverCreate); 158 .ve 159 160 Then, your solver can be chosen with the procedural interface via 161 $ PCSetType(pc,"my_solver") 162 or at runtime via the option 163 $ -pc_type my_solver 164 165 Level: advanced 166 167 Notes: ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR}, or ${any environmental variable} 168 occuring in pathname will be replaced with appropriate values. 169 If your function is not being put into a shared library then use PCRegister() instead 170 171 .keywords: PC, register 172 173 .seealso: PCRegisterAll(), PCRegisterDestroy() 174 M*/ 175 #if defined(PETSC_USE_DYNAMIC_LIBRARIES) 176 #define PCRegisterDynamic(a,b,c,d) PCRegister(a,b,c,0) 177 #else 178 #define PCRegisterDynamic(a,b,c,d) PCRegister(a,b,c,d) 179 #endif 180 181 extern PetscErrorCode PCReset(PC); 182 extern PetscErrorCode PCDestroy(PC); 183 extern PetscErrorCode PCSetFromOptions(PC); 184 extern PetscErrorCode PCGetType(PC,const PCType*); 185 186 extern PetscErrorCode PCFactorGetMatrix(PC,Mat*); 187 extern PetscErrorCode PCSetModifySubMatrices(PC,PetscErrorCode(*)(PC,PetscInt,const IS[],const IS[],Mat[],void*),void*); 188 extern PetscErrorCode PCModifySubMatrices(PC,PetscInt,const IS[],const IS[],Mat[],void*); 189 190 extern PetscErrorCode PCSetOperators(PC,Mat,Mat,MatStructure); 191 extern PetscErrorCode PCGetOperators(PC,Mat*,Mat*,MatStructure*); 192 extern PetscErrorCode PCGetOperatorsSet(PC,PetscBool *,PetscBool *); 193 194 extern PetscErrorCode PCView(PC,PetscViewer); 195 196 extern PetscErrorCode PCSetOptionsPrefix(PC,const char[]); 197 extern PetscErrorCode PCAppendOptionsPrefix(PC,const char[]); 198 extern PetscErrorCode PCGetOptionsPrefix(PC,const char*[]); 199 200 extern PetscErrorCode PCComputeExplicitOperator(PC,Mat*); 201 202 /* 203 These are used to provide extra scaling of preconditioned 204 operator for time-stepping schemes like in SUNDIALS 205 */ 206 extern PetscErrorCode PCGetDiagonalScale(PC,PetscBool *); 207 extern PetscErrorCode PCDiagonalScaleLeft(PC,Vec,Vec); 208 extern PetscErrorCode PCDiagonalScaleRight(PC,Vec,Vec); 209 extern PetscErrorCode PCSetDiagonalScale(PC,Vec); 210 211 /* ------------- options specific to particular preconditioners --------- */ 212 213 extern PetscErrorCode PCJacobiSetUseRowMax(PC); 214 extern PetscErrorCode PCJacobiSetUseRowSum(PC); 215 extern PetscErrorCode PCJacobiSetUseAbs(PC); 216 extern PetscErrorCode PCSORSetSymmetric(PC,MatSORType); 217 extern PetscErrorCode PCSORSetOmega(PC,PetscReal); 218 extern PetscErrorCode PCSORSetIterations(PC,PetscInt,PetscInt); 219 220 extern PetscErrorCode PCEisenstatSetOmega(PC,PetscReal); 221 extern PetscErrorCode PCEisenstatNoDiagonalScaling(PC); 222 223 #define USE_PRECONDITIONER_MATRIX 0 224 #define USE_TRUE_MATRIX 1 225 extern PetscErrorCode PCBJacobiSetUseTrueLocal(PC); 226 extern PetscErrorCode PCBJacobiSetTotalBlocks(PC,PetscInt,const PetscInt[]); 227 extern PetscErrorCode PCBJacobiSetLocalBlocks(PC,PetscInt,const PetscInt[]); 228 229 extern PetscErrorCode PCKSPSetUseTrue(PC); 230 231 extern PetscErrorCode PCShellSetApply(PC,PetscErrorCode (*)(PC,Vec,Vec)); 232 extern PetscErrorCode PCShellSetApplyBA(PC,PetscErrorCode (*)(PC,PCSide,Vec,Vec,Vec)); 233 extern PetscErrorCode PCShellSetApplyTranspose(PC,PetscErrorCode (*)(PC,Vec,Vec)); 234 extern PetscErrorCode PCShellSetSetUp(PC,PetscErrorCode (*)(PC)); 235 extern PetscErrorCode PCShellSetApplyRichardson(PC,PetscErrorCode (*)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool ,PetscInt*,PCRichardsonConvergedReason*)); 236 extern PetscErrorCode PCShellSetView(PC,PetscErrorCode (*)(PC,PetscViewer)); 237 extern PetscErrorCode PCShellSetDestroy(PC,PetscErrorCode (*)(PC)); 238 extern PetscErrorCode PCShellGetContext(PC,void**); 239 extern PetscErrorCode PCShellSetContext(PC,void*); 240 extern PetscErrorCode PCShellSetName(PC,const char[]); 241 extern PetscErrorCode PCShellGetName(PC,char*[]); 242 243 extern PetscErrorCode PCFactorSetZeroPivot(PC,PetscReal); 244 245 extern PetscErrorCode PCFactorSetShiftType(PC,MatFactorShiftType); 246 extern PetscErrorCode PCFactorSetShiftAmount(PC,PetscReal); 247 248 extern PetscErrorCode PCFactorSetMatSolverPackage(PC,const MatSolverPackage); 249 extern PetscErrorCode PCFactorGetMatSolverPackage(PC,const MatSolverPackage*); 250 extern PetscErrorCode PCFactorSetUpMatSolverPackage(PC); 251 252 extern PetscErrorCode PCFactorSetFill(PC,PetscReal); 253 extern PetscErrorCode PCFactorSetColumnPivot(PC,PetscReal); 254 extern PetscErrorCode PCFactorReorderForNonzeroDiagonal(PC,PetscReal); 255 256 extern PetscErrorCode PCFactorSetMatOrderingType(PC,const MatOrderingType); 257 extern PetscErrorCode PCFactorSetReuseOrdering(PC,PetscBool ); 258 extern PetscErrorCode PCFactorSetReuseFill(PC,PetscBool ); 259 extern PetscErrorCode PCFactorSetUseInPlace(PC); 260 extern PetscErrorCode PCFactorSetAllowDiagonalFill(PC); 261 extern PetscErrorCode PCFactorSetPivotInBlocks(PC,PetscBool ); 262 263 extern PetscErrorCode PCFactorSetLevels(PC,PetscInt); 264 extern PetscErrorCode PCFactorSetDropTolerance(PC,PetscReal,PetscReal,PetscInt); 265 266 extern PetscErrorCode PCASMSetLocalSubdomains(PC,PetscInt,IS[],IS[]); 267 extern PetscErrorCode PCASMSetTotalSubdomains(PC,PetscInt,IS[],IS[]); 268 extern PetscErrorCode PCASMSetOverlap(PC,PetscInt); 269 extern PetscErrorCode PCASMSetSortIndices(PC,PetscBool ); 270 271 /*E 272 PCASMType - Type of additive Schwarz method to use 273 274 $ PC_ASM_BASIC - symmetric version where residuals from the ghost points are used 275 $ and computed values in ghost regions are added together. Classical 276 $ standard additive Schwarz 277 $ PC_ASM_RESTRICT - residuals from ghost points are used but computed values in ghost 278 $ region are discarded. Default 279 $ PC_ASM_INTERPOLATE - residuals from ghost points are not used, computed values in ghost 280 $ region are added back in 281 $ PC_ASM_NONE - ghost point residuals are not used, computed ghost values are discarded 282 $ not very good. 283 284 Level: beginner 285 286 .seealso: PCASMSetType() 287 E*/ 288 typedef enum {PC_ASM_BASIC = 3,PC_ASM_RESTRICT = 1,PC_ASM_INTERPOLATE = 2,PC_ASM_NONE = 0} PCASMType; 289 extern const char *PCASMTypes[]; 290 291 extern PetscErrorCode PCASMSetType(PC,PCASMType); 292 extern PetscErrorCode PCASMCreateSubdomains(Mat,PetscInt,IS*[]); 293 extern PetscErrorCode PCASMDestroySubdomains(PetscInt,IS[],IS[]); 294 extern PetscErrorCode PCASMCreateSubdomains2D(PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,IS**,IS**); 295 extern PetscErrorCode PCASMGetLocalSubdomains(PC,PetscInt*,IS*[],IS*[]); 296 extern PetscErrorCode PCASMGetLocalSubmatrices(PC,PetscInt*,Mat*[]); 297 298 /*E 299 PCGASMType - Type of generalized additive Schwarz method to use (differs from ASM in allowing multiple processors per domain) 300 301 $ PC_GASM_BASIC - symmetric version where residuals from the ghost points are used 302 $ and computed values in ghost regions are added together. Classical 303 $ standard additive Schwarz 304 $ PC_GASM_RESTRICT - residuals from ghost points are used but computed values in ghost 305 $ region are discarded. Default 306 $ PC_GASM_INTERPOLATE - residuals from ghost points are not used, computed values in ghost 307 $ region are added back in 308 $ PC_GASM_NONE - ghost point residuals are not used, computed ghost values are discarded 309 $ not very good. 310 311 Level: beginner 312 313 .seealso: PCGASMSetType() 314 E*/ 315 typedef enum {PC_GASM_BASIC = 3,PC_GASM_RESTRICT = 1,PC_GASM_INTERPOLATE = 2,PC_GASM_NONE = 0} PCGASMType; 316 extern const char *PCGASMTypes[]; 317 318 extern PetscErrorCode PCGASMSetLocalSubdomains(PC,PetscInt,IS[],IS[]); 319 extern PetscErrorCode PCGASMSetTotalSubdomains(PC,PetscInt); 320 extern PetscErrorCode PCGASMSetOverlap(PC,PetscInt); 321 extern PetscErrorCode PCGASMSetSortIndices(PC,PetscBool ); 322 323 extern PetscErrorCode PCGASMSetType(PC,PCGASMType); 324 extern PetscErrorCode PCGASMCreateSubdomains(Mat,PetscInt,IS*[]); 325 extern PetscErrorCode PCGASMDestroySubdomains(PetscInt,IS[],IS[]); 326 extern PetscErrorCode PCGASMCreateSubdomains2D(PC,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,IS**,IS**); 327 extern PetscErrorCode PCGASMGetLocalSubdomains(PC,PetscInt*,IS*[],IS*[]); 328 extern PetscErrorCode PCGASMGetLocalSubmatrices(PC,PetscInt*,Mat*[]); 329 330 /*E 331 PCCompositeType - Determines how two or more preconditioner are composed 332 333 $ PC_COMPOSITE_ADDITIVE - results from application of all preconditioners are added together 334 $ PC_COMPOSITE_MULTIPLICATIVE - preconditioners are applied sequentially to the residual freshly 335 $ computed after the previous preconditioner application 336 $ PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE - preconditioners are applied sequentially to the residual freshly 337 $ computed from first preconditioner to last and then back (Use only for symmetric matrices and preconditions) 338 $ PC_COMPOSITE_SPECIAL - This is very special for a matrix of the form alpha I + R + S 339 $ where first preconditioner is built from alpha I + S and second from 340 $ alpha I + R 341 342 Level: beginner 343 344 .seealso: PCCompositeSetType() 345 E*/ 346 typedef enum {PC_COMPOSITE_ADDITIVE,PC_COMPOSITE_MULTIPLICATIVE,PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE,PC_COMPOSITE_SPECIAL,PC_COMPOSITE_SCHUR} PCCompositeType; 347 extern const char *PCCompositeTypes[]; 348 349 extern PetscErrorCode PCCompositeSetUseTrue(PC); 350 extern PetscErrorCode PCCompositeSetType(PC,PCCompositeType); 351 extern PetscErrorCode PCCompositeAddPC(PC,PCType); 352 extern PetscErrorCode PCCompositeGetPC(PC,PetscInt,PC *); 353 extern PetscErrorCode PCCompositeSpecialSetAlpha(PC,PetscScalar); 354 355 extern PetscErrorCode PCRedundantSetNumber(PC,PetscInt); 356 extern PetscErrorCode PCRedundantSetScatter(PC,VecScatter,VecScatter); 357 extern PetscErrorCode PCRedundantGetOperators(PC,Mat*,Mat*); 358 359 extern PetscErrorCode PCSPAISetEpsilon(PC,double); 360 extern PetscErrorCode PCSPAISetNBSteps(PC,PetscInt); 361 extern PetscErrorCode PCSPAISetMax(PC,PetscInt); 362 extern PetscErrorCode PCSPAISetMaxNew(PC,PetscInt); 363 extern PetscErrorCode PCSPAISetBlockSize(PC,PetscInt); 364 extern PetscErrorCode PCSPAISetCacheSize(PC,PetscInt); 365 extern PetscErrorCode PCSPAISetVerbose(PC,PetscInt); 366 extern PetscErrorCode PCSPAISetSp(PC,PetscInt); 367 368 extern PetscErrorCode PCHYPRESetType(PC,const char[]); 369 extern PetscErrorCode PCHYPREGetType(PC,const char*[]); 370 extern PetscErrorCode PCBJacobiGetLocalBlocks(PC,PetscInt*,const PetscInt*[]); 371 extern PetscErrorCode PCBJacobiGetTotalBlocks(PC,PetscInt*,const PetscInt*[]); 372 373 extern PetscErrorCode PCFieldSplitSetFields(PC,const char[],PetscInt,const PetscInt*); 374 extern PetscErrorCode PCFieldSplitSetType(PC,PCCompositeType); 375 extern PetscErrorCode PCFieldSplitSetBlockSize(PC,PetscInt); 376 extern PetscErrorCode PCFieldSplitSetIS(PC,const char[],IS); 377 378 /*E 379 PCFieldSplitSchurPreType - Determines how to precondition Schur complement 380 381 Level: intermediate 382 383 .seealso: PCFieldSplitSchurPrecondition() 384 E*/ 385 typedef enum {PC_FIELDSPLIT_SCHUR_PRE_SELF,PC_FIELDSPLIT_SCHUR_PRE_DIAG,PC_FIELDSPLIT_SCHUR_PRE_USER} PCFieldSplitSchurPreType; 386 extern const char *const PCFieldSplitSchurPreTypes[]; 387 388 extern PetscErrorCode PCFieldSplitSchurPrecondition(PC,PCFieldSplitSchurPreType,Mat); 389 extern PetscErrorCode PCFieldSplitGetSchurBlocks(PC,Mat*,Mat*,Mat*,Mat*); 390 391 extern PetscErrorCode PCGalerkinSetRestriction(PC,Mat); 392 extern PetscErrorCode PCGalerkinSetInterpolation(PC,Mat); 393 394 extern PetscErrorCode PCSetCoordinates(PC,PetscInt,PetscReal*); 395 extern PetscErrorCode PCSASetVectors(PC,PetscInt,PetscReal *); 396 397 extern PetscErrorCode PCPythonSetType(PC,const char[]); 398 399 extern PetscErrorCode PCSetDM(PC,DM); 400 extern PetscErrorCode PCGetDM(PC,DM*); 401 402 extern PetscErrorCode PCBiCGStabCUSPSetTolerance(PC,PetscReal); 403 extern PetscErrorCode PCBiCGStabCUSPSetIterations(PC,PetscInt); 404 extern PetscErrorCode PCBiCGStabCUSPSetUseVerboseMonitor(PC,PetscBool); 405 406 extern PetscErrorCode PCAINVCUSPSetDropTolerance(PC,PetscReal); 407 extern PetscErrorCode PCAINVCUSPUseScaling(PC,PetscBool); 408 extern PetscErrorCode PCAINVCUSPSetNonzeros(PC,PetscInt); 409 extern PetscErrorCode PCAINVCUSPSetLinParameter(PC,PetscInt); 410 PETSC_EXTERN_CXX_END 411 412 #endif /* __PETSCPC_H */ 413