Actual source code: petscksp.h
1: /*
2: Defines the interface functions for the Krylov subspace accelerators.
3: */
4: #ifndef __PETSCKSP_H
6: #include petscpc.h
11: /*S
12: KSP - Abstract PETSc object that manages all Krylov methods
14: Level: beginner
16: Concepts: Krylov methods
18: .seealso: KSPCreate(), KSPSetType(), KSPType, SNES, TS, PC, KSP
19: S*/
20: typedef struct _p_KSP* KSP;
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()
27: Level: beginner
29: .seealso: KSPSetType(), KSP
30: E*/
31: #define KSPType char*
32: #define KSPRICHARDSON "richardson"
33: #define KSPCHEBYCHEV "chebychev"
34: #define KSPCG "cg"
35: #define KSPCGNE "cgne"
36: #define KSPNASH "nash"
37: #define KSPSTCG "stcg"
38: #define KSPGLTR "gltr"
39: #define KSPGMRES "gmres"
40: #define KSPFGMRES "fgmres"
41: #define KSPLGMRES "lgmres"
42: #define KSPDGMRES "dgmres"
43: #define KSPTCQMR "tcqmr"
44: #define KSPBCGS "bcgs"
45: #define KSPIBCGS "ibcgs"
46: #define KSPBCGSL "bcgsl"
47: #define KSPCGS "cgs"
48: #define KSPTFQMR "tfqmr"
49: #define KSPCR "cr"
50: #define KSPLSQR "lsqr"
51: #define KSPPREONLY "preonly"
52: #define KSPQCG "qcg"
53: #define KSPBICG "bicg"
54: #define KSPMINRES "minres"
55: #define KSPSYMMLQ "symmlq"
56: #define KSPLCD "lcd"
57: #define KSPPYTHON "python"
58: #define KSPBROYDEN "broyden"
59: #define KSPGCR "gcr"
60: #define KSPNGMRES "ngmres"
61: #define KSPSPECEST "specest"
63: /* Logging support */
81: /*MC
82: KSPRegisterDynamic - Adds a method to the Krylov subspace solver package.
84: Synopsis:
85: PetscErrorCode KSPRegisterDynamic(const char *name_solver,const char *path,const char *name_create,PetscErrorCode (*routine_create)(KSP))
87: Not Collective
89: Input Parameters:
90: + name_solver - name of a new user-defined solver
91: . path - path (either absolute or relative) the library containing this solver
92: . name_create - name of routine to create method context
93: - routine_create - routine to create method context
95: Notes:
96: KSPRegisterDynamic() may be called multiple times to add several user-defined solvers.
98: If dynamic libraries are used, then the fourth input argument (routine_create)
99: is ignored.
101: Sample usage:
102: .vb
103: KSPRegisterDynamic("my_solver",/home/username/my_lib/lib/libO/solaris/mylib.a,
104: "MySolverCreate",MySolverCreate);
105: .ve
107: Then, your solver can be chosen with the procedural interface via
108: $ KSPSetType(ksp,"my_solver")
109: or at runtime via the option
110: $ -ksp_type my_solver
112: Level: advanced
114: Notes: Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR},
115: and others of the form ${any_environmental_variable} occuring in pathname will be
116: replaced with appropriate values.
117: If your function is not being put into a shared library then use KSPRegister() instead
119: .keywords: KSP, register
121: .seealso: KSPRegisterAll(), KSPRegisterDestroy()
123: M*/
124: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
125: #define KSPRegisterDynamic(a,b,c,d) KSPRegister(a,b,c,0)
126: #else
127: #define KSPRegisterDynamic(a,b,c,d) KSPRegister(a,b,c,d)
128: #endif
162: /* not sure where to put this */
199: /*E
200: KSPGMRESCGSRefinementType - How the classical (unmodified) Gram-Schmidt is performed.
202: Level: advanced
204: .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), KSPGMRESGetOrthogonalization(),
205: KSPGMRESSetCGSRefinementType(), KSPGMRESGetCGSRefinementType(), KSPGMRESModifiedGramSchmidtOrthogonalization()
207: E*/
208: typedef enum {KSP_GMRES_CGS_REFINE_NEVER, KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS} KSPGMRESCGSRefinementType;
210: /*MC
211: KSP_GMRES_CGS_REFINE_NEVER - Just do the classical (unmodified) Gram-Schmidt process
213: Level: advanced
215: Note: Possible unstable, but the fastest to compute
217: .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), KSPGMRESGetOrthogonalization(),
218: KSPGMRESSetCGSRefinementType(), KSPGMRESGetCGSRefinementType(), KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS,
219: KSPGMRESModifiedGramSchmidtOrthogonalization()
220: M*/
222: /*MC
223: KSP_GMRES_CGS_REFINE_IFNEEDED - Do the classical (unmodified) Gram-Schmidt process and one step of
224: iterative refinement if an estimate of the orthogonality of the resulting vectors indicates
225: poor orthogonality.
227: Level: advanced
229: Note: This is slower than KSP_GMRES_CGS_REFINE_NEVER because it requires an extra norm computation to
230: estimate the orthogonality but is more stable.
232: .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), KSPGMRESGetOrthogonalization(),
233: KSPGMRESSetCGSRefinementType(), KSPGMRESGetCGSRefinementType(), KSP_GMRES_CGS_REFINE_NEVER, KSP_GMRES_CGS_REFINE_ALWAYS,
234: KSPGMRESModifiedGramSchmidtOrthogonalization()
235: M*/
237: /*MC
238: KSP_GMRES_CGS_REFINE_NEVER - Do two steps of the classical (unmodified) Gram-Schmidt process.
240: Level: advanced
242: Note: This is roughly twice the cost of KSP_GMRES_CGS_REFINE_NEVER because it performs the process twice
243: but it saves the extra norm calculation needed by KSP_GMRES_CGS_REFINE_IFNEEDED.
245: You should only use this if you absolutely know that the iterative refinement is needed.
247: .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), KSPGMRESGetOrthogonalization(),
248: KSPGMRESSetCGSRefinementType(), KSPGMRESGetCGSRefinementType(), KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS,
249: KSPGMRESModifiedGramSchmidtOrthogonalization()
250: M*/
304: /*E
305: KSPNormType - Norm that is passed in the Krylov convergence
306: test routines.
308: Level: advanced
310: Each solver only supports a subset of these and some may support different ones
311: depending on left or right preconditioning, see KSPSetPCSide()
313: Notes: this must match finclude/petscksp.h
315: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPSetNormType(),
316: KSPSetConvergenceTest(), KSPSetPCSide()
317: E*/
318: typedef enum {KSP_NORM_DEFAULT = -1,KSP_NORM_NONE = 0,KSP_NORM_PRECONDITIONED = 1,KSP_NORM_UNPRECONDITIONED = 2,KSP_NORM_NATURAL = 3} KSPNormType;
319: #define KSP_NORM_MAX (KSP_NORM_NATURAL + 1)
322: /*MC
323: KSP_NORM_NONE - Do not compute a norm during the Krylov process. This will
324: possibly save some computation but means the convergence test cannot
325: be based on a norm of a residual etc.
327: Level: advanced
329: Note: Some Krylov methods need to compute a residual norm and then this option is ignored
331: .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_PRECONDITIONED, KSP_NORM_UNPRECONDITIONED, KSP_NORM_NATURAL
332: M*/
334: /*MC
335: KSP_NORM_PRECONDITIONED - Compute the norm of the preconditioned residual and pass that to the
336: convergence test routine.
338: Level: advanced
340: .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_NONE, KSP_NORM_UNPRECONDITIONED, KSP_NORM_NATURAL, KSPSetConvergenceTest()
341: M*/
343: /*MC
344: KSP_NORM_UNPRECONDITIONED - Compute the norm of the true residual (b - A*x) and pass that to the
345: convergence test routine.
347: Level: advanced
349: .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_NONE, KSP_NORM_PRECONDITIONED, KSP_NORM_NATURAL, KSPSetConvergenceTest()
350: M*/
352: /*MC
353: KSP_NORM_NATURAL - Compute the 'natural norm' of residual sqrt((b - A*x)*B*(b - A*x)) and pass that to the
354: convergence test routine. This is only supported by KSPCG, KSPCR, KSPCGNE, KSPCGS
356: Level: advanced
358: .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_NONE, KSP_NORM_PRECONDITIONED, KSP_NORM_UNPRECONDITIONED, KSPSetConvergenceTest()
359: M*/
367: /*E
368: KSPConvergedReason - reason a Krylov method was said to
369: have converged or diverged
371: Level: beginner
373: Notes: See KSPGetConvergedReason() for explanation of each value
375: Developer notes: this must match finclude/petscksp.h
377: The string versions of these are KSPConvergedReasons; if you change
378: any of the values here also change them that array of names.
380: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPSetTolerances()
381: E*/
382: typedef enum {/* converged */
383: KSP_CONVERGED_RTOL_NORMAL = 1,
384: KSP_CONVERGED_ATOL_NORMAL = 9,
385: KSP_CONVERGED_RTOL = 2,
386: KSP_CONVERGED_ATOL = 3,
387: KSP_CONVERGED_ITS = 4,
388: KSP_CONVERGED_CG_NEG_CURVE = 5,
389: KSP_CONVERGED_CG_CONSTRAINED = 6,
390: KSP_CONVERGED_STEP_LENGTH = 7,
391: KSP_CONVERGED_HAPPY_BREAKDOWN = 8,
392: /* diverged */
393: KSP_DIVERGED_NULL = -2,
394: KSP_DIVERGED_ITS = -3,
395: KSP_DIVERGED_DTOL = -4,
396: KSP_DIVERGED_BREAKDOWN = -5,
397: KSP_DIVERGED_BREAKDOWN_BICG = -6,
398: KSP_DIVERGED_NONSYMMETRIC = -7,
399: KSP_DIVERGED_INDEFINITE_PC = -8,
400: KSP_DIVERGED_NAN = -9,
401: KSP_DIVERGED_INDEFINITE_MAT = -10,
402:
403: KSP_CONVERGED_ITERATING = 0} KSPConvergedReason;
406: /*MC
407: KSP_CONVERGED_RTOL - norm(r) <= rtol*norm(b)
409: Level: beginner
411: See KSPNormType and KSPSetNormType() for possible norms that may be used. By default
412: for left preconditioning it is the 2-norm of the preconditioned residual, and the
413: 2-norm of the residual for right preconditioning
415: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
417: M*/
419: /*MC
420: KSP_CONVERGED_ATOL - norm(r) <= atol
422: Level: beginner
424: See KSPNormType and KSPSetNormType() for possible norms that may be used. By default
425: for left preconditioning it is the 2-norm of the preconditioned residual, and the
426: 2-norm of the residual for right preconditioning
428: Level: beginner
430: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
432: M*/
434: /*MC
435: KSP_DIVERGED_DTOL - norm(r) >= dtol*norm(b)
437: Level: beginner
439: See KSPNormType and KSPSetNormType() for possible norms that may be used. By default
440: for left preconditioning it is the 2-norm of the preconditioned residual, and the
441: 2-norm of the residual for right preconditioning
443: Level: beginner
445: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
447: M*/
449: /*MC
450: KSP_DIVERGED_ITS - Ran out of iterations before any convergence criteria was
451: reached
453: Level: beginner
455: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
457: M*/
459: /*MC
460: KSP_CONVERGED_ITS - Used by the KSPPREONLY solver after the single iteration of
461: the preconditioner is applied. Also used when the KSPSkipConverged() convergence
462: test routine is set in KSP.
465: Level: beginner
468: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
470: M*/
472: /*MC
473: KSP_DIVERGED_BREAKDOWN - A breakdown in the Krylov method was detected so the
474: method could not continue to enlarge the Krylov space. Could be due to a singlular matrix or
475: preconditioner.
477: Level: beginner
479: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
481: M*/
483: /*MC
484: KSP_DIVERGED_BREAKDOWN_BICG - A breakdown in the KSPBICG method was detected so the
485: method could not continue to enlarge the Krylov space.
488: Level: beginner
491: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
493: M*/
495: /*MC
496: KSP_DIVERGED_NONSYMMETRIC - It appears the operator or preconditioner is not
497: symmetric and this Krylov method (KSPCG, KSPMINRES, KSPCR) requires symmetry
499: Level: beginner
501: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
503: M*/
505: /*MC
506: KSP_DIVERGED_INDEFINITE_PC - It appears the preconditioner is indefinite (has both
507: positive and negative eigenvalues) and this Krylov method (KSPCG) requires it to
508: be positive definite
510: Level: beginner
512: Notes: This can happen with the PCICC preconditioner, use -pc_factor_shift_positive_definite to force
513: the PCICC preconditioner to generate a positive definite preconditioner
515: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
517: M*/
519: /*MC
520: KSP_CONVERGED_ITERATING - This flag is returned if you call KSPGetConvergedReason()
521: while the KSPSolve() is still running.
523: Level: beginner
525: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
527: M*/
542: /*E
543: KSPCGType - Determines what type of CG to use
545: Level: beginner
547: .seealso: KSPCGSetType()
548: E*/
549: typedef enum {KSP_CG_SYMMETRIC=0,KSP_CG_HERMITIAN=1} KSPCGType;
587: /* see src/ksp/ksp/interface/iguess.c */
588: typedef struct _p_KSPFischerGuess {PetscInt method,curl,maxl,refcnt;PetscBool monitor;Mat mat; KSP ksp;}* KSPFischerGuess;
616: #endif