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
7: PETSC_EXTERN_CXX_BEGIN
9: extern PetscErrorCode KSPInitializePackage(const char[]);
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: /*J
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: J*/
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 */
64: extern PetscClassId KSP_CLASSID;
66: extern PetscErrorCode KSPCreate(MPI_Comm,KSP *);
67: extern PetscErrorCode KSPSetType(KSP,const KSPType);
68: extern PetscErrorCode KSPSetUp(KSP);
69: extern PetscErrorCode KSPSetUpOnBlocks(KSP);
70: extern PetscErrorCode KSPSolve(KSP,Vec,Vec);
71: extern PetscErrorCode KSPSolveTranspose(KSP,Vec,Vec);
72: extern PetscErrorCode KSPReset(KSP);
73: extern PetscErrorCode KSPDestroy(KSP*);
75: extern PetscFList KSPList;
76: extern PetscBool KSPRegisterAllCalled;
77: extern PetscErrorCode KSPRegisterAll(const char[]);
78: extern PetscErrorCode KSPRegisterDestroy(void);
79: extern PetscErrorCode KSPRegister(const char[],const char[],const char[],PetscErrorCode (*)(KSP));
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
130: extern PetscErrorCode KSPGetType(KSP,const KSPType *);
131: extern PetscErrorCode KSPSetPCSide(KSP,PCSide);
132: extern PetscErrorCode KSPGetPCSide(KSP,PCSide*);
133: extern PetscErrorCode KSPGetTolerances(KSP,PetscReal*,PetscReal*,PetscReal*,PetscInt*);
134: extern PetscErrorCode KSPSetTolerances(KSP,PetscReal,PetscReal,PetscReal,PetscInt);
135: extern PetscErrorCode KSPSetInitialGuessNonzero(KSP,PetscBool );
136: extern PetscErrorCode KSPGetInitialGuessNonzero(KSP,PetscBool *);
137: extern PetscErrorCode KSPSetInitialGuessKnoll(KSP,PetscBool );
138: extern PetscErrorCode KSPGetInitialGuessKnoll(KSP,PetscBool *);
139: extern PetscErrorCode KSPSetErrorIfNotConverged(KSP,PetscBool );
140: extern PetscErrorCode KSPGetErrorIfNotConverged(KSP,PetscBool *);
141: extern PetscErrorCode KSPGetComputeEigenvalues(KSP,PetscBool *);
142: extern PetscErrorCode KSPSetComputeEigenvalues(KSP,PetscBool );
143: extern PetscErrorCode KSPGetComputeSingularValues(KSP,PetscBool *);
144: extern PetscErrorCode KSPSetComputeSingularValues(KSP,PetscBool );
145: extern PetscErrorCode KSPGetRhs(KSP,Vec *);
146: extern PetscErrorCode KSPGetSolution(KSP,Vec *);
147: extern PetscErrorCode KSPGetResidualNorm(KSP,PetscReal*);
148: extern PetscErrorCode KSPGetIterationNumber(KSP,PetscInt*);
149: extern PetscErrorCode KSPSetNullSpace(KSP,MatNullSpace);
150: extern PetscErrorCode KSPGetNullSpace(KSP,MatNullSpace*);
151: extern PetscErrorCode KSPGetVecs(KSP,PetscInt,Vec**,PetscInt,Vec**);
153: extern PetscErrorCode KSPSetPC(KSP,PC);
154: extern PetscErrorCode KSPGetPC(KSP,PC*);
156: extern PetscErrorCode KSPMonitor(KSP,PetscInt,PetscReal);
157: extern PetscErrorCode KSPMonitorSet(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscReal,void*),void *,PetscErrorCode (*)(void**));
158: extern PetscErrorCode KSPMonitorCancel(KSP);
159: extern PetscErrorCode KSPGetMonitorContext(KSP,void **);
160: extern PetscErrorCode KSPGetResidualHistory(KSP,PetscReal*[],PetscInt *);
161: extern PetscErrorCode KSPSetResidualHistory(KSP,PetscReal[],PetscInt,PetscBool );
163: /* not sure where to put this */
164: extern PetscErrorCode PCKSPGetKSP(PC,KSP*);
165: extern PetscErrorCode PCBJacobiGetSubKSP(PC,PetscInt*,PetscInt*,KSP*[]);
166: extern PetscErrorCode PCASMGetSubKSP(PC,PetscInt*,PetscInt*,KSP*[]);
167: extern PetscErrorCode PCGASMGetSubKSP(PC,PetscInt*,PetscInt*,KSP*[]);
168: extern PetscErrorCode PCFieldSplitGetSubKSP(PC,PetscInt*,KSP*[]);
170: extern PetscErrorCode PCGalerkinGetKSP(PC,KSP *);
172: extern PetscErrorCode KSPBuildSolution(KSP,Vec,Vec *);
173: extern PetscErrorCode KSPBuildResidual(KSP,Vec,Vec,Vec *);
175: extern PetscErrorCode KSPRichardsonSetScale(KSP,PetscReal);
176: extern PetscErrorCode KSPRichardsonSetSelfScale(KSP,PetscBool );
177: extern PetscErrorCode KSPChebychevSetEigenvalues(KSP,PetscReal,PetscReal);
178: extern PetscErrorCode KSPChebychevSetEstimateEigenvalues(KSP,PetscReal,PetscReal,PetscReal,PetscReal);
179: extern PetscErrorCode KSPComputeExtremeSingularValues(KSP,PetscReal*,PetscReal*);
180: extern PetscErrorCode KSPComputeEigenvalues(KSP,PetscInt,PetscReal*,PetscReal*,PetscInt *);
181: extern PetscErrorCode KSPComputeEigenvaluesExplicitly(KSP,PetscInt,PetscReal*,PetscReal*);
183: extern PetscErrorCode KSPGMRESSetRestart(KSP, PetscInt);
184: extern PetscErrorCode KSPGMRESGetRestart(KSP, PetscInt*);
185: extern PetscErrorCode KSPGMRESSetHapTol(KSP,PetscReal);
187: extern PetscErrorCode KSPGMRESSetPreAllocateVectors(KSP);
188: extern PetscErrorCode KSPGMRESSetOrthogonalization(KSP,PetscErrorCode (*)(KSP,PetscInt));
189: extern PetscErrorCode KSPGMRESGetOrthogonalization(KSP,PetscErrorCode (**)(KSP,PetscInt));
190: extern PetscErrorCode KSPGMRESModifiedGramSchmidtOrthogonalization(KSP,PetscInt);
191: extern PetscErrorCode KSPGMRESClassicalGramSchmidtOrthogonalization(KSP,PetscInt);
193: extern PetscErrorCode KSPLGMRESSetAugDim(KSP,PetscInt);
194: extern PetscErrorCode KSPLGMRESSetConstant(KSP);
196: extern PetscErrorCode KSPGCRSetRestart(KSP,PetscInt);
197: extern PetscErrorCode KSPGCRGetRestart(KSP,PetscInt*);
198: extern PetscErrorCode KSPGCRSetModifyPC(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscReal,void*),void*,PetscErrorCode(*)(void*));
200: /*E
201: KSPGMRESCGSRefinementType - How the classical (unmodified) Gram-Schmidt is performed.
203: Level: advanced
205: .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), KSPGMRESGetOrthogonalization(),
206: KSPGMRESSetCGSRefinementType(), KSPGMRESGetCGSRefinementType(), KSPGMRESModifiedGramSchmidtOrthogonalization()
208: E*/
209: typedef enum {KSP_GMRES_CGS_REFINE_NEVER, KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS} KSPGMRESCGSRefinementType;
210: extern const char *KSPGMRESCGSRefinementTypes[];
211: /*MC
212: KSP_GMRES_CGS_REFINE_NEVER - Just do the classical (unmodified) Gram-Schmidt process
214: Level: advanced
216: Note: Possible unstable, but the fastest to compute
218: .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), KSPGMRESGetOrthogonalization(),
219: KSPGMRESSetCGSRefinementType(), KSPGMRESGetCGSRefinementType(), KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS,
220: KSPGMRESModifiedGramSchmidtOrthogonalization()
221: M*/
223: /*MC
224: KSP_GMRES_CGS_REFINE_IFNEEDED - Do the classical (unmodified) Gram-Schmidt process and one step of
225: iterative refinement if an estimate of the orthogonality of the resulting vectors indicates
226: poor orthogonality.
228: Level: advanced
230: Note: This is slower than KSP_GMRES_CGS_REFINE_NEVER because it requires an extra norm computation to
231: estimate the orthogonality but is more stable.
233: .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), KSPGMRESGetOrthogonalization(),
234: KSPGMRESSetCGSRefinementType(), KSPGMRESGetCGSRefinementType(), KSP_GMRES_CGS_REFINE_NEVER, KSP_GMRES_CGS_REFINE_ALWAYS,
235: KSPGMRESModifiedGramSchmidtOrthogonalization()
236: M*/
238: /*MC
239: KSP_GMRES_CGS_REFINE_NEVER - Do two steps of the classical (unmodified) Gram-Schmidt process.
241: Level: advanced
243: Note: This is roughly twice the cost of KSP_GMRES_CGS_REFINE_NEVER because it performs the process twice
244: but it saves the extra norm calculation needed by KSP_GMRES_CGS_REFINE_IFNEEDED.
246: You should only use this if you absolutely know that the iterative refinement is needed.
248: .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(), KSPGMRESGetOrthogonalization(),
249: KSPGMRESSetCGSRefinementType(), KSPGMRESGetCGSRefinementType(), KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS,
250: KSPGMRESModifiedGramSchmidtOrthogonalization()
251: M*/
253: extern PetscErrorCode KSPGMRESSetCGSRefinementType(KSP,KSPGMRESCGSRefinementType);
254: extern PetscErrorCode KSPGMRESGetCGSRefinementType(KSP,KSPGMRESCGSRefinementType*);
256: extern PetscErrorCode KSPFGMRESModifyPCNoChange(KSP,PetscInt,PetscInt,PetscReal,void*);
257: extern PetscErrorCode KSPFGMRESModifyPCKSP(KSP,PetscInt,PetscInt,PetscReal,void*);
258: extern PetscErrorCode KSPFGMRESSetModifyPC(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscInt,PetscReal,void*),void*,PetscErrorCode(*)(void*));
260: extern PetscErrorCode KSPQCGSetTrustRegionRadius(KSP,PetscReal);
261: extern PetscErrorCode KSPQCGGetQuadratic(KSP,PetscReal*);
262: extern PetscErrorCode KSPQCGGetTrialStepNorm(KSP,PetscReal*);
264: extern PetscErrorCode KSPBCGSLSetXRes(KSP,PetscReal);
265: extern PetscErrorCode KSPBCGSLSetPol(KSP,PetscBool );
266: extern PetscErrorCode KSPBCGSLSetEll(KSP,PetscInt);
268: extern PetscErrorCode KSPSetFromOptions(KSP);
269: extern PetscErrorCode KSPAddOptionsChecker(PetscErrorCode (*)(KSP));
271: extern PetscErrorCode KSPMonitorSingularValue(KSP,PetscInt,PetscReal,void *);
272: extern PetscErrorCode KSPMonitorDefault(KSP,PetscInt,PetscReal,void *);
273: extern PetscErrorCode KSPMonitorDefaultLSQR(KSP,PetscInt,PetscReal,void *);
274: extern PetscErrorCode KSPMonitorRange(KSP,PetscInt,PetscReal,void *);
275: extern PetscErrorCode KSPMonitorTrueResidualNorm(KSP,PetscInt,PetscReal,void *);
276: extern PetscErrorCode KSPMonitorDefaultShort(KSP,PetscInt,PetscReal,void *);
277: extern PetscErrorCode KSPMonitorSolution(KSP,PetscInt,PetscReal,void *);
278: extern PetscErrorCode KSPGMRESMonitorKrylov(KSP,PetscInt,PetscReal,void *);
280: extern PetscErrorCode KSPUnwindPreconditioner(KSP,Vec,Vec);
281: extern PetscErrorCode KSPDefaultBuildSolution(KSP,Vec,Vec*);
282: extern PetscErrorCode KSPDefaultBuildResidual(KSP,Vec,Vec,Vec *);
283: extern PetscErrorCode KSPInitialResidual(KSP,Vec,Vec,Vec,Vec,Vec);
285: extern PetscErrorCode KSPSetOperators(KSP,Mat,Mat,MatStructure);
286: extern PetscErrorCode KSPGetOperators(KSP,Mat*,Mat*,MatStructure*);
287: extern PetscErrorCode KSPGetOperatorsSet(KSP,PetscBool *,PetscBool *);
288: extern PetscErrorCode KSPSetOptionsPrefix(KSP,const char[]);
289: extern PetscErrorCode KSPAppendOptionsPrefix(KSP,const char[]);
290: extern PetscErrorCode KSPGetOptionsPrefix(KSP,const char*[]);
292: extern PetscErrorCode KSPSetDiagonalScale(KSP,PetscBool );
293: extern PetscErrorCode KSPGetDiagonalScale(KSP,PetscBool *);
294: extern PetscErrorCode KSPSetDiagonalScaleFix(KSP,PetscBool );
295: extern PetscErrorCode KSPGetDiagonalScaleFix(KSP,PetscBool *);
297: extern PetscErrorCode KSPView(KSP,PetscViewer);
299: extern PetscErrorCode KSPLSQRSetStandardErrorVec(KSP,Vec);
300: extern PetscErrorCode KSPLSQRGetStandardErrorVec(KSP,Vec*);
302: extern PetscErrorCode PCRedundantGetKSP(PC,KSP*);
303: extern PetscErrorCode PCRedistributeGetKSP(PC,KSP*);
305: /*E
306: KSPNormType - Norm that is passed in the Krylov convergence
307: test routines.
309: Level: advanced
311: Each solver only supports a subset of these and some may support different ones
312: depending on left or right preconditioning, see KSPSetPCSide()
314: Notes: this must match finclude/petscksp.h
316: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPSetNormType(),
317: KSPSetConvergenceTest(), KSPSetPCSide()
318: E*/
319: typedef enum {KSP_NORM_DEFAULT = -1,KSP_NORM_NONE = 0,KSP_NORM_PRECONDITIONED = 1,KSP_NORM_UNPRECONDITIONED = 2,KSP_NORM_NATURAL = 3} KSPNormType;
320: #define KSP_NORM_MAX (KSP_NORM_NATURAL + 1)
321: extern const char *const*const KSPNormTypes;
323: /*MC
324: KSP_NORM_NONE - Do not compute a norm during the Krylov process. This will
325: possibly save some computation but means the convergence test cannot
326: be based on a norm of a residual etc.
328: Level: advanced
330: Note: Some Krylov methods need to compute a residual norm and then this option is ignored
332: .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_PRECONDITIONED, KSP_NORM_UNPRECONDITIONED, KSP_NORM_NATURAL
333: M*/
335: /*MC
336: KSP_NORM_PRECONDITIONED - Compute the norm of the preconditioned residual and pass that to the
337: convergence test routine.
339: Level: advanced
341: .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_NONE, KSP_NORM_UNPRECONDITIONED, KSP_NORM_NATURAL, KSPSetConvergenceTest()
342: M*/
344: /*MC
345: KSP_NORM_UNPRECONDITIONED - Compute the norm of the true residual (b - A*x) and pass that to the
346: convergence test routine.
348: Level: advanced
350: .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_NONE, KSP_NORM_PRECONDITIONED, KSP_NORM_NATURAL, KSPSetConvergenceTest()
351: M*/
353: /*MC
354: KSP_NORM_NATURAL - Compute the 'natural norm' of residual sqrt((b - A*x)*B*(b - A*x)) and pass that to the
355: convergence test routine. This is only supported by KSPCG, KSPCR, KSPCGNE, KSPCGS
357: Level: advanced
359: .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_NONE, KSP_NORM_PRECONDITIONED, KSP_NORM_UNPRECONDITIONED, KSPSetConvergenceTest()
360: M*/
362: extern PetscErrorCode KSPSetNormType(KSP,KSPNormType);
363: extern PetscErrorCode KSPGetNormType(KSP,KSPNormType*);
364: extern PetscErrorCode KSPSetSupportedNorm(KSP ksp,KSPNormType,PCSide,PetscInt);
365: extern PetscErrorCode KSPSetCheckNormIteration(KSP,PetscInt);
366: extern PetscErrorCode KSPSetLagNorm(KSP,PetscBool );
368: /*E
369: KSPConvergedReason - reason a Krylov method was said to
370: have converged or diverged
372: Level: beginner
374: Notes: See KSPGetConvergedReason() for explanation of each value
376: Developer notes: this must match finclude/petscksp.h
378: The string versions of these are KSPConvergedReasons; if you change
379: any of the values here also change them that array of names.
381: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPSetTolerances()
382: E*/
383: typedef enum {/* converged */
384: KSP_CONVERGED_RTOL_NORMAL = 1,
385: KSP_CONVERGED_ATOL_NORMAL = 9,
386: KSP_CONVERGED_RTOL = 2,
387: KSP_CONVERGED_ATOL = 3,
388: KSP_CONVERGED_ITS = 4,
389: KSP_CONVERGED_CG_NEG_CURVE = 5,
390: KSP_CONVERGED_CG_CONSTRAINED = 6,
391: KSP_CONVERGED_STEP_LENGTH = 7,
392: KSP_CONVERGED_HAPPY_BREAKDOWN = 8,
393: /* diverged */
394: KSP_DIVERGED_NULL = -2,
395: KSP_DIVERGED_ITS = -3,
396: KSP_DIVERGED_DTOL = -4,
397: KSP_DIVERGED_BREAKDOWN = -5,
398: KSP_DIVERGED_BREAKDOWN_BICG = -6,
399: KSP_DIVERGED_NONSYMMETRIC = -7,
400: KSP_DIVERGED_INDEFINITE_PC = -8,
401: KSP_DIVERGED_NAN = -9,
402: KSP_DIVERGED_INDEFINITE_MAT = -10,
403:
404: KSP_CONVERGED_ITERATING = 0} KSPConvergedReason;
405: extern const char *const*KSPConvergedReasons;
407: /*MC
408: KSP_CONVERGED_RTOL - norm(r) <= rtol*norm(b)
410: Level: beginner
412: See KSPNormType and KSPSetNormType() for possible norms that may be used. By default
413: for left preconditioning it is the 2-norm of the preconditioned residual, and the
414: 2-norm of the residual for right preconditioning
416: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
418: M*/
420: /*MC
421: KSP_CONVERGED_ATOL - norm(r) <= atol
423: Level: beginner
425: See KSPNormType and KSPSetNormType() for possible norms that may be used. By default
426: for left preconditioning it is the 2-norm of the preconditioned residual, and the
427: 2-norm of the residual for right preconditioning
429: Level: beginner
431: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
433: M*/
435: /*MC
436: KSP_DIVERGED_DTOL - norm(r) >= dtol*norm(b)
438: Level: beginner
440: See KSPNormType and KSPSetNormType() for possible norms that may be used. By default
441: for left preconditioning it is the 2-norm of the preconditioned residual, and the
442: 2-norm of the residual for right preconditioning
444: Level: beginner
446: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
448: M*/
450: /*MC
451: KSP_DIVERGED_ITS - Ran out of iterations before any convergence criteria was
452: reached
454: Level: beginner
456: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
458: M*/
460: /*MC
461: KSP_CONVERGED_ITS - Used by the KSPPREONLY solver after the single iteration of
462: the preconditioner is applied. Also used when the KSPSkipConverged() convergence
463: test routine is set in KSP.
466: Level: beginner
469: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
471: M*/
473: /*MC
474: KSP_DIVERGED_BREAKDOWN - A breakdown in the Krylov method was detected so the
475: method could not continue to enlarge the Krylov space. Could be due to a singlular matrix or
476: preconditioner.
478: Level: beginner
480: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
482: M*/
484: /*MC
485: KSP_DIVERGED_BREAKDOWN_BICG - A breakdown in the KSPBICG method was detected so the
486: method could not continue to enlarge the Krylov space.
489: Level: beginner
492: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
494: M*/
496: /*MC
497: KSP_DIVERGED_NONSYMMETRIC - It appears the operator or preconditioner is not
498: symmetric and this Krylov method (KSPCG, KSPMINRES, KSPCR) requires symmetry
500: Level: beginner
502: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
504: M*/
506: /*MC
507: KSP_DIVERGED_INDEFINITE_PC - It appears the preconditioner is indefinite (has both
508: positive and negative eigenvalues) and this Krylov method (KSPCG) requires it to
509: be positive definite
511: Level: beginner
513: Notes: This can happen with the PCICC preconditioner, use -pc_factor_shift_positive_definite to force
514: the PCICC preconditioner to generate a positive definite preconditioner
516: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
518: M*/
520: /*MC
521: KSP_CONVERGED_ITERATING - This flag is returned if you call KSPGetConvergedReason()
522: while the KSPSolve() is still running.
524: Level: beginner
526: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
528: M*/
530: extern PetscErrorCode KSPSetConvergenceTest(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*),void *,PetscErrorCode (*)(void*));
531: extern PetscErrorCode KSPGetConvergenceContext(KSP,void **);
532: extern PetscErrorCode KSPDefaultConverged(KSP,PetscInt,PetscReal,KSPConvergedReason*,void *);
533: extern PetscErrorCode KSPConvergedLSQR(KSP,PetscInt,PetscReal,KSPConvergedReason*,void *);
534: extern PetscErrorCode KSPDefaultConvergedDestroy(void *);
535: extern PetscErrorCode KSPDefaultConvergedCreate(void **);
536: extern PetscErrorCode KSPDefaultConvergedSetUIRNorm(KSP);
537: extern PetscErrorCode KSPDefaultConvergedSetUMIRNorm(KSP);
538: extern PetscErrorCode KSPSkipConverged(KSP,PetscInt,PetscReal,KSPConvergedReason*,void *);
539: extern PetscErrorCode KSPGetConvergedReason(KSP,KSPConvergedReason *);
541: extern PetscErrorCode KSPComputeExplicitOperator(KSP,Mat *);
543: /*E
544: KSPCGType - Determines what type of CG to use
546: Level: beginner
548: .seealso: KSPCGSetType()
549: E*/
550: typedef enum {KSP_CG_SYMMETRIC=0,KSP_CG_HERMITIAN=1} KSPCGType;
551: extern const char *KSPCGTypes[];
553: extern PetscErrorCode KSPCGSetType(KSP,KSPCGType);
554: extern PetscErrorCode KSPCGUseSingleReduction(KSP,PetscBool );
556: extern PetscErrorCode KSPNASHSetRadius(KSP,PetscReal);
557: extern PetscErrorCode KSPNASHGetNormD(KSP,PetscReal *);
558: extern PetscErrorCode KSPNASHGetObjFcn(KSP,PetscReal *);
560: extern PetscErrorCode KSPSTCGSetRadius(KSP,PetscReal);
561: extern PetscErrorCode KSPSTCGGetNormD(KSP,PetscReal *);
562: extern PetscErrorCode KSPSTCGGetObjFcn(KSP,PetscReal *);
564: extern PetscErrorCode KSPGLTRSetRadius(KSP,PetscReal);
565: extern PetscErrorCode KSPGLTRGetNormD(KSP,PetscReal *);
566: extern PetscErrorCode KSPGLTRGetObjFcn(KSP,PetscReal *);
567: extern PetscErrorCode KSPGLTRGetMinEig(KSP,PetscReal *);
568: extern PetscErrorCode KSPGLTRGetLambda(KSP,PetscReal *);
570: extern PetscErrorCode KSPPythonSetType(KSP,const char[]);
572: extern PetscErrorCode PCPreSolve(PC,KSP);
573: extern PetscErrorCode PCPostSolve(PC,KSP);
575: extern PetscErrorCode KSPMonitorLGCreate(const char[],const char[],int,int,int,int,PetscDrawLG*);
576: extern PetscErrorCode KSPMonitorLG(KSP,PetscInt,PetscReal,void*);
577: extern PetscErrorCode KSPMonitorLGDestroy(PetscDrawLG*);
578: extern PetscErrorCode KSPMonitorLGTrueResidualNormCreate(MPI_Comm,const char[],const char[],int,int,int,int,PetscDrawLG*);
579: extern PetscErrorCode KSPMonitorLGTrueResidualNorm(KSP,PetscInt,PetscReal,void*);
580: extern PetscErrorCode KSPMonitorLGTrueResidualNormDestroy(PetscDrawLG*);
581: extern PetscErrorCode KSPMonitorLGRangeCreate(const char[],const char[],int,int,int,int,PetscDrawLG*);
582: extern PetscErrorCode KSPMonitorLGRange(KSP,PetscInt,PetscReal,void*);
583: extern PetscErrorCode KSPMonitorLGRangeDestroy(PetscDrawLG*);
585: extern PetscErrorCode PCShellSetPreSolve(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec));
586: extern PetscErrorCode PCShellSetPostSolve(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec));
588: /* see src/ksp/ksp/interface/iguess.c */
589: typedef struct _p_KSPFischerGuess {PetscInt method,curl,maxl,refcnt;PetscBool monitor;Mat mat; KSP ksp;}* KSPFischerGuess;
591: extern PetscErrorCode KSPFischerGuessCreate(KSP,PetscInt,PetscInt,KSPFischerGuess*);
592: extern PetscErrorCode KSPFischerGuessDestroy(KSPFischerGuess*);
593: extern PetscErrorCode KSPFischerGuessReset(KSPFischerGuess);
594: extern PetscErrorCode KSPFischerGuessUpdate(KSPFischerGuess,Vec);
595: extern PetscErrorCode KSPFischerGuessFormGuess(KSPFischerGuess,Vec,Vec);
596: extern PetscErrorCode KSPFischerGuessSetFromOptions(KSPFischerGuess);
598: extern PetscErrorCode KSPSetUseFischerGuess(KSP,PetscInt,PetscInt);
599: extern PetscErrorCode KSPSetFischerGuess(KSP,KSPFischerGuess);
600: extern PetscErrorCode KSPGetFischerGuess(KSP,KSPFischerGuess*);
602: extern PetscErrorCode MatCreateSchurComplement(Mat,Mat,Mat,Mat,Mat,Mat*);
603: extern PetscErrorCode MatSchurComplementGetKSP(Mat,KSP*);
604: extern PetscErrorCode MatSchurComplementUpdate(Mat,Mat,Mat,Mat,Mat,Mat,MatStructure);
605: extern PetscErrorCode MatSchurComplementGetSubmatrices(Mat,Mat*,Mat*,Mat*,Mat*,Mat*);
606: extern PetscErrorCode MatGetSchurComplement(Mat,IS,IS,IS,IS,MatReuse,Mat *,MatReuse,Mat *);
608: extern PetscErrorCode MatGetSchurComplement_Basic(Mat mat,IS isrow0,IS iscol0,IS isrow1,IS iscol1,MatReuse mreuse,Mat *newmat,MatReuse preuse,Mat *newpmat);
610: extern PetscErrorCode KSPSetDM(KSP,DM);
611: extern PetscErrorCode KSPSetDMActive(KSP,PetscBool );
612: extern PetscErrorCode KSPGetDM(KSP,DM*);
613: extern PetscErrorCode KSPSetApplicationContext(KSP,void*);
614: extern PetscErrorCode KSPGetApplicationContext(KSP,void*);
616: PETSC_EXTERN_CXX_END
617: #endif