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 KSPSPECEST "specest"
62: /* Logging support */
63: extern PetscClassId KSP_CLASSID;
65: extern PetscErrorCode KSPCreate(MPI_Comm,KSP *);
66: extern PetscErrorCode KSPSetType(KSP,const KSPType);
67: extern PetscErrorCode KSPSetUp(KSP);
68: extern PetscErrorCode KSPSetUpOnBlocks(KSP);
69: extern PetscErrorCode KSPSolve(KSP,Vec,Vec);
70: extern PetscErrorCode KSPSolveTranspose(KSP,Vec,Vec);
71: extern PetscErrorCode KSPReset(KSP);
72: extern PetscErrorCode KSPDestroy(KSP*);
74: extern PetscFList KSPList;
75: extern PetscBool KSPRegisterAllCalled;
76: extern PetscErrorCode KSPRegisterAll(const char[]);
77: extern PetscErrorCode KSPRegisterDestroy(void);
78: extern PetscErrorCode KSPRegister(const char[],const char[],const char[],PetscErrorCode (*)(KSP));
80: /*MC
81: KSPRegisterDynamic - Adds a method to the Krylov subspace solver package.
83: Synopsis:
84: PetscErrorCode KSPRegisterDynamic(const char *name_solver,const char *path,const char *name_create,PetscErrorCode (*routine_create)(KSP))
86: Not Collective
88: Input Parameters:
89: + name_solver - name of a new user-defined solver
90: . path - path (either absolute or relative) the library containing this solver
91: . name_create - name of routine to create method context
92: - routine_create - routine to create method context
94: Notes:
95: KSPRegisterDynamic() may be called multiple times to add several user-defined solvers.
97: If dynamic libraries are used, then the fourth input argument (routine_create)
98: is ignored.
100: Sample usage:
101: .vb
102: KSPRegisterDynamic("my_solver",/home/username/my_lib/lib/libO/solaris/mylib.a,
103: "MySolverCreate",MySolverCreate);
104: .ve
106: Then, your solver can be chosen with the procedural interface via
107: $ KSPSetType(ksp,"my_solver")
108: or at runtime via the option
109: $ -ksp_type my_solver
111: Level: advanced
113: Notes: Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR},
114: and others of the form ${any_environmental_variable} occuring in pathname will be
115: replaced with appropriate values.
116: If your function is not being put into a shared library then use KSPRegister() instead
118: .keywords: KSP, register
120: .seealso: KSPRegisterAll(), KSPRegisterDestroy()
122: M*/
123: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
124: #define KSPRegisterDynamic(a,b,c,d) KSPRegister(a,b,c,0)
125: #else
126: #define KSPRegisterDynamic(a,b,c,d) KSPRegister(a,b,c,d)
127: #endif
129: extern PetscErrorCode KSPGetType(KSP,const KSPType *);
130: extern PetscErrorCode KSPSetPCSide(KSP,PCSide);
131: extern PetscErrorCode KSPGetPCSide(KSP,PCSide*);
132: extern PetscErrorCode KSPGetTolerances(KSP,PetscReal*,PetscReal*,PetscReal*,PetscInt*);
133: extern PetscErrorCode KSPSetTolerances(KSP,PetscReal,PetscReal,PetscReal,PetscInt);
134: extern PetscErrorCode KSPSetInitialGuessNonzero(KSP,PetscBool );
135: extern PetscErrorCode KSPGetInitialGuessNonzero(KSP,PetscBool *);
136: extern PetscErrorCode KSPSetInitialGuessKnoll(KSP,PetscBool );
137: extern PetscErrorCode KSPGetInitialGuessKnoll(KSP,PetscBool *);
138: extern PetscErrorCode KSPSetErrorIfNotConverged(KSP,PetscBool );
139: extern PetscErrorCode KSPGetErrorIfNotConverged(KSP,PetscBool *);
140: extern PetscErrorCode KSPGetComputeEigenvalues(KSP,PetscBool *);
141: extern PetscErrorCode KSPSetComputeEigenvalues(KSP,PetscBool );
142: extern PetscErrorCode KSPGetComputeSingularValues(KSP,PetscBool *);
143: extern PetscErrorCode KSPSetComputeSingularValues(KSP,PetscBool );
144: extern PetscErrorCode KSPGetRhs(KSP,Vec *);
145: extern PetscErrorCode KSPGetSolution(KSP,Vec *);
146: extern PetscErrorCode KSPGetResidualNorm(KSP,PetscReal*);
147: extern PetscErrorCode KSPGetIterationNumber(KSP,PetscInt*);
148: extern PetscErrorCode KSPSetNullSpace(KSP,MatNullSpace);
149: extern PetscErrorCode KSPGetNullSpace(KSP,MatNullSpace*);
150: extern PetscErrorCode KSPGetVecs(KSP,PetscInt,Vec**,PetscInt,Vec**);
152: extern PetscErrorCode KSPSetPC(KSP,PC);
153: extern PetscErrorCode KSPGetPC(KSP,PC*);
155: extern PetscErrorCode KSPMonitor(KSP,PetscInt,PetscReal);
156: extern PetscErrorCode KSPMonitorSet(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscReal,void*),void *,PetscErrorCode (*)(void**));
157: extern PetscErrorCode KSPMonitorCancel(KSP);
158: extern PetscErrorCode KSPGetMonitorContext(KSP,void **);
159: extern PetscErrorCode KSPGetResidualHistory(KSP,PetscReal*[],PetscInt *);
160: extern PetscErrorCode KSPSetResidualHistory(KSP,PetscReal[],PetscInt,PetscBool );
162: /* not sure where to put this */
163: extern PetscErrorCode PCKSPGetKSP(PC,KSP*);
164: extern PetscErrorCode PCBJacobiGetSubKSP(PC,PetscInt*,PetscInt*,KSP*[]);
165: extern PetscErrorCode PCASMGetSubKSP(PC,PetscInt*,PetscInt*,KSP*[]);
166: extern PetscErrorCode PCGASMGetSubKSP(PC,PetscInt*,PetscInt*,KSP*[]);
167: extern PetscErrorCode PCFieldSplitGetSubKSP(PC,PetscInt*,KSP*[]);
169: extern PetscErrorCode PCGalerkinGetKSP(PC,KSP *);
171: extern PetscErrorCode KSPBuildSolution(KSP,Vec,Vec *);
172: extern PetscErrorCode KSPBuildResidual(KSP,Vec,Vec,Vec *);
174: extern PetscErrorCode KSPRichardsonSetScale(KSP,PetscReal);
175: extern PetscErrorCode KSPRichardsonSetSelfScale(KSP,PetscBool );
176: extern PetscErrorCode KSPChebychevSetEigenvalues(KSP,PetscReal,PetscReal);
177: extern PetscErrorCode KSPChebychevSetEstimateEigenvalues(KSP,PetscReal,PetscReal,PetscReal,PetscReal);
178: extern PetscErrorCode KSPComputeExtremeSingularValues(KSP,PetscReal*,PetscReal*);
179: extern PetscErrorCode KSPComputeEigenvalues(KSP,PetscInt,PetscReal*,PetscReal*,PetscInt *);
180: extern PetscErrorCode KSPComputeEigenvaluesExplicitly(KSP,PetscInt,PetscReal*,PetscReal*);
182: extern PetscErrorCode KSPGMRESSetRestart(KSP, PetscInt);
183: extern PetscErrorCode KSPGMRESGetRestart(KSP, PetscInt*);
184: extern PetscErrorCode KSPGMRESSetHapTol(KSP,PetscReal);
186: extern PetscErrorCode KSPGMRESSetPreAllocateVectors(KSP);
187: extern PetscErrorCode KSPGMRESSetOrthogonalization(KSP,PetscErrorCode (*)(KSP,PetscInt));
188: extern PetscErrorCode KSPGMRESGetOrthogonalization(KSP,PetscErrorCode (**)(KSP,PetscInt));
189: extern PetscErrorCode KSPGMRESModifiedGramSchmidtOrthogonalization(KSP,PetscInt);
190: extern PetscErrorCode KSPGMRESClassicalGramSchmidtOrthogonalization(KSP,PetscInt);
192: extern PetscErrorCode KSPLGMRESSetAugDim(KSP,PetscInt);
193: extern PetscErrorCode KSPLGMRESSetConstant(KSP);
195: extern PetscErrorCode KSPGCRSetRestart(KSP,PetscInt);
196: extern PetscErrorCode KSPGCRGetRestart(KSP,PetscInt*);
197: extern PetscErrorCode KSPGCRSetModifyPC(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscReal,void*),void*,PetscErrorCode(*)(void*));
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;
209: extern const char *KSPGMRESCGSRefinementTypes[];
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*/
252: extern PetscErrorCode KSPGMRESSetCGSRefinementType(KSP,KSPGMRESCGSRefinementType);
253: extern PetscErrorCode KSPGMRESGetCGSRefinementType(KSP,KSPGMRESCGSRefinementType*);
255: extern PetscErrorCode KSPFGMRESModifyPCNoChange(KSP,PetscInt,PetscInt,PetscReal,void*);
256: extern PetscErrorCode KSPFGMRESModifyPCKSP(KSP,PetscInt,PetscInt,PetscReal,void*);
257: extern PetscErrorCode KSPFGMRESSetModifyPC(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscInt,PetscReal,void*),void*,PetscErrorCode(*)(void*));
259: extern PetscErrorCode KSPQCGSetTrustRegionRadius(KSP,PetscReal);
260: extern PetscErrorCode KSPQCGGetQuadratic(KSP,PetscReal*);
261: extern PetscErrorCode KSPQCGGetTrialStepNorm(KSP,PetscReal*);
263: extern PetscErrorCode KSPBCGSLSetXRes(KSP,PetscReal);
264: extern PetscErrorCode KSPBCGSLSetPol(KSP,PetscBool );
265: extern PetscErrorCode KSPBCGSLSetEll(KSP,PetscInt);
267: extern PetscErrorCode KSPSetFromOptions(KSP);
268: extern PetscErrorCode KSPAddOptionsChecker(PetscErrorCode (*)(KSP));
270: extern PetscErrorCode KSPMonitorSingularValue(KSP,PetscInt,PetscReal,void *);
271: extern PetscErrorCode KSPMonitorDefault(KSP,PetscInt,PetscReal,void *);
272: extern PetscErrorCode KSPLSQRMonitorDefault(KSP,PetscInt,PetscReal,void *);
273: extern PetscErrorCode KSPMonitorRange(KSP,PetscInt,PetscReal,void *);
274: extern PetscErrorCode KSPMonitorTrueResidualNorm(KSP,PetscInt,PetscReal,void *);
275: extern PetscErrorCode KSPMonitorDefaultShort(KSP,PetscInt,PetscReal,void *);
276: extern PetscErrorCode KSPMonitorSolution(KSP,PetscInt,PetscReal,void *);
277: extern PetscErrorCode KSPGMRESMonitorKrylov(KSP,PetscInt,PetscReal,void *);
279: extern PetscErrorCode KSPUnwindPreconditioner(KSP,Vec,Vec);
280: extern PetscErrorCode KSPDefaultBuildSolution(KSP,Vec,Vec*);
281: extern PetscErrorCode KSPDefaultBuildResidual(KSP,Vec,Vec,Vec *);
282: extern PetscErrorCode KSPInitialResidual(KSP,Vec,Vec,Vec,Vec,Vec);
284: extern PetscErrorCode KSPSetOperators(KSP,Mat,Mat,MatStructure);
285: extern PetscErrorCode KSPGetOperators(KSP,Mat*,Mat*,MatStructure*);
286: extern PetscErrorCode KSPGetOperatorsSet(KSP,PetscBool *,PetscBool *);
287: extern PetscErrorCode KSPSetOptionsPrefix(KSP,const char[]);
288: extern PetscErrorCode KSPAppendOptionsPrefix(KSP,const char[]);
289: extern PetscErrorCode KSPGetOptionsPrefix(KSP,const char*[]);
291: extern PetscErrorCode KSPSetDiagonalScale(KSP,PetscBool );
292: extern PetscErrorCode KSPGetDiagonalScale(KSP,PetscBool *);
293: extern PetscErrorCode KSPSetDiagonalScaleFix(KSP,PetscBool );
294: extern PetscErrorCode KSPGetDiagonalScaleFix(KSP,PetscBool *);
296: extern PetscErrorCode KSPView(KSP,PetscViewer);
298: extern PetscErrorCode KSPLSQRSetStandardErrorVec(KSP,Vec);
299: extern PetscErrorCode KSPLSQRGetStandardErrorVec(KSP,Vec*);
301: extern PetscErrorCode PCRedundantGetKSP(PC,KSP*);
302: extern PetscErrorCode PCRedistributeGetKSP(PC,KSP*);
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)
320: extern const char *const*const KSPNormTypes;
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*/
361: extern PetscErrorCode KSPSetNormType(KSP,KSPNormType);
362: extern PetscErrorCode KSPGetNormType(KSP,KSPNormType*);
363: extern PetscErrorCode KSPSetSupportedNorm(KSP ksp,KSPNormType,PCSide,PetscInt);
364: extern PetscErrorCode KSPSetCheckNormIteration(KSP,PetscInt);
365: extern PetscErrorCode KSPSetLagNorm(KSP,PetscBool );
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;
404: extern const char *const*KSPConvergedReasons;
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*/
529: extern PetscErrorCode KSPSetConvergenceTest(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*),void *,PetscErrorCode (*)(void*));
530: extern PetscErrorCode KSPGetConvergenceContext(KSP,void **);
531: extern PetscErrorCode KSPDefaultConverged(KSP,PetscInt,PetscReal,KSPConvergedReason*,void *);
532: extern PetscErrorCode KSPConvergedLSQR(KSP,PetscInt,PetscReal,KSPConvergedReason*,void *);
533: extern PetscErrorCode KSPDefaultConvergedDestroy(void *);
534: extern PetscErrorCode KSPDefaultConvergedCreate(void **);
535: extern PetscErrorCode KSPDefaultConvergedSetUIRNorm(KSP);
536: extern PetscErrorCode KSPDefaultConvergedSetUMIRNorm(KSP);
537: extern PetscErrorCode KSPSkipConverged(KSP,PetscInt,PetscReal,KSPConvergedReason*,void *);
538: extern PetscErrorCode KSPGetConvergedReason(KSP,KSPConvergedReason *);
540: extern PetscErrorCode KSPComputeExplicitOperator(KSP,Mat *);
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;
550: extern const char *KSPCGTypes[];
552: extern PetscErrorCode KSPCGSetType(KSP,KSPCGType);
553: extern PetscErrorCode KSPCGUseSingleReduction(KSP,PetscBool );
555: extern PetscErrorCode KSPNASHSetRadius(KSP,PetscReal);
556: extern PetscErrorCode KSPNASHGetNormD(KSP,PetscReal *);
557: extern PetscErrorCode KSPNASHGetObjFcn(KSP,PetscReal *);
559: extern PetscErrorCode KSPSTCGSetRadius(KSP,PetscReal);
560: extern PetscErrorCode KSPSTCGGetNormD(KSP,PetscReal *);
561: extern PetscErrorCode KSPSTCGGetObjFcn(KSP,PetscReal *);
563: extern PetscErrorCode KSPGLTRSetRadius(KSP,PetscReal);
564: extern PetscErrorCode KSPGLTRGetNormD(KSP,PetscReal *);
565: extern PetscErrorCode KSPGLTRGetObjFcn(KSP,PetscReal *);
566: extern PetscErrorCode KSPGLTRGetMinEig(KSP,PetscReal *);
567: extern PetscErrorCode KSPGLTRGetLambda(KSP,PetscReal *);
569: extern PetscErrorCode KSPPythonSetType(KSP,const char[]);
571: extern PetscErrorCode PCPreSolve(PC,KSP);
572: extern PetscErrorCode PCPostSolve(PC,KSP);
574: extern PetscErrorCode KSPMonitorLGCreate(const char[],const char[],int,int,int,int,PetscDrawLG*);
575: extern PetscErrorCode KSPMonitorLG(KSP,PetscInt,PetscReal,void*);
576: extern PetscErrorCode KSPMonitorLGDestroy(PetscDrawLG*);
577: extern PetscErrorCode KSPMonitorLGTrueResidualNormCreate(MPI_Comm,const char[],const char[],int,int,int,int,PetscDrawLG*);
578: extern PetscErrorCode KSPMonitorLGTrueResidualNorm(KSP,PetscInt,PetscReal,void*);
579: extern PetscErrorCode KSPMonitorLGTrueResidualNormDestroy(PetscDrawLG*);
580: extern PetscErrorCode KSPMonitorLGRangeCreate(const char[],const char[],int,int,int,int,PetscDrawLG*);
581: extern PetscErrorCode KSPMonitorLGRange(KSP,PetscInt,PetscReal,void*);
582: extern PetscErrorCode KSPMonitorLGRangeDestroy(PetscDrawLG*);
584: extern PetscErrorCode PCShellSetPreSolve(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec));
585: extern PetscErrorCode PCShellSetPostSolve(PC,PetscErrorCode (*)(PC,KSP,Vec,Vec));
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;
590: extern PetscErrorCode KSPFischerGuessCreate(KSP,PetscInt,PetscInt,KSPFischerGuess*);
591: extern PetscErrorCode KSPFischerGuessDestroy(KSPFischerGuess*);
592: extern PetscErrorCode KSPFischerGuessReset(KSPFischerGuess);
593: extern PetscErrorCode KSPFischerGuessUpdate(KSPFischerGuess,Vec);
594: extern PetscErrorCode KSPFischerGuessFormGuess(KSPFischerGuess,Vec,Vec);
595: extern PetscErrorCode KSPFischerGuessSetFromOptions(KSPFischerGuess);
597: extern PetscErrorCode KSPSetUseFischerGuess(KSP,PetscInt,PetscInt);
598: extern PetscErrorCode KSPSetFischerGuess(KSP,KSPFischerGuess);
599: extern PetscErrorCode KSPGetFischerGuess(KSP,KSPFischerGuess*);
601: extern PetscErrorCode MatCreateSchurComplement(Mat,Mat,Mat,Mat,Mat,Mat*);
602: extern PetscErrorCode MatSchurComplementGetKSP(Mat,KSP*);
603: extern PetscErrorCode MatSchurComplementUpdate(Mat,Mat,Mat,Mat,Mat,Mat,MatStructure);
604: extern PetscErrorCode MatSchurComplementGetSubmatrices(Mat,Mat*,Mat*,Mat*,Mat*,Mat*);
605: extern PetscErrorCode MatGetSchurComplement(Mat,IS,IS,IS,IS,MatReuse,Mat *,MatReuse,Mat *);
607: extern PetscErrorCode MatGetSchurComplement_Basic(Mat mat,IS isrow0,IS iscol0,IS isrow1,IS iscol1,MatReuse mreuse,Mat *newmat,MatReuse preuse,Mat *newpmat);
609: extern PetscErrorCode KSPSetDM(KSP,DM);
610: extern PetscErrorCode KSPSetDMActive(KSP,PetscBool );
611: extern PetscErrorCode KSPGetDM(KSP,DM*);
612: extern PetscErrorCode KSPSetApplicationContext(KSP,void*);
613: extern PetscErrorCode KSPGetApplicationContext(KSP,void*);
615: PETSC_EXTERN_CXX_END
616: #endif