1 /* 2 Defines the interface functions for the Krylov subspace accelerators. 3 */ 4 #ifndef __PETSCKSP_H 5 #define __PETSCKSP_H 6 #include "petscpc.h" 7 PETSC_EXTERN_CXX_BEGIN 8 9 EXTERN PetscErrorCode KSPInitializePackage(const char[]); 10 11 /*S 12 KSP - Abstract PETSc object that manages all Krylov methods 13 14 Level: beginner 15 16 Concepts: Krylov methods 17 18 .seealso: KSPCreate(), KSPSetType(), KSPType, SNES, TS, PC, KSP 19 S*/ 20 typedef struct _p_KSP* KSP; 21 22 /*E 23 KSPType - String with the name of a PETSc Krylov method or the creation function 24 with an optional dynamic library name, for example 25 http://www.mcs.anl.gov/petsc/lib.a:mykspcreate() 26 27 Level: beginner 28 29 .seealso: KSPSetType(), KSP 30 E*/ 31 #define KSPRICHARDSON "richardson" 32 #define KSPCHEBYCHEV "chebychev" 33 #define KSPCG "cg" 34 #define KSPCGNE "cgne" 35 #define KSPGMRES "gmres" 36 #define KSPTCQMR "tcqmr" 37 #define KSPBCGS "bcgs" 38 #define KSPCGS "cgs" 39 #define KSPTFQMR "tfqmr" 40 #define KSPCR "cr" 41 #define KSPLSQR "lsqr" 42 #define KSPPREONLY "preonly" 43 #define KSPQCG "qcg" 44 #define KSPBICG "bicg" 45 #define KSPFGMRES "fgmres" 46 #define KSPMINRES "minres" 47 #define KSPSYMMLQ "symmlq" 48 #define KSPLGMRES "lgmres" 49 #define KSPType char* 50 51 /* Logging support */ 52 extern PetscCookie KSP_COOKIE; 53 extern PetscEvent KSP_GMRESOrthogonalization, KSP_SetUp, KSP_Solve; 54 55 EXTERN PetscErrorCode KSPCreate(MPI_Comm,KSP *); 56 EXTERN PetscErrorCode KSPSetType(KSP,const KSPType); 57 EXTERN PetscErrorCode KSPSetUp(KSP); 58 EXTERN PetscErrorCode KSPSetUpOnBlocks(KSP); 59 EXTERN PetscErrorCode KSPSolve(KSP,Vec,Vec); 60 EXTERN PetscErrorCode KSPSolveTranspose(KSP,Vec,Vec); 61 EXTERN PetscErrorCode KSPDestroy(KSP); 62 63 extern PetscFList KSPList; 64 EXTERN PetscErrorCode KSPRegisterAll(const char[]); 65 EXTERN PetscErrorCode KSPRegisterDestroy(void); 66 67 EXTERN PetscErrorCode KSPRegister(const char[],const char[],const char[],PetscErrorCode (*)(KSP)); 68 69 /*MC 70 KSPRegisterDynamic - Adds a method to the Krylov subspace solver package. 71 72 Synopsis: 73 int KSPRegisterDynamic(char *name_solver,char *path,char *name_create,PetscErrorCode (*routine_create)(KSP)) 74 75 Not Collective 76 77 Input Parameters: 78 + name_solver - name of a new user-defined solver 79 . path - path (either absolute or relative) the library containing this solver 80 . name_create - name of routine to create method context 81 - routine_create - routine to create method context 82 83 Notes: 84 KSPRegisterDynamic() may be called multiple times to add several user-defined solvers. 85 86 If dynamic libraries are used, then the fourth input argument (routine_create) 87 is ignored. 88 89 Sample usage: 90 .vb 91 KSPRegisterDynamic("my_solver",/home/username/my_lib/lib/libO/solaris/mylib.a, 92 "MySolverCreate",MySolverCreate); 93 .ve 94 95 Then, your solver can be chosen with the procedural interface via 96 $ KSPSetType(ksp,"my_solver") 97 or at runtime via the option 98 $ -ksp_type my_solver 99 100 Level: advanced 101 102 Notes: Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR}, ${BOPT}, 103 and others of the form ${any_environmental_variable} occuring in pathname will be 104 replaced with appropriate values. 105 If your function is not being put into a shared library then use KSPRegister() instead 106 107 .keywords: KSP, register 108 109 .seealso: KSPRegisterAll(), KSPRegisterDestroy() 110 111 M*/ 112 #if defined(PETSC_USE_DYNAMIC_LIBRARIES) 113 #define KSPRegisterDynamic(a,b,c,d) KSPRegister(a,b,c,0) 114 #else 115 #define KSPRegisterDynamic(a,b,c,d) KSPRegister(a,b,c,d) 116 #endif 117 118 EXTERN PetscErrorCode KSPGetType(KSP,KSPType *); 119 EXTERN PetscErrorCode KSPSetPreconditionerSide(KSP,PCSide); 120 EXTERN PetscErrorCode KSPGetPreconditionerSide(KSP,PCSide*); 121 EXTERN PetscErrorCode KSPGetTolerances(KSP,PetscReal*,PetscReal*,PetscReal*,int*); 122 EXTERN PetscErrorCode KSPSetTolerances(KSP,PetscReal,PetscReal,PetscReal,int); 123 EXTERN PetscErrorCode KSPSetInitialGuessNonzero(KSP,PetscTruth); 124 EXTERN PetscErrorCode KSPGetInitialGuessNonzero(KSP,PetscTruth *); 125 EXTERN PetscErrorCode KSPSetInitialGuessKnoll(KSP,PetscTruth); 126 EXTERN PetscErrorCode KSPGetInitialGuessKnoll(KSP,PetscTruth*); 127 EXTERN PetscErrorCode KSPSetComputeEigenvalues(KSP,PetscTruth); 128 EXTERN PetscErrorCode KSPSetComputeSingularValues(KSP,PetscTruth); 129 EXTERN PetscErrorCode KSPGetRhs(KSP,Vec *); 130 EXTERN PetscErrorCode KSPGetSolution(KSP,Vec *); 131 EXTERN PetscErrorCode KSPGetResidualNorm(KSP,PetscReal*); 132 EXTERN PetscErrorCode KSPGetIterationNumber(KSP,int*); 133 EXTERN PetscErrorCode KSPSetNullSpace(KSP,MatNullSpace); 134 EXTERN PetscErrorCode KSPGetNullSpace(KSP,MatNullSpace*); 135 136 EXTERN PetscErrorCode KSPSetPC(KSP,PC); 137 EXTERN PetscErrorCode KSPGetPC(KSP,PC*); 138 139 EXTERN PetscErrorCode KSPSetMonitor(KSP,PetscErrorCode (*)(KSP,int,PetscReal,void*),void *,PetscErrorCode (*)(void*)); 140 EXTERN PetscErrorCode KSPClearMonitor(KSP); 141 EXTERN PetscErrorCode KSPGetMonitorContext(KSP,void **); 142 EXTERN PetscErrorCode KSPGetResidualHistory(KSP,PetscReal*[],int *); 143 EXTERN PetscErrorCode KSPSetResidualHistory(KSP,PetscReal[],int,PetscTruth); 144 145 /* not sure where to put this */ 146 EXTERN PetscErrorCode PCKSPGetKSP(PC,KSP*); 147 EXTERN PetscErrorCode PCBJacobiGetSubKSP(PC,int*,int*,KSP*[]); 148 EXTERN PetscErrorCode PCASMGetSubKSP(PC,int*,int*,KSP*[]); 149 150 EXTERN PetscErrorCode KSPBuildSolution(KSP,Vec,Vec *); 151 EXTERN PetscErrorCode KSPBuildResidual(KSP,Vec,Vec,Vec *); 152 153 EXTERN PetscErrorCode KSPRichardsonSetScale(KSP,PetscReal); 154 EXTERN PetscErrorCode KSPChebychevSetEigenvalues(KSP,PetscReal,PetscReal); 155 EXTERN PetscErrorCode KSPComputeExtremeSingularValues(KSP,PetscReal*,PetscReal*); 156 EXTERN PetscErrorCode KSPComputeEigenvalues(KSP,int,PetscReal*,PetscReal*,int *); 157 EXTERN PetscErrorCode KSPComputeEigenvaluesExplicitly(KSP,int,PetscReal*,PetscReal*); 158 159 EXTERN PetscErrorCode KSPGMRESSetRestart(KSP, int); 160 EXTERN PetscErrorCode KSPGMRESSetHapTol(KSP,PetscReal); 161 162 EXTERN PetscErrorCode KSPGMRESSetPreAllocateVectors(KSP); 163 EXTERN PetscErrorCode KSPGMRESSetOrthogonalization(KSP,PetscErrorCode (*)(KSP,int)); 164 EXTERN PetscErrorCode KSPGMRESModifiedGramSchmidtOrthogonalization(KSP,int); 165 EXTERN PetscErrorCode KSPGMRESClassicalGramSchmidtOrthogonalization(KSP,int); 166 167 EXTERN PetscErrorCode KSPLGMRESSetAugDim(KSP,int); 168 EXTERN PetscErrorCode KSPLGMRESSetConstant(KSP); 169 170 /*E 171 KSPGMRESCGSRefinementType - How the classical (unmodified) Gram-Schmidt is performed. 172 173 Level: advanced 174 175 .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), 176 KSPGMRESSetCGSRefinementType(), KSPGMRESModifiedGramSchmidtOrthogonalization() 177 178 E*/ 179 typedef enum {KSP_GMRES_CGS_REFINE_NEVER, KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS} KSPGMRESCGSRefinementType; 180 181 /*M 182 KSP_GMRES_CGS_REFINE_NEVER - Just do the classical (unmodified) Gram-Schmidt process 183 184 Level: advanced 185 186 Note: Possible unstable, but the fastest to compute 187 188 .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), 189 KSPGMRESSetCGSRefinementType(), KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS, 190 KSPGMRESModifiedGramSchmidtOrthogonalization() 191 M*/ 192 193 /*M 194 KSP_GMRES_CGS_REFINE_IFNEEDED - Do the classical (unmodified) Gram-Schmidt process and one step of 195 iterative refinement if an estimate of the orthogonality of the resulting vectors indicates 196 poor orthogonality. 197 198 Level: advanced 199 200 Note: This is slower than KSP_GMRES_CGS_REFINE_NEVER because it requires an extra norm computation to 201 estimate the orthogonality but is more stable. 202 203 .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), 204 KSPGMRESSetCGSRefinementType(), KSP_GMRES_CGS_REFINE_NEVER, KSP_GMRES_CGS_REFINE_ALWAYS, 205 KSPGMRESModifiedGramSchmidtOrthogonalization() 206 M*/ 207 208 /*M 209 KSP_GMRES_CGS_REFINE_NEVER - Do two steps of the classical (unmodified) Gram-Schmidt process. 210 211 Level: advanced 212 213 Note: This is roughly twice the cost of KSP_GMRES_CGS_REFINE_NEVER because it performs the process twice 214 but it saves the extra norm calculation needed by KSP_GMRES_CGS_REFINE_IFNEEDED. 215 216 You should only use this if you absolutely know that the iterative refinement is needed. 217 218 .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), 219 KSPGMRESSetCGSRefinementType(), KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS, 220 KSPGMRESModifiedGramSchmidtOrthogonalization() 221 M*/ 222 223 EXTERN PetscErrorCode KSPGMRESSetCGSRefinementType(KSP,KSPGMRESCGSRefinementType); 224 225 EXTERN PetscErrorCode KSPFGMRESModifyPCNoChange(KSP,int,int,PetscReal,void*); 226 EXTERN PetscErrorCode KSPFGMRESModifyPCKSP(KSP,int,int,PetscReal,void*); 227 EXTERN PetscErrorCode KSPFGMRESSetModifyPC(KSP,PetscErrorCode (*)(KSP,int,int,PetscReal,void*),void*,PetscErrorCode(*)(void*)); 228 229 EXTERN PetscErrorCode KSPQCGSetTrustRegionRadius(KSP,PetscReal); 230 EXTERN PetscErrorCode KSPQCGGetQuadratic(KSP,PetscReal*); 231 EXTERN PetscErrorCode KSPQCGGetTrialStepNorm(KSP,PetscReal*); 232 233 EXTERN PetscErrorCode KSPSetFromOptions(KSP); 234 EXTERN PetscErrorCode KSPAddOptionsChecker(PetscErrorCode (*)(KSP)); 235 236 EXTERN PetscErrorCode KSPSingularValueMonitor(KSP,int,PetscReal,void *); 237 EXTERN PetscErrorCode KSPDefaultMonitor(KSP,int,PetscReal,void *); 238 EXTERN PetscErrorCode KSPTrueMonitor(KSP,int,PetscReal,void *); 239 EXTERN PetscErrorCode KSPDefaultSMonitor(KSP,int,PetscReal,void *); 240 EXTERN PetscErrorCode KSPVecViewMonitor(KSP,int,PetscReal,void *); 241 EXTERN PetscErrorCode KSPGMRESKrylovMonitor(KSP,int,PetscReal,void *); 242 243 EXTERN PetscErrorCode KSPUnwindPreconditioner(KSP,Vec,Vec); 244 EXTERN PetscErrorCode KSPDefaultBuildSolution(KSP,Vec,Vec*); 245 EXTERN PetscErrorCode KSPDefaultBuildResidual(KSP,Vec,Vec,Vec *); 246 247 EXTERN PetscErrorCode KSPSetOperators(KSP,Mat,Mat,MatStructure); 248 EXTERN PetscErrorCode KSPSetOptionsPrefix(KSP,const char[]); 249 EXTERN PetscErrorCode KSPAppendOptionsPrefix(KSP,const char[]); 250 EXTERN PetscErrorCode KSPGetOptionsPrefix(KSP,char*[]); 251 252 EXTERN PetscErrorCode KSPSetDiagonalScale(KSP,PetscTruth); 253 EXTERN PetscErrorCode KSPGetDiagonalScale(KSP,PetscTruth*); 254 EXTERN PetscErrorCode KSPSetDiagonalScaleFix(KSP,PetscTruth); 255 EXTERN PetscErrorCode KSPGetDiagonalScaleFix(KSP,PetscTruth*); 256 257 EXTERN PetscErrorCode KSPView(KSP,PetscViewer); 258 259 /*E 260 KSPNormType - Norm that is passed in the Krylov convergence 261 test routines. 262 263 Level: advanced 264 265 Notes: this must match finclude/petscksp.h 266 267 .seealso: KSPSolve(), KSPGetConvergedReason(), KSPSetNormType(), 268 KSPSetConvergenceTest() 269 E*/ 270 typedef enum {KSP_NO_NORM = 0, 271 KSP_PRECONDITIONED_NORM = 1, 272 KSP_UNPRECONDITIONED_NORM = 2, 273 KSP_NATURAL_NORM = 3} KSPNormType; 274 275 /*M 276 KSP_NO_NORM - Do not compute a norm during the Krylov process. This will 277 possibly save some computation but means the convergence test cannot 278 be based on a norm of a residual etc. 279 280 Level: advanced 281 282 Note: Some Krylov methods need to compute a residual norm and then this is ignored 283 284 .seealso: KSPNormType, KSPSetNormType(), KSP_PRECONDITIONED_NORM, KSP_UNPRECONDITIONED_NORM, KSP_NATURAL_NORM 285 M*/ 286 287 /*M 288 KSP_PRECONDITIONED_NORM - Compute the norm of the preconditioned residual and pass that to the 289 convergence test routine. 290 291 Level: advanced 292 293 .seealso: KSPNormType, KSPSetNormType(), KSP_NO_NORM, KSP_UNPRECONDITIONED_NORM, KSP_NATURAL_NORM, KSPSetConvergenceTest() 294 M*/ 295 296 /*M 297 KSP_UNPRECONDITIONED_NORM - Compute the norm of the true residual (b - A*x) and pass that to the 298 convergence test routine. 299 300 Level: advanced 301 302 .seealso: KSPNormType, KSPSetNormType(), KSP_NO_NORM, KSP_PRECONDITIONED_NORM, KSP_NATURAL_NORM, KSPSetConvergenceTest() 303 M*/ 304 305 /*M 306 KSP_NATURAL_NORM - Compute the 'natural norm' of residual sqrt((b - A*x)*B*(b - A*x)) and pass that to the 307 convergence test routine. 308 309 Level: advanced 310 311 .seealso: KSPNormType, KSPSetNormType(), KSP_NO_NORM, KSP_PRECONDITIONED_NORM, KSP_UNPRECONDITIONED_NORM, KSPSetConvergenceTest() 312 M*/ 313 314 EXTERN PetscErrorCode KSPSetNormType(KSP,KSPNormType); 315 316 /*E 317 KSPConvergedReason - reason a Krylov method was said to 318 have converged or diverged 319 320 Level: beginner 321 322 Notes: this must match finclude/petscksp.h 323 324 Developer note: The string versions of these are in 325 src/ksp/ksp/interface/itfunc.c called convergedreasons. 326 If these enums are changed you much change those. 327 328 .seealso: KSPSolve(), KSPGetConvergedReason() 329 E*/ 330 typedef enum {/* converged */ 331 KSP_CONVERGED_RTOL = 2, 332 KSP_CONVERGED_ATOL = 3, 333 KSP_CONVERGED_ITS = 4, 334 KSP_CONVERGED_QCG_NEG_CURVE = 5, 335 KSP_CONVERGED_QCG_CONSTRAINED = 6, 336 KSP_CONVERGED_STEP_LENGTH = 7, 337 /* diverged */ 338 KSP_DIVERGED_NULL = -2, 339 KSP_DIVERGED_ITS = -3, 340 KSP_DIVERGED_DTOL = -4, 341 KSP_DIVERGED_BREAKDOWN = -5, 342 KSP_DIVERGED_BREAKDOWN_BICG = -6, 343 KSP_DIVERGED_NONSYMMETRIC = -7, 344 KSP_DIVERGED_INDEFINITE_PC = -8, 345 346 KSP_CONVERGED_ITERATING = 0} KSPConvergedReason; 347 348 EXTERN PetscErrorCode KSPSetConvergenceTest(KSP,PetscErrorCode (*)(KSP,int,PetscReal,KSPConvergedReason*,void*),void *); 349 EXTERN PetscErrorCode KSPGetConvergenceContext(KSP,void **); 350 EXTERN PetscErrorCode KSPDefaultConverged(KSP,int,PetscReal,KSPConvergedReason*,void *); 351 EXTERN PetscErrorCode KSPSkipConverged(KSP,int,PetscReal,KSPConvergedReason*,void *); 352 EXTERN PetscErrorCode KSPGetConvergedReason(KSP,KSPConvergedReason *); 353 354 EXTERN PetscErrorCode KSPComputeExplicitOperator(KSP,Mat *); 355 356 /*E 357 KSPCGType - Determines what type of CG to use 358 359 Level: beginner 360 361 .seealso: KSPCGSetType() 362 E*/ 363 typedef enum {KSP_CG_SYMMETRIC=1,KSP_CG_HERMITIAN=2} KSPCGType; 364 365 EXTERN PetscErrorCode KSPCGSetType(KSP,KSPCGType); 366 367 EXTERN PetscErrorCode PCPreSolve(PC,KSP); 368 EXTERN PetscErrorCode PCPostSolve(PC,KSP); 369 370 EXTERN PetscErrorCode KSPLGMonitorCreate(const char[],const char[],int,int,int,int,PetscDrawLG*); 371 EXTERN PetscErrorCode KSPLGMonitor(KSP,int,PetscReal,void*); 372 EXTERN PetscErrorCode KSPLGMonitorDestroy(PetscDrawLG); 373 EXTERN PetscErrorCode KSPLGTrueMonitorCreate(MPI_Comm,const char[],const char[],int,int,int,int,PetscDrawLG*); 374 EXTERN PetscErrorCode KSPLGTrueMonitor(KSP,int,PetscReal,void*); 375 EXTERN PetscErrorCode KSPLGTrueMonitorDestroy(PetscDrawLG); 376 377 PETSC_EXTERN_CXX_END 378 #endif 379