Actual source code: petscksp.h
1: /* $Id: petscksp.h,v 1.107 2001/08/06 21:16:38 bsmith Exp $ */
2: /*
3: Defines the interface functions for the Krylov subspace accelerators.
4: */
5: #ifndef __PETSCKSP_H
7: #include petscpc.h
9: /*S
10: KSP - Abstract PETSc object that manages all Krylov methods
12: Level: beginner
14: Concepts: Krylov methods
16: .seealso: KSPCreate(), KSPSetType(), KSPType, SNES, TS, PC, SLES
17: S*/
18: typedef struct _p_KSP* KSP;
20: /*E
21: KSPType - String with the name of a PETSc Krylov method or the creation function
22: with an optional dynamic library name, for example
23: http://www.mcs.anl.gov/petsc/lib.a:mykspcreate()
25: Level: beginner
27: .seealso: KSPSetType(), KSP
28: E*/
29: #define KSPRICHARDSON "richardson"
30: #define KSPCHEBYCHEV "chebychev"
31: #define KSPCG "cg"
32: #define KSPGMRES "gmres"
33: #define KSPTCQMR "tcqmr"
34: #define KSPBCGS "bcgs"
35: #define KSPCGS "cgs"
36: #define KSPTFQMR "tfqmr"
37: #define KSPCR "cr"
38: #define KSPLSQR "lsqr"
39: #define KSPPREONLY "preonly"
40: #define KSPQCG "qcg"
41: #define KSPBICG "bicg"
42: #define KSPFGMRES "fgmres"
43: #define KSPMINRES "minres"
44: #define KSPSYMMLQ "symmlq"
45: typedef char * KSPType;
47: /* Logging support */
48: extern int KSP_COOKIE;
49: extern int KSP_GMRESOrthogonalization;
51: EXTERN int KSPCreate(MPI_Comm,KSP *);
52: EXTERN int KSPSetType(KSP,KSPType);
53: EXTERN int KSPSetUp(KSP);
54: EXTERN int KSPSolve(KSP,int *);
55: EXTERN int KSPSolveTranspose(KSP,int *);
56: EXTERN int KSPDestroy(KSP);
58: extern PetscFList KSPList;
59: EXTERN int KSPRegisterAll(char *);
60: EXTERN int KSPRegisterDestroy(void);
62: EXTERN int KSPRegister(char*,char*,char*,int(*)(KSP));
63: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
64: #define KSPRegisterDynamic(a,b,c,d) KSPRegister(a,b,c,0)
65: #else
66: #define KSPRegisterDynamic(a,b,c,d) KSPRegister(a,b,c,d)
67: #endif
69: EXTERN int KSPGetType(KSP,KSPType *);
70: EXTERN int KSPSetPreconditionerSide(KSP,PCSide);
71: EXTERN int KSPGetPreconditionerSide(KSP,PCSide*);
72: EXTERN int KSPGetTolerances(KSP,PetscReal*,PetscReal*,PetscReal*,int*);
73: EXTERN int KSPSetTolerances(KSP,PetscReal,PetscReal,PetscReal,int);
74: EXTERN int KSPSetInitialGuessNonzero(KSP,PetscTruth);
75: EXTERN int KSPGetInitialGuessNonzero(KSP,PetscTruth *);
76: EXTERN int KSPSetInitialGuessKnoll(KSP,PetscTruth);
77: EXTERN int KSPGetInitialGuessKnoll(KSP,PetscTruth*);
78: EXTERN int KSPSetComputeEigenvalues(KSP,PetscTruth);
79: EXTERN int KSPSetComputeSingularValues(KSP,PetscTruth);
80: EXTERN int KSPSetRhs(KSP,Vec);
81: EXTERN int KSPGetRhs(KSP,Vec *);
82: EXTERN int KSPSetSolution(KSP,Vec);
83: EXTERN int KSPGetSolution(KSP,Vec *);
84: EXTERN int KSPGetResidualNorm(KSP,PetscReal*);
85: EXTERN int KSPGetIterationNumber(KSP,int*);
87: EXTERN int KSPSetPC(KSP,PC);
88: EXTERN int KSPGetPC(KSP,PC*);
90: EXTERN int KSPSetMonitor(KSP,int (*)(KSP,int,PetscReal,void*),void *,int (*)(void*));
91: EXTERN int KSPClearMonitor(KSP);
92: EXTERN int KSPGetMonitorContext(KSP,void **);
93: EXTERN int KSPGetResidualHistory(KSP,PetscReal **,int *);
94: EXTERN int KSPSetResidualHistory(KSP,PetscReal *,int,PetscTruth);
97: EXTERN int KSPBuildSolution(KSP,Vec,Vec *);
98: EXTERN int KSPBuildResidual(KSP,Vec,Vec,Vec *);
100: EXTERN int KSPRichardsonSetScale(KSP,PetscReal);
101: EXTERN int KSPChebychevSetEigenvalues(KSP,PetscReal,PetscReal);
102: EXTERN int KSPComputeExtremeSingularValues(KSP,PetscReal*,PetscReal*);
103: EXTERN int KSPComputeEigenvalues(KSP,int,PetscReal*,PetscReal*,int *);
104: EXTERN int KSPComputeEigenvaluesExplicitly(KSP,int,PetscReal*,PetscReal*);
106: #define KSPGMRESSetRestart(ksp,r) PetscTryMethod((ksp),KSPGMRESSetRestart_C,(KSP,int),((ksp),(r)))
107: #define KSPGMRESSetHapTol(ksp,tol) PetscTryMethod((ksp),KSPGMRESSetHapTol_C,(KSP,PetscReal),((ksp),(tol)))
109: EXTERN int KSPGMRESSetPreAllocateVectors(KSP);
110: EXTERN int KSPGMRESSetOrthogonalization(KSP,int (*)(KSP,int));
111: EXTERN int KSPGMRESUnmodifiedGramSchmidtOrthogonalization(KSP,int);
112: EXTERN int KSPGMRESModifiedGramSchmidtOrthogonalization(KSP,int);
113: EXTERN int KSPGMRESIROrthogonalization(KSP,int);
115: EXTERN int KSPFGMRESModifyPCNoChange(KSP,int,int,PetscReal,void*);
116: EXTERN int KSPFGMRESModifyPCSLES(KSP,int,int,PetscReal,void*);
117: EXTERN int KSPFGMRESSetModifyPC(KSP,int (*)(KSP,int,int,PetscReal,void*),void*,int(*)(void*));
119: EXTERN int KSPQCGSetTrustRegionRadius(KSP,PetscReal);
120: EXTERN int KSPQCGGetQuadratic(KSP,PetscReal*);
121: EXTERN int KSPQCGGetTrialStepNorm(KSP,PetscReal*);
123: EXTERN int KSPSetFromOptions(KSP);
124: EXTERN int KSPAddOptionsChecker(int (*)(KSP));
126: EXTERN int KSPSingularValueMonitor(KSP,int,PetscReal,void *);
127: EXTERN int KSPDefaultMonitor(KSP,int,PetscReal,void *);
128: EXTERN int KSPTrueMonitor(KSP,int,PetscReal,void *);
129: EXTERN int KSPDefaultSMonitor(KSP,int,PetscReal,void *);
130: EXTERN int KSPVecViewMonitor(KSP,int,PetscReal,void *);
131: EXTERN int KSPGMRESKrylovMonitor(KSP,int,PetscReal,void *);
133: EXTERN int KSPUnwindPreconditioner(KSP,Vec,Vec);
134: EXTERN int KSPDefaultBuildSolution(KSP,Vec,Vec*);
135: EXTERN int KSPDefaultBuildResidual(KSP,Vec,Vec,Vec *);
137: EXTERN int KSPSetOptionsPrefix(KSP,char*);
138: EXTERN int KSPAppendOptionsPrefix(KSP,char*);
139: EXTERN int KSPGetOptionsPrefix(KSP,char**);
141: EXTERN int KSPView(KSP,PetscViewer);
143: /*E
144: KSPNormType - Norm that is passed in the Krylov convergence
145: test routines.
147: Level: advanced
149: Notes: this must match finclude/petscksp.h
151: .seealso: SLESSolve(), KSPSolve(), KSPGetConvergedReason(), KSPSetNormType(),
152: KSPSetConvergenceTest()
153: E*/
154: typedef enum {KSP_NO_NORM = 0,
155: KSP_PRECONDITIONED_NORM = 1,
156: KSP_UNPRECONDITIONED_NORM = 2,
157: KSP_NATURAL_NORM = 3} KSPNormType;
158: EXTERN int KSPSetNormType(KSP,KSPNormType);
160: /*E
161: KSPConvergedReason - reason a Krylov method was said to
162: have converged or diverged
164: Level: beginner
166: Notes: this must match finclude/petscksp.h
168: .seealso: SLESSolve(), KSPSolve(), KSPGetConvergedReason()
169: E*/
170: typedef enum {/* converged */
171: KSP_CONVERGED_RTOL = 2,
172: KSP_CONVERGED_ATOL = 3,
173: KSP_CONVERGED_ITS = 4,
174: KSP_CONVERGED_QCG_NEG_CURVE = 5,
175: KSP_CONVERGED_QCG_CONSTRAINED = 6,
176: KSP_CONVERGED_STEP_LENGTH = 7,
177: /* diverged */
178: KSP_DIVERGED_ITS = -3,
179: KSP_DIVERGED_DTOL = -4,
180: KSP_DIVERGED_BREAKDOWN = -5,
181: KSP_DIVERGED_BREAKDOWN_BICG = -6,
182: KSP_DIVERGED_NONSYMMETRIC = -7,
183: KSP_DIVERGED_INDEFINITE_PC = -8,
184:
185: KSP_CONVERGED_ITERATING = 0} KSPConvergedReason;
187: EXTERN int KSPSetConvergenceTest(KSP,int (*)(KSP,int,PetscReal,KSPConvergedReason*,void*),void *);
188: EXTERN int KSPGetConvergenceContext(KSP,void **);
189: EXTERN int KSPDefaultConverged(KSP,int,PetscReal,KSPConvergedReason*,void *);
190: EXTERN int KSPSkipConverged(KSP,int,PetscReal,KSPConvergedReason*,void *);
191: EXTERN int KSPGetConvergedReason(KSP,KSPConvergedReason *);
193: EXTERN int KSPComputeExplicitOperator(KSP,Mat *);
195: /*E
196: KSPCGType - Determines what type of CG to use
198: Level: beginner
200: .seealso: KSPCGSetType()
201: E*/
202: typedef enum {KSP_CG_SYMMETRIC=1,KSP_CG_HERMITIAN=2} KSPCGType;
204: EXTERN int KSPCGSetType(KSP,KSPCGType);
206: EXTERN int PCPreSolve(PC,KSP);
207: EXTERN int PCPostSolve(PC,KSP);
209: EXTERN int KSPLGMonitorCreate(char*,char*,int,int,int,int,PetscDrawLG*);
210: EXTERN int KSPLGMonitor(KSP,int,PetscReal,void*);
211: EXTERN int KSPLGMonitorDestroy(PetscDrawLG);
212: EXTERN int KSPLGTrueMonitorCreate(MPI_Comm,char*,char*,int,int,int,int,PetscDrawLG*);
213: EXTERN int KSPLGTrueMonitor(KSP,int,PetscReal,void*);
214: EXTERN int KSPLGTrueMonitorDestroy(PetscDrawLG);
216: #endif