1 /* 2 Preconditioner module. 3 */ 4 #if !defined(__PETSCPC_H) 5 #define __PETSCPC_H 6 #include <petscdm.h> 7 8 PETSC_EXTERN PetscErrorCode PCInitializePackage(const char[]); 9 10 /* 11 PCList contains the list of preconditioners currently registered 12 These are added with the PCRegisterDynamic() macro 13 */ 14 PETSC_EXTERN PetscFunctionList PCList; 15 16 /*S 17 PC - Abstract PETSc object that manages all preconditioners 18 19 Level: beginner 20 21 Concepts: preconditioners 22 23 .seealso: PCCreate(), PCSetType(), PCType (for list of available types) 24 S*/ 25 typedef struct _p_PC* PC; 26 27 /*J 28 PCType - String with the name of a PETSc preconditioner method or the creation function 29 with an optional dynamic library name, for example 30 http://www.mcs.anl.gov/petsc/lib.a:mypccreate() 31 32 Level: beginner 33 34 Notes: Click on the links below to see details on a particular solver 35 36 .seealso: PCSetType(), PC, PCCreate() 37 J*/ 38 typedef const char* PCType; 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 PCGASM "gasm" 51 #define PCKSP "ksp" 52 #define PCCOMPOSITE "composite" 53 #define PCREDUNDANT "redundant" 54 #define PCSPAI "spai" 55 #define PCNN "nn" 56 #define PCCHOLESKY "cholesky" 57 #define PCPBJACOBI "pbjacobi" 58 #define PCMAT "mat" 59 #define PCHYPRE "hypre" 60 #define PCPARMS "parms" 61 #define PCFIELDSPLIT "fieldsplit" 62 #define PCTFS "tfs" 63 #define PCML "ml" 64 #define PCGALERKIN "galerkin" 65 #define PCEXOTIC "exotic" 66 #define PCHMPI "hmpi" 67 #define PCSUPPORTGRAPH "supportgraph" 68 #define PCASA "asa" 69 #define PCCP "cp" 70 #define PCBFBT "bfbt" 71 #define PCLSC "lsc" 72 #define PCPYTHON "python" 73 #define PCPFMG "pfmg" 74 #define PCSYSPFMG "syspfmg" 75 #define PCREDISTRIBUTE "redistribute" 76 #define PCSVD "svd" 77 #define PCGAMG "gamg" 78 #define PCSACUSP "sacusp" /* these four run on NVIDIA GPUs using CUSP */ 79 #define PCSACUSPPOLY "sacusppoly" 80 #define PCBICGSTABCUSP "bicgstabcusp" 81 #define PCAINVCUSP "ainvcusp" 82 #define PCBDDC "bddc" 83 84 /* Logging support */ 85 PETSC_EXTERN PetscClassId PC_CLASSID; 86 87 /*E 88 PCSide - If the preconditioner is to be applied to the left, right 89 or symmetrically around the operator. 90 91 Level: beginner 92 93 .seealso: 94 E*/ 95 typedef enum { PC_SIDE_DEFAULT=-1,PC_LEFT,PC_RIGHT,PC_SYMMETRIC} PCSide; 96 #define PC_SIDE_MAX (PC_SYMMETRIC + 1) 97 PETSC_EXTERN const char *const *const PCSides; 98 99 PETSC_EXTERN PetscErrorCode PCCreate(MPI_Comm,PC*); 100 PETSC_EXTERN PetscErrorCode PCSetType(PC,PCType); 101 PETSC_EXTERN PetscErrorCode PCSetUp(PC); 102 PETSC_EXTERN PetscErrorCode PCSetUpOnBlocks(PC); 103 PETSC_EXTERN PetscErrorCode PCApply(PC,Vec,Vec); 104 PETSC_EXTERN PetscErrorCode PCApplySymmetricLeft(PC,Vec,Vec); 105 PETSC_EXTERN PetscErrorCode PCApplySymmetricRight(PC,Vec,Vec); 106 PETSC_EXTERN PetscErrorCode PCApplyBAorAB(PC,PCSide,Vec,Vec,Vec); 107 PETSC_EXTERN PetscErrorCode PCApplyTranspose(PC,Vec,Vec); 108 PETSC_EXTERN PetscErrorCode PCApplyTransposeExists(PC,PetscBool *); 109 PETSC_EXTERN PetscErrorCode PCApplyBAorABTranspose(PC,PCSide,Vec,Vec,Vec); 110 111 #define PC_FILE_CLASSID 1211222 112 113 /*E 114 PCRichardsonConvergedReason - reason a PCApplyRichardson method terminates 115 116 Level: advanced 117 118 Notes: this must match finclude/petscpc.h and the KSPConvergedReason values in petscksp.h 119 120 .seealso: PCApplyRichardson() 121 E*/ 122 typedef enum { 123 PCRICHARDSON_CONVERGED_RTOL = 2, 124 PCRICHARDSON_CONVERGED_ATOL = 3, 125 PCRICHARDSON_CONVERGED_ITS = 4, 126 PCRICHARDSON_DIVERGED_DTOL = -4} PCRichardsonConvergedReason; 127 128 PETSC_EXTERN PetscErrorCode PCApplyRichardson(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool ,PetscInt*,PCRichardsonConvergedReason*); 129 PETSC_EXTERN PetscErrorCode PCApplyRichardsonExists(PC,PetscBool *); 130 PETSC_EXTERN PetscErrorCode PCSetInitialGuessNonzero(PC,PetscBool ); 131 132 PETSC_EXTERN PetscErrorCode PCRegisterDestroy(void); 133 PETSC_EXTERN PetscErrorCode PCRegisterAll(const char[]); 134 PETSC_EXTERN PetscBool PCRegisterAllCalled; 135 136 PETSC_EXTERN PetscErrorCode PCRegister(const char[],const char[],const char[],PetscErrorCode(*)(PC)); 137 138 /*MC 139 PCRegisterDynamic - Adds a method to the preconditioner package. 140 141 Synopsis: 142 #include "petscpc.h" 143 PetscErrorCode PCRegisterDynamic(const char *name_solver,const char *path,const char *name_create,PetscErrorCode (*routine_create)(PC)) 144 145 Not collective 146 147 Input Parameters: 148 + name_solver - name of a new user-defined solver 149 . path - path (either absolute or relative) the library containing this solver 150 . name_create - name of routine to create method context 151 - routine_create - routine to create method context 152 153 Notes: 154 PCRegisterDynamic() may be called multiple times to add several user-defined preconditioners. 155 156 If dynamic libraries are used, then the fourth input argument (routine_create) 157 is ignored. 158 159 Sample usage: 160 .vb 161 PCRegisterDynamic("my_solver","/home/username/my_lib/lib/libO/solaris/mylib", 162 "MySolverCreate",MySolverCreate); 163 .ve 164 165 Then, your solver can be chosen with the procedural interface via 166 $ PCSetType(pc,"my_solver") 167 or at runtime via the option 168 $ -pc_type my_solver 169 170 Level: advanced 171 172 Notes: ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR}, or ${any environmental variable} 173 occuring in pathname will be replaced with appropriate values. 174 If your function is not being put into a shared library then use PCRegister() instead 175 176 .keywords: PC, register 177 178 .seealso: PCRegisterAll(), PCRegisterDestroy() 179 M*/ 180 #if defined(PETSC_USE_DYNAMIC_LIBRARIES) 181 #define PCRegisterDynamic(a,b,c,d) PCRegister(a,b,c,0) 182 #else 183 #define PCRegisterDynamic(a,b,c,d) PCRegister(a,b,c,d) 184 #endif 185 186 PETSC_EXTERN PetscErrorCode PCReset(PC); 187 PETSC_EXTERN PetscErrorCode PCDestroy(PC*); 188 PETSC_EXTERN PetscErrorCode PCSetFromOptions(PC); 189 PETSC_EXTERN PetscErrorCode PCGetType(PC,PCType*); 190 191 PETSC_EXTERN PetscErrorCode PCFactorGetMatrix(PC,Mat*); 192 PETSC_EXTERN PetscErrorCode PCSetModifySubMatrices(PC,PetscErrorCode(*)(PC,PetscInt,const IS[],const IS[],Mat[],void*),void*); 193 PETSC_EXTERN PetscErrorCode PCModifySubMatrices(PC,PetscInt,const IS[],const IS[],Mat[],void*); 194 195 PETSC_EXTERN PetscErrorCode PCSetOperators(PC,Mat,Mat,MatStructure); 196 PETSC_EXTERN PetscErrorCode PCGetOperators(PC,Mat*,Mat*,MatStructure*); 197 PETSC_EXTERN PetscErrorCode PCGetOperatorsSet(PC,PetscBool *,PetscBool *); 198 199 PETSC_EXTERN PetscErrorCode PCView(PC,PetscViewer); 200 PETSC_EXTERN PetscErrorCode PCLoad(PC,PetscViewer); 201 202 PETSC_EXTERN PetscErrorCode PCSetOptionsPrefix(PC,const char[]); 203 PETSC_EXTERN PetscErrorCode PCAppendOptionsPrefix(PC,const char[]); 204 PETSC_EXTERN PetscErrorCode PCGetOptionsPrefix(PC,const char*[]); 205 206 PETSC_EXTERN PetscErrorCode PCComputeExplicitOperator(PC,Mat*); 207 208 /* 209 These are used to provide extra scaling of preconditioned 210 operator for time-stepping schemes like in SUNDIALS 211 */ 212 PETSC_EXTERN PetscErrorCode PCGetDiagonalScale(PC,PetscBool *); 213 PETSC_EXTERN PetscErrorCode PCDiagonalScaleLeft(PC,Vec,Vec); 214 PETSC_EXTERN PetscErrorCode PCDiagonalScaleRight(PC,Vec,Vec); 215 PETSC_EXTERN PetscErrorCode PCSetDiagonalScale(PC,Vec); 216 217 /* ------------- options specific to particular preconditioners --------- */ 218 219 PETSC_EXTERN PetscErrorCode PCJacobiSetUseRowMax(PC); 220 PETSC_EXTERN PetscErrorCode PCJacobiSetUseRowSum(PC); 221 PETSC_EXTERN PetscErrorCode PCJacobiSetUseAbs(PC); 222 PETSC_EXTERN PetscErrorCode PCSORSetSymmetric(PC,MatSORType); 223 PETSC_EXTERN PetscErrorCode PCSORSetOmega(PC,PetscReal); 224 PETSC_EXTERN PetscErrorCode PCSORSetIterations(PC,PetscInt,PetscInt); 225 226 PETSC_EXTERN PetscErrorCode PCEisenstatSetOmega(PC,PetscReal); 227 PETSC_EXTERN PetscErrorCode PCEisenstatNoDiagonalScaling(PC); 228 229 #define USE_PRECONDITIONER_MATRIX 0 230 #define USE_TRUE_MATRIX 1 231 PETSC_EXTERN PetscErrorCode PCBJacobiSetUseTrueLocal(PC); 232 PETSC_EXTERN PetscErrorCode PCBJacobiSetTotalBlocks(PC,PetscInt,const PetscInt[]); 233 PETSC_EXTERN PetscErrorCode PCBJacobiSetLocalBlocks(PC,PetscInt,const PetscInt[]); 234 235 PETSC_EXTERN PetscErrorCode PCKSPSetUseTrue(PC); 236 237 PETSC_EXTERN PetscErrorCode PCShellSetApply(PC,PetscErrorCode (*)(PC,Vec,Vec)); 238 PETSC_EXTERN PetscErrorCode PCShellSetApplyBA(PC,PetscErrorCode (*)(PC,PCSide,Vec,Vec,Vec)); 239 PETSC_EXTERN PetscErrorCode PCShellSetApplyTranspose(PC,PetscErrorCode (*)(PC,Vec,Vec)); 240 PETSC_EXTERN PetscErrorCode PCShellSetSetUp(PC,PetscErrorCode (*)(PC)); 241 PETSC_EXTERN PetscErrorCode PCShellSetApplyRichardson(PC,PetscErrorCode (*)(PC,Vec,Vec,Vec,PetscReal,PetscReal,PetscReal,PetscInt,PetscBool ,PetscInt*,PCRichardsonConvergedReason*)); 242 PETSC_EXTERN PetscErrorCode PCShellSetView(PC,PetscErrorCode (*)(PC,PetscViewer)); 243 PETSC_EXTERN PetscErrorCode PCShellSetDestroy(PC,PetscErrorCode (*)(PC)); 244 PETSC_EXTERN PetscErrorCode PCShellGetContext(PC,void**); 245 PETSC_EXTERN PetscErrorCode PCShellSetContext(PC,void*); 246 PETSC_EXTERN PetscErrorCode PCShellSetName(PC,const char[]); 247 PETSC_EXTERN PetscErrorCode PCShellGetName(PC,const char*[]); 248 249 PETSC_EXTERN PetscErrorCode PCFactorSetZeroPivot(PC,PetscReal); 250 251 PETSC_EXTERN PetscErrorCode PCFactorSetShiftType(PC,MatFactorShiftType); 252 PETSC_EXTERN PetscErrorCode PCFactorSetShiftAmount(PC,PetscReal); 253 254 PETSC_EXTERN PetscErrorCode PCFactorSetMatSolverPackage(PC,const MatSolverPackage); 255 PETSC_EXTERN PetscErrorCode PCFactorGetMatSolverPackage(PC,const MatSolverPackage*); 256 PETSC_EXTERN PetscErrorCode PCFactorSetUpMatSolverPackage(PC); 257 258 PETSC_EXTERN PetscErrorCode PCFactorSetFill(PC,PetscReal); 259 PETSC_EXTERN PetscErrorCode PCFactorSetColumnPivot(PC,PetscReal); 260 PETSC_EXTERN PetscErrorCode PCFactorReorderForNonzeroDiagonal(PC,PetscReal); 261 262 PETSC_EXTERN PetscErrorCode PCFactorSetMatOrderingType(PC,MatOrderingType); 263 PETSC_EXTERN PetscErrorCode PCFactorSetReuseOrdering(PC,PetscBool ); 264 PETSC_EXTERN PetscErrorCode PCFactorSetReuseFill(PC,PetscBool ); 265 PETSC_EXTERN PetscErrorCode PCFactorSetUseInPlace(PC); 266 PETSC_EXTERN PetscErrorCode PCFactorSetAllowDiagonalFill(PC); 267 PETSC_EXTERN PetscErrorCode PCFactorSetPivotInBlocks(PC,PetscBool ); 268 269 PETSC_EXTERN PetscErrorCode PCFactorSetLevels(PC,PetscInt); 270 PETSC_EXTERN PetscErrorCode PCFactorSetDropTolerance(PC,PetscReal,PetscReal,PetscInt); 271 272 PETSC_EXTERN PetscErrorCode PCASMSetLocalSubdomains(PC,PetscInt,IS[],IS[]); 273 PETSC_EXTERN PetscErrorCode PCASMSetTotalSubdomains(PC,PetscInt,IS[],IS[]); 274 PETSC_EXTERN PetscErrorCode PCASMSetOverlap(PC,PetscInt); 275 PETSC_EXTERN PetscErrorCode PCASMSetSortIndices(PC,PetscBool ); 276 277 /*E 278 PCASMType - Type of additive Schwarz method to use 279 280 $ PC_ASM_BASIC - Symmetric version where residuals from the ghost points are used 281 $ and computed values in ghost regions are added together. 282 $ Classical standard additive Schwarz. 283 $ PC_ASM_RESTRICT - Residuals from ghost points are used but computed values in ghost 284 $ region are discarded. 285 $ Default. 286 $ PC_ASM_INTERPOLATE - Residuals from ghost points are not used, computed values in ghost 287 $ region are added back in. 288 $ PC_ASM_NONE - Residuals from ghost points are not used, computed ghost values are 289 $ discarded. 290 $ Not very good. 291 292 Level: beginner 293 294 .seealso: PCASMSetType() 295 E*/ 296 typedef enum {PC_ASM_BASIC = 3,PC_ASM_RESTRICT = 1,PC_ASM_INTERPOLATE = 2,PC_ASM_NONE = 0} PCASMType; 297 PETSC_EXTERN const char *const PCASMTypes[]; 298 299 PETSC_EXTERN PetscErrorCode PCASMSetType(PC,PCASMType); 300 PETSC_EXTERN PetscErrorCode PCASMCreateSubdomains(Mat,PetscInt,IS*[]); 301 PETSC_EXTERN PetscErrorCode PCASMDestroySubdomains(PetscInt,IS[],IS[]); 302 PETSC_EXTERN PetscErrorCode PCASMCreateSubdomains2D(PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,IS**,IS**); 303 PETSC_EXTERN PetscErrorCode PCASMGetLocalSubdomains(PC,PetscInt*,IS*[],IS*[]); 304 PETSC_EXTERN PetscErrorCode PCASMGetLocalSubmatrices(PC,PetscInt*,Mat*[]); 305 306 /*E 307 PCGASMType - Type of generalized additive Schwarz method to use (differs from ASM in allowing multiple processors per subdomain). 308 309 Each subdomain has nested inner and outer parts. The inner subdomains are assumed to form a non-overlapping covering of the computational 310 domain, while the outer subdomains contain the inner subdomains and overlap with each other. This preconditioner will compute 311 a subdomain correction over each *outer* subdomain from a residual computed there, but its different variants will differ in 312 (a) how the outer subdomain residual is computed, and (b) how the outer subdomain correction is computed. 313 314 $ PC_GASM_BASIC - Symmetric version where the full from the outer subdomain is used, and the resulting correction is applied 315 $ over the outer subdomains. As a result, points in the overlap will receive the sum of the corrections 316 $ from neighboring subdomains. 317 $ Classical standard additive Schwarz. 318 $ PC_GASM_RESTRICT - Residual from the outer subdomain is used but the correction is restricted to the inner subdomain only 319 $ (i.e., zeroed out over the overlap portion of the outer subdomain before being applied). As a result, 320 $ each point will receive a correction only from the unique inner subdomain containing it (nonoverlapping covering 321 $ assumption). 322 $ Default. 323 $ PC_GASM_INTERPOLATE - Residual is zeroed out over the overlap portion of the outer subdomain, but the resulting correction is 324 $ applied over the outer subdomain. As a result, points in the overlap will receive the sum of the corrections 325 $ from neighboring subdomains. 326 $ 327 $ PC_GASM_NONE - Residuals and corrections are zeroed out outside the local subdomains. 328 $ Not very good. 329 330 Level: beginner 331 332 .seealso: PCGASMSetType() 333 E*/ 334 typedef enum {PC_GASM_BASIC = 3,PC_GASM_RESTRICT = 1,PC_GASM_INTERPOLATE = 2,PC_GASM_NONE = 0} PCGASMType; 335 PETSC_EXTERN const char *const PCGASMTypes[]; 336 337 PETSC_EXTERN PetscErrorCode PCGASMSetSubdomains(PC,PetscInt,IS[],IS[]); 338 PETSC_EXTERN PetscErrorCode PCGASMSetTotalSubdomains(PC,PetscInt,PetscBool); 339 PETSC_EXTERN PetscErrorCode PCGASMSetOverlap(PC,PetscInt); 340 PETSC_EXTERN PetscErrorCode PCGASMSetSortIndices(PC,PetscBool ); 341 342 PETSC_EXTERN PetscErrorCode PCGASMSetType(PC,PCGASMType); 343 PETSC_EXTERN PetscErrorCode PCGASMCreateLocalSubdomains(Mat,PetscInt,PetscInt,IS*[],IS*[]); 344 PETSC_EXTERN PetscErrorCode PCGASMDestroySubdomains(PetscInt,IS[],IS[]); 345 PETSC_EXTERN PetscErrorCode PCGASMCreateSubdomains2D(PC,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt*,IS**,IS**); 346 PETSC_EXTERN PetscErrorCode PCGASMGetSubdomains(PC,PetscInt*,IS*[],IS*[]); 347 PETSC_EXTERN PetscErrorCode PCGASMGetSubmatrices(PC,PetscInt*,Mat*[]); 348 349 /*E 350 PCCompositeType - Determines how two or more preconditioner are composed 351 352 $ PC_COMPOSITE_ADDITIVE - results from application of all preconditioners are added together 353 $ PC_COMPOSITE_MULTIPLICATIVE - preconditioners are applied sequentially to the residual freshly 354 $ computed after the previous preconditioner application 355 $ PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE - preconditioners are applied sequentially to the residual freshly 356 $ computed from first preconditioner to last and then back (Use only for symmetric matrices and preconditions) 357 $ PC_COMPOSITE_SPECIAL - This is very special for a matrix of the form alpha I + R + S 358 $ where first preconditioner is built from alpha I + S and second from 359 $ alpha I + R 360 361 Level: beginner 362 363 .seealso: PCCompositeSetType() 364 E*/ 365 typedef enum {PC_COMPOSITE_ADDITIVE,PC_COMPOSITE_MULTIPLICATIVE,PC_COMPOSITE_SYMMETRIC_MULTIPLICATIVE,PC_COMPOSITE_SPECIAL,PC_COMPOSITE_SCHUR} PCCompositeType; 366 PETSC_EXTERN const char *const PCCompositeTypes[]; 367 368 PETSC_EXTERN PetscErrorCode PCCompositeSetUseTrue(PC); 369 PETSC_EXTERN PetscErrorCode PCCompositeSetType(PC,PCCompositeType); 370 PETSC_EXTERN PetscErrorCode PCCompositeAddPC(PC,PCType); 371 PETSC_EXTERN PetscErrorCode PCCompositeGetPC(PC,PetscInt,PC *); 372 PETSC_EXTERN PetscErrorCode PCCompositeSpecialSetAlpha(PC,PetscScalar); 373 374 PETSC_EXTERN PetscErrorCode PCRedundantSetNumber(PC,PetscInt); 375 PETSC_EXTERN PetscErrorCode PCRedundantSetScatter(PC,VecScatter,VecScatter); 376 PETSC_EXTERN PetscErrorCode PCRedundantGetOperators(PC,Mat*,Mat*); 377 378 PETSC_EXTERN PetscErrorCode PCSPAISetEpsilon(PC,double); 379 PETSC_EXTERN PetscErrorCode PCSPAISetNBSteps(PC,PetscInt); 380 PETSC_EXTERN PetscErrorCode PCSPAISetMax(PC,PetscInt); 381 PETSC_EXTERN PetscErrorCode PCSPAISetMaxNew(PC,PetscInt); 382 PETSC_EXTERN PetscErrorCode PCSPAISetBlockSize(PC,PetscInt); 383 PETSC_EXTERN PetscErrorCode PCSPAISetCacheSize(PC,PetscInt); 384 PETSC_EXTERN PetscErrorCode PCSPAISetVerbose(PC,PetscInt); 385 PETSC_EXTERN PetscErrorCode PCSPAISetSp(PC,PetscInt); 386 387 PETSC_EXTERN PetscErrorCode PCHYPRESetType(PC,const char[]); 388 PETSC_EXTERN PetscErrorCode PCHYPREGetType(PC,const char*[]); 389 PETSC_EXTERN PetscErrorCode PCBJacobiGetLocalBlocks(PC,PetscInt*,const PetscInt*[]); 390 PETSC_EXTERN PetscErrorCode PCBJacobiGetTotalBlocks(PC,PetscInt*,const PetscInt*[]); 391 392 PETSC_EXTERN PetscErrorCode PCFieldSplitSetFields(PC,const char[],PetscInt,const PetscInt*,const PetscInt*); 393 PETSC_EXTERN PetscErrorCode PCFieldSplitGetType(PC,PCCompositeType*); 394 PETSC_EXTERN PetscErrorCode PCFieldSplitSetType(PC,PCCompositeType); 395 PETSC_EXTERN PetscErrorCode PCFieldSplitSetBlockSize(PC,PetscInt); 396 PETSC_EXTERN PetscErrorCode PCFieldSplitSetIS(PC,const char[],IS); 397 PETSC_EXTERN PetscErrorCode PCFieldSplitGetIS(PC,const char[],IS*); 398 399 /*E 400 PCFieldSplitSchurPreType - Determines how to precondition Schur complement 401 402 Level: intermediate 403 404 .seealso: PCFieldSplitSchurPrecondition() 405 E*/ 406 typedef enum {PC_FIELDSPLIT_SCHUR_PRE_SELF,PC_FIELDSPLIT_SCHUR_PRE_A11,PC_FIELDSPLIT_SCHUR_PRE_USER} PCFieldSplitSchurPreType; 407 PETSC_EXTERN const char *const PCFieldSplitSchurPreTypes[]; 408 409 /*E 410 PCFieldSplitSchurFactType - determines which off-diagonal parts of the approximate block factorization to use 411 412 Level: intermediate 413 414 .seealso: PCFieldSplitSetSchurFactType() 415 E*/ 416 typedef enum { 417 PC_FIELDSPLIT_SCHUR_FACT_DIAG, 418 PC_FIELDSPLIT_SCHUR_FACT_LOWER, 419 PC_FIELDSPLIT_SCHUR_FACT_UPPER, 420 PC_FIELDSPLIT_SCHUR_FACT_FULL 421 } PCFieldSplitSchurFactType; 422 PETSC_EXTERN const char *const PCFieldSplitSchurFactTypes[]; 423 424 PETSC_EXTERN PetscErrorCode PCFieldSplitSchurPrecondition(PC,PCFieldSplitSchurPreType,Mat); 425 PETSC_EXTERN PetscErrorCode PCFieldSplitSetSchurFactType(PC,PCFieldSplitSchurFactType); 426 PETSC_EXTERN PetscErrorCode PCFieldSplitGetSchurBlocks(PC,Mat*,Mat*,Mat*,Mat*); 427 428 PETSC_EXTERN PetscErrorCode PCGalerkinSetRestriction(PC,Mat); 429 PETSC_EXTERN PetscErrorCode PCGalerkinSetInterpolation(PC,Mat); 430 431 PETSC_EXTERN PetscErrorCode PCSetCoordinates(PC,PetscInt,PetscInt,PetscReal*); 432 PETSC_EXTERN PetscErrorCode PCSASetVectors(PC,PetscInt,PetscReal *); 433 434 PETSC_EXTERN PetscErrorCode PCPythonSetType(PC,const char[]); 435 436 PETSC_EXTERN PetscErrorCode PCSetDM(PC,DM); 437 PETSC_EXTERN PetscErrorCode PCGetDM(PC,DM*); 438 439 PETSC_EXTERN PetscErrorCode PCSetApplicationContext(PC,void*); 440 PETSC_EXTERN PetscErrorCode PCGetApplicationContext(PC,void*); 441 442 PETSC_EXTERN PetscErrorCode PCBiCGStabCUSPSetTolerance(PC,PetscReal); 443 PETSC_EXTERN PetscErrorCode PCBiCGStabCUSPSetIterations(PC,PetscInt); 444 PETSC_EXTERN PetscErrorCode PCBiCGStabCUSPSetUseVerboseMonitor(PC,PetscBool); 445 446 PETSC_EXTERN PetscErrorCode PCAINVCUSPSetDropTolerance(PC,PetscReal); 447 PETSC_EXTERN PetscErrorCode PCAINVCUSPUseScaling(PC,PetscBool); 448 PETSC_EXTERN PetscErrorCode PCAINVCUSPSetNonzeros(PC,PetscInt); 449 PETSC_EXTERN PetscErrorCode PCAINVCUSPSetLinParameter(PC,PetscInt); 450 /*E 451 PCPARMSGlobalType - Determines the global preconditioner method in PARMS 452 453 Level: intermediate 454 455 .seealso: PCPARMSSetGlobal() 456 E*/ 457 typedef enum {PC_PARMS_GLOBAL_RAS,PC_PARMS_GLOBAL_SCHUR,PC_PARMS_GLOBAL_BJ} PCPARMSGlobalType; 458 PETSC_EXTERN const char *const PCPARMSGlobalTypes[]; 459 /*E 460 PCPARMSLocalType - Determines the local preconditioner method in PARMS 461 462 Level: intermediate 463 464 .seealso: PCPARMSSetLocal() 465 E*/ 466 typedef enum {PC_PARMS_LOCAL_ILU0,PC_PARMS_LOCAL_ILUK,PC_PARMS_LOCAL_ILUT,PC_PARMS_LOCAL_ARMS} PCPARMSLocalType; 467 PETSC_EXTERN const char *const PCPARMSLocalTypes[]; 468 469 PETSC_EXTERN PetscErrorCode PCPARMSSetGlobal(PC pc,PCPARMSGlobalType type); 470 PETSC_EXTERN PetscErrorCode PCPARMSSetLocal(PC pc,PCPARMSLocalType type); 471 PETSC_EXTERN PetscErrorCode PCPARMSSetSolveTolerances(PC pc,PetscReal tol,PetscInt maxits); 472 PETSC_EXTERN PetscErrorCode PCPARMSSetSolveRestart(PC pc,PetscInt restart); 473 PETSC_EXTERN PetscErrorCode PCPARMSSetNonsymPerm(PC pc,PetscBool nonsym); 474 PETSC_EXTERN PetscErrorCode PCPARMSSetFill(PC pc,PetscInt lfil0,PetscInt lfil1,PetscInt lfil2); 475 476 PETSC_EXTERN PetscErrorCode PCGAMGSetProcEqLim(PC,PetscInt); 477 PETSC_EXTERN PetscErrorCode PCGAMGSetRepartitioning(PC,PetscBool); 478 PETSC_EXTERN PetscErrorCode PCGAMGSetUseASMAggs(PC,PetscBool); 479 PETSC_EXTERN PetscErrorCode PCGAMGSetSolverType(PC,char[],PetscInt); 480 PETSC_EXTERN PetscErrorCode PCGAMGSetThreshold(PC,PetscReal); 481 PETSC_EXTERN PetscErrorCode PCGAMGSetCoarseEqLim(PC,PetscInt); 482 PETSC_EXTERN PetscErrorCode PCGAMGSetNlevels(PC,PetscInt); 483 typedef const char* PCGAMGType; 484 PETSC_EXTERN PetscErrorCode PCGAMGSetType( PC,PCGAMGType ); 485 PETSC_EXTERN PetscErrorCode PCGAMGSetNSmooths(PC pc, PetscInt n); 486 PETSC_EXTERN PetscErrorCode PCGAMGSetSymGraph(PC pc, PetscBool n); 487 PETSC_EXTERN PetscErrorCode PCGAMGSetSquareGraph(PC,PetscBool); 488 PETSC_EXTERN PetscErrorCode PCGAMGSetReuseProl(PC,PetscBool); 489 490 #if defined(PETSC_HAVE_PCBDDC) 491 /* Enum defining how to treat the coarse problem */ 492 typedef enum {SEQUENTIAL_BDDC,REPLICATED_BDDC,PARALLEL_BDDC,MULTILEVEL_BDDC} CoarseProblemType; 493 PETSC_EXTERN PetscErrorCode PCBDDCSetCoarseningRatio(PC,PetscInt); 494 PETSC_EXTERN PetscErrorCode PCBDDCSetMaxLevels(PC,PetscInt); 495 PETSC_EXTERN PetscErrorCode PCBDDCSetNullSpace(PC,MatNullSpace); 496 PETSC_EXTERN PetscErrorCode PCBDDCSetDirichletBoundaries(PC,IS); 497 PETSC_EXTERN PetscErrorCode PCBDDCGetDirichletBoundaries(PC,IS*); 498 PETSC_EXTERN PetscErrorCode PCBDDCSetNeumannBoundaries(PC,IS); 499 PETSC_EXTERN PetscErrorCode PCBDDCGetNeumannBoundaries(PC,IS*); 500 PETSC_EXTERN PetscErrorCode PCBDDCSetCoarseProblemType(PC,CoarseProblemType); 501 PETSC_EXTERN PetscErrorCode PCBDDCSetDofsSplitting(PC,PetscInt,IS[]); 502 PETSC_EXTERN PetscErrorCode PCBDDCSetLocalAdjacencyGraph(PC,PetscInt,const PetscInt[],const PetscInt[],PetscCopyMode); 503 PETSC_EXTERN PetscErrorCode PCBDDCCreateFETIDPOperators(PC,Mat*,PC*); 504 PETSC_EXTERN PetscErrorCode PCBDDCMatFETIDPGetRHS(Mat,Vec,Vec); 505 PETSC_EXTERN PetscErrorCode PCBDDCMatFETIDPGetSolution(Mat,Vec,Vec); 506 #endif 507 508 PETSC_EXTERN PetscErrorCode PCISSetUseStiffnessScaling(PC,PetscBool); 509 PETSC_EXTERN PetscErrorCode PCISSetSubdomainScalingFactor(PC,PetscScalar); 510 PETSC_EXTERN PetscErrorCode PCISSetSubdomainDiagonalScaling(PC,Vec); 511 512 #endif /* __PETSCPC_H */ 513