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
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: /*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 const 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 KSPTCQMR "tcqmr"
43: #define KSPBCGS "bcgs"
44: #define KSPBCGSL "bcgsl"
45: #define KSPCGS "cgs"
46: #define KSPTFQMR "tfqmr"
47: #define KSPCR "cr"
48: #define KSPLSQR "lsqr"
49: #define KSPPREONLY "preonly"
50: #define KSPQCG "qcg"
51: #define KSPBICG "bicg"
52: #define KSPMINRES "minres"
53: #define KSPSYMMLQ "symmlq"
54: #define KSPLCD "lcd"
56: /* Logging support */
59: EXTERN PetscErrorCode KSPCreate(MPI_Comm,KSP *);
60: EXTERN PetscErrorCode KSPSetType(KSP,KSPType);
61: EXTERN PetscErrorCode KSPSetUp(KSP);
62: EXTERN PetscErrorCode KSPSetUpOnBlocks(KSP);
63: EXTERN PetscErrorCode KSPSolve(KSP,Vec,Vec);
64: EXTERN PetscErrorCode KSPSolveTranspose(KSP,Vec,Vec);
65: EXTERN PetscErrorCode KSPDestroy(KSP);
68: EXTERN PetscErrorCode KSPRegisterAll(const char[]);
69: EXTERN PetscErrorCode KSPRegisterDestroy(void);
70: EXTERN PetscErrorCode KSPRegister(const char[],const char[],const char[],PetscErrorCode (*)(KSP));
72: /*MC
73: KSPRegisterDynamic - Adds a method to the Krylov subspace solver package.
75: Synopsis:
76: PetscErrorCode KSPRegisterDynamic(char *name_solver,char *path,char *name_create,PetscErrorCode (*routine_create)(KSP))
78: Not Collective
80: Input Parameters:
81: + name_solver - name of a new user-defined solver
82: . path - path (either absolute or relative) the library containing this solver
83: . name_create - name of routine to create method context
84: - routine_create - routine to create method context
86: Notes:
87: KSPRegisterDynamic() may be called multiple times to add several user-defined solvers.
89: If dynamic libraries are used, then the fourth input argument (routine_create)
90: is ignored.
92: Sample usage:
93: .vb
94: KSPRegisterDynamic("my_solver",/home/username/my_lib/lib/libO/solaris/mylib.a,
95: "MySolverCreate",MySolverCreate);
96: .ve
98: Then, your solver can be chosen with the procedural interface via
99: $ KSPSetType(ksp,"my_solver")
100: or at runtime via the option
101: $ -ksp_type my_solver
103: Level: advanced
105: Notes: Environmental variables such as ${PETSC_ARCH}, ${PETSC_DIR}, ${PETSC_LIB_DIR},
106: and others of the form ${any_environmental_variable} occuring in pathname will be
107: replaced with appropriate values.
108: If your function is not being put into a shared library then use KSPRegister() instead
110: .keywords: KSP, register
112: .seealso: KSPRegisterAll(), KSPRegisterDestroy()
114: M*/
115: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
116: #define KSPRegisterDynamic(a,b,c,d) KSPRegister(a,b,c,0)
117: #else
118: #define KSPRegisterDynamic(a,b,c,d) KSPRegister(a,b,c,d)
119: #endif
121: EXTERN PetscErrorCode KSPGetType(KSP,KSPType *);
122: EXTERN PetscErrorCode KSPSetPreconditionerSide(KSP,PCSide);
123: EXTERN PetscErrorCode KSPGetPreconditionerSide(KSP,PCSide*);
124: EXTERN PetscErrorCode KSPGetTolerances(KSP,PetscReal*,PetscReal*,PetscReal*,PetscInt*);
125: EXTERN PetscErrorCode KSPSetTolerances(KSP,PetscReal,PetscReal,PetscReal,PetscInt);
126: EXTERN PetscErrorCode KSPSetInitialGuessNonzero(KSP,PetscTruth);
127: EXTERN PetscErrorCode KSPGetInitialGuessNonzero(KSP,PetscTruth *);
128: EXTERN PetscErrorCode KSPSetInitialGuessKnoll(KSP,PetscTruth);
129: EXTERN PetscErrorCode KSPGetInitialGuessKnoll(KSP,PetscTruth*);
130: EXTERN PetscErrorCode KSPGetComputeEigenvalues(KSP,PetscTruth*);
131: EXTERN PetscErrorCode KSPSetComputeEigenvalues(KSP,PetscTruth);
132: EXTERN PetscErrorCode KSPGetComputeSingularValues(KSP,PetscTruth*);
133: EXTERN PetscErrorCode KSPSetComputeSingularValues(KSP,PetscTruth);
134: EXTERN PetscErrorCode KSPGetRhs(KSP,Vec *);
135: EXTERN PetscErrorCode KSPGetSolution(KSP,Vec *);
136: EXTERN PetscErrorCode KSPGetResidualNorm(KSP,PetscReal*);
137: EXTERN PetscErrorCode KSPGetIterationNumber(KSP,PetscInt*);
138: EXTERN PetscErrorCode KSPSetNullSpace(KSP,MatNullSpace);
139: EXTERN PetscErrorCode KSPGetNullSpace(KSP,MatNullSpace*);
140: EXTERN PetscErrorCode KSPGetVecs(KSP,PetscInt,Vec**,PetscInt,Vec**);
142: EXTERN PetscErrorCode KSPSetPC(KSP,PC);
143: EXTERN PetscErrorCode KSPGetPC(KSP,PC*);
145: EXTERN PetscErrorCode KSPMonitorSet(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscReal,void*),void *,PetscErrorCode (*)(void*));
146: EXTERN PetscErrorCode KSPMonitorCancel(KSP);
147: EXTERN PetscErrorCode KSPGetMonitorContext(KSP,void **);
148: EXTERN PetscErrorCode KSPGetResidualHistory(KSP,PetscReal*[],PetscInt *);
149: EXTERN PetscErrorCode KSPSetResidualHistory(KSP,PetscReal[],PetscInt,PetscTruth);
151: /* not sure where to put this */
152: EXTERN PetscErrorCode PCKSPGetKSP(PC,KSP*);
153: EXTERN PetscErrorCode PCBJacobiGetSubKSP(PC,PetscInt*,PetscInt*,KSP*[]);
154: EXTERN PetscErrorCode PCASMGetSubKSP(PC,PetscInt*,PetscInt*,KSP*[]);
155: EXTERN PetscErrorCode PCFieldSplitGetSubKSP(PC,PetscInt*,KSP*[]);
157: EXTERN PetscErrorCode PCGalerkinGetKSP(PC,KSP *);
159: EXTERN PetscErrorCode KSPBuildSolution(KSP,Vec,Vec *);
160: EXTERN PetscErrorCode KSPBuildResidual(KSP,Vec,Vec,Vec *);
162: EXTERN PetscErrorCode KSPRichardsonSetScale(KSP,PetscReal);
163: EXTERN PetscErrorCode KSPChebychevSetEigenvalues(KSP,PetscReal,PetscReal);
164: EXTERN PetscErrorCode KSPComputeExtremeSingularValues(KSP,PetscReal*,PetscReal*);
165: EXTERN PetscErrorCode KSPComputeEigenvalues(KSP,PetscInt,PetscReal*,PetscReal*,PetscInt *);
166: EXTERN PetscErrorCode KSPComputeEigenvaluesExplicitly(KSP,PetscInt,PetscReal*,PetscReal*);
168: EXTERN PetscErrorCode KSPGMRESSetRestart(KSP, PetscInt);
169: EXTERN PetscErrorCode KSPGMRESSetHapTol(KSP,PetscReal);
171: EXTERN PetscErrorCode KSPGMRESSetPreAllocateVectors(KSP);
172: EXTERN PetscErrorCode KSPGMRESSetOrthogonalization(KSP,PetscErrorCode (*)(KSP,PetscInt));
173: EXTERN PetscErrorCode KSPGMRESModifiedGramSchmidtOrthogonalization(KSP,PetscInt);
174: EXTERN PetscErrorCode KSPGMRESClassicalGramSchmidtOrthogonalization(KSP,PetscInt);
176: EXTERN PetscErrorCode KSPLGMRESSetAugDim(KSP,PetscInt);
177: EXTERN PetscErrorCode KSPLGMRESSetConstant(KSP);
179: /*E
180: KSPGMRESCGSRefinementType - How the classical (unmodified) Gram-Schmidt is performed.
182: Level: advanced
184: .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(),
185: KSPGMRESSetCGSRefinementType(), KSPGMRESModifiedGramSchmidtOrthogonalization()
187: E*/
188: typedef enum {KSP_GMRES_CGS_REFINE_NEVER, KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS} KSPGMRESCGSRefinementType;
190: /*MC
191: KSP_GMRES_CGS_REFINE_NEVER - Just do the classical (unmodified) Gram-Schmidt process
193: Level: advanced
195: Note: Possible unstable, but the fastest to compute
197: .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(),
198: KSPGMRESSetCGSRefinementType(), KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS,
199: KSPGMRESModifiedGramSchmidtOrthogonalization()
200: M*/
202: /*MC
203: KSP_GMRES_CGS_REFINE_IFNEEDED - Do the classical (unmodified) Gram-Schmidt process and one step of
204: iterative refinement if an estimate of the orthogonality of the resulting vectors indicates
205: poor orthogonality.
207: Level: advanced
209: Note: This is slower than KSP_GMRES_CGS_REFINE_NEVER because it requires an extra norm computation to
210: estimate the orthogonality but is more stable.
212: .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(),
213: KSPGMRESSetCGSRefinementType(), KSP_GMRES_CGS_REFINE_NEVER, KSP_GMRES_CGS_REFINE_ALWAYS,
214: KSPGMRESModifiedGramSchmidtOrthogonalization()
215: M*/
217: /*MC
218: KSP_GMRES_CGS_REFINE_NEVER - Do two steps of the classical (unmodified) Gram-Schmidt process.
220: Level: advanced
222: Note: This is roughly twice the cost of KSP_GMRES_CGS_REFINE_NEVER because it performs the process twice
223: but it saves the extra norm calculation needed by KSP_GMRES_CGS_REFINE_IFNEEDED.
225: You should only use this if you absolutely know that the iterative refinement is needed.
227: .seealso: KSPGMRESClassicalGramSchmidtOrthogonalization(), KSPGMRESSetOrthogonalization(),
228: KSPGMRESSetCGSRefinementType(), KSP_GMRES_CGS_REFINE_IFNEEDED, KSP_GMRES_CGS_REFINE_ALWAYS,
229: KSPGMRESModifiedGramSchmidtOrthogonalization()
230: M*/
232: EXTERN PetscErrorCode KSPGMRESSetCGSRefinementType(KSP,KSPGMRESCGSRefinementType);
234: EXTERN PetscErrorCode KSPFGMRESModifyPCNoChange(KSP,PetscInt,PetscInt,PetscReal,void*);
235: EXTERN PetscErrorCode KSPFGMRESModifyPCKSP(KSP,PetscInt,PetscInt,PetscReal,void*);
236: EXTERN PetscErrorCode KSPFGMRESSetModifyPC(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscInt,PetscReal,void*),void*,PetscErrorCode(*)(void*));
238: EXTERN PetscErrorCode KSPQCGSetTrustRegionRadius(KSP,PetscReal);
239: EXTERN PetscErrorCode KSPQCGGetQuadratic(KSP,PetscReal*);
240: EXTERN PetscErrorCode KSPQCGGetTrialStepNorm(KSP,PetscReal*);
242: EXTERN PetscErrorCode KSPBCGSLSetXRes(KSP,PetscReal);
243: EXTERN PetscErrorCode KSPBCGSLSetPol(KSP,PetscTruth);
244: EXTERN PetscErrorCode KSPBCGSLSetEll(KSP,int);
246: EXTERN PetscErrorCode KSPSetFromOptions(KSP);
247: EXTERN PetscErrorCode KSPAddOptionsChecker(PetscErrorCode (*)(KSP));
249: EXTERN PetscErrorCode KSPMonitorSingularValue(KSP,PetscInt,PetscReal,void *);
250: EXTERN PetscErrorCode KSPMonitorDefault(KSP,PetscInt,PetscReal,void *);
251: EXTERN PetscErrorCode KSPMonitorTrueResidualNorm(KSP,PetscInt,PetscReal,void *);
252: EXTERN PetscErrorCode KSPMonitorDefaultShort(KSP,PetscInt,PetscReal,void *);
253: EXTERN PetscErrorCode KSPMonitorSolution(KSP,PetscInt,PetscReal,void *);
254: EXTERN PetscErrorCode KSPGMRESMonitorKrylov(KSP,PetscInt,PetscReal,void *);
256: EXTERN PetscErrorCode KSPUnwindPreconditioner(KSP,Vec,Vec);
257: EXTERN PetscErrorCode KSPDefaultBuildSolution(KSP,Vec,Vec*);
258: EXTERN PetscErrorCode KSPDefaultBuildResidual(KSP,Vec,Vec,Vec *);
259: EXTERN PetscErrorCode KSPInitialResidual(KSP,Vec,Vec,Vec,Vec,Vec);
261: EXTERN PetscErrorCode KSPSetOperators(KSP,Mat,Mat,MatStructure);
262: EXTERN PetscErrorCode KSPGetOperators(KSP,Mat*,Mat*,MatStructure*);
263: EXTERN PetscErrorCode KSPGetOperatorsSet(KSP,PetscTruth*,PetscTruth*);
264: EXTERN PetscErrorCode KSPSetOptionsPrefix(KSP,const char[]);
265: EXTERN PetscErrorCode KSPAppendOptionsPrefix(KSP,const char[]);
266: EXTERN PetscErrorCode KSPGetOptionsPrefix(KSP,const char*[]);
268: EXTERN PetscErrorCode KSPSetDiagonalScale(KSP,PetscTruth);
269: EXTERN PetscErrorCode KSPGetDiagonalScale(KSP,PetscTruth*);
270: EXTERN PetscErrorCode KSPSetDiagonalScaleFix(KSP,PetscTruth);
271: EXTERN PetscErrorCode KSPGetDiagonalScaleFix(KSP,PetscTruth*);
273: EXTERN PetscErrorCode KSPView(KSP,PetscViewer);
275: EXTERN PetscErrorCode KSPLSQRSetStandardErrorVec(KSP,Vec);
276: EXTERN PetscErrorCode KSPLSQRGetStandardErrorVec(KSP,Vec*);
278: /*E
279: KSPNormType - Norm that is passed in the Krylov convergence
280: test routines.
282: Level: advanced
284: Notes: this must match finclude/petscksp.h
286: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPSetNormType(),
287: KSPSetConvergenceTest()
288: E*/
289: typedef enum {KSP_NORM_NO = 0,KSP_NORM_PRECONDITIONED = 1,KSP_NORM_UNPRECONDITIONED = 2,KSP_NORM_NATURAL = 3} KSPNormType;
291: /*MC
292: KSP_NORM_NO - Do not compute a norm during the Krylov process. This will
293: possibly save some computation but means the convergence test cannot
294: be based on a norm of a residual etc.
296: Level: advanced
298: Note: Some Krylov methods need to compute a residual norm and then this option is ignored
300: .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_PRECONDITIONED, KSP_NORM_UNPRECONDITIONED, KSP_NORM_NATURAL
301: M*/
303: /*MC
304: KSP_NORM_PRECONDITIONED - Compute the norm of the preconditioned residual and pass that to the
305: convergence test routine.
307: Level: advanced
309: .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_NO, KSP_NORM_UNPRECONDITIONED, KSP_NORM_NATURAL, KSPSetConvergenceTest()
310: M*/
312: /*MC
313: KSP_NORM_UNPRECONDITIONED - Compute the norm of the true residual (b - A*x) and pass that to the
314: convergence test routine.
316: Level: advanced
318: .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_NO, KSP_NORM_PRECONDITIONED, KSP_NORM_NATURAL, KSPSetConvergenceTest()
319: M*/
321: /*MC
322: KSP_NORM_NATURAL - Compute the 'natural norm' of residual sqrt((b - A*x)*B*(b - A*x)) and pass that to the
323: convergence test routine.
325: Level: advanced
327: .seealso: KSPNormType, KSPSetNormType(), KSP_NORM_NO, KSP_NORM_PRECONDITIONED, KSP_NORM_UNPRECONDITIONED, KSPSetConvergenceTest()
328: M*/
330: EXTERN PetscErrorCode KSPSetNormType(KSP,KSPNormType);
331: EXTERN PetscErrorCode KSPGetNormType(KSP,KSPNormType*);
332: EXTERN PetscErrorCode KSPSetCheckNormIteration(KSP,PetscInt);
334: /*E
335: KSPConvergedReason - reason a Krylov method was said to
336: have converged or diverged
338: Level: beginner
340: Notes: this must match finclude/petscksp.h
342: Developer note: The string versions of these are in
343: src/ksp/ksp/interface/itfunc.c called convergedreasons.
344: If these enums are changed you must change those.
346: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPSetTolerances()
347: E*/
348: typedef enum {/* converged */
349: KSP_CONVERGED_RTOL = 2,
350: KSP_CONVERGED_ATOL = 3,
351: KSP_CONVERGED_ITS = 4,
352: KSP_CONVERGED_CG_NEG_CURVE = 5,
353: KSP_CONVERGED_CG_CONSTRAINED = 6,
354: KSP_CONVERGED_STEP_LENGTH = 7,
355: KSP_CONVERGED_HAPPY_BREAKDOWN = 8,
356: /* diverged */
357: KSP_DIVERGED_NULL = -2,
358: KSP_DIVERGED_ITS = -3,
359: KSP_DIVERGED_DTOL = -4,
360: KSP_DIVERGED_BREAKDOWN = -5,
361: KSP_DIVERGED_BREAKDOWN_BICG = -6,
362: KSP_DIVERGED_NONSYMMETRIC = -7,
363: KSP_DIVERGED_INDEFINITE_PC = -8,
364: KSP_DIVERGED_NAN = -9,
365: KSP_DIVERGED_INDEFINITE_MAT = -10,
366:
367: KSP_CONVERGED_ITERATING = 0} KSPConvergedReason;
370: /*MC
371: KSP_CONVERGED_RTOL - norm(r) <= rtol*norm(b)
373: Level: beginner
375: See KSPNormType and KSPSetNormType() for possible norms that may be used. By default
376: for left preconditioning it is the 2-norm of the preconditioned residual, and the
377: 2-norm of the residual for right preconditioning
379: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
381: M*/
383: /*MC
384: KSP_CONVERGED_ATOL - norm(r) <= atol
386: Level: beginner
388: See KSPNormType and KSPSetNormType() for possible norms that may be used. By default
389: for left preconditioning it is the 2-norm of the preconditioned residual, and the
390: 2-norm of the residual for right preconditioning
392: Level: beginner
394: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
396: M*/
398: /*MC
399: KSP_DIVERGED_DTOL - norm(r) >= dtol*norm(b)
401: Level: beginner
403: See KSPNormType and KSPSetNormType() for possible norms that may be used. By default
404: for left preconditioning it is the 2-norm of the preconditioned residual, and the
405: 2-norm of the residual for right preconditioning
407: Level: beginner
409: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
411: M*/
413: /*MC
414: KSP_DIVERGED_ITS - Ran out of iterations before any convergence criteria was
415: reached
417: Level: beginner
419: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
421: M*/
423: /*MC
424: KSP_CONVERGED_ITS - Used by the KSPPREONLY solver after the single iteration of
425: the preconditioner is applied. Also used when the KSPSkipConverged() convergence
426: test rutine is set in KSP.
429: Level: beginner
432: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
434: M*/
436: /*MC
437: KSP_DIVERGED_BREAKDOWN - A breakdown in the Krylov method was detected so the
438: method could not continue to enlarge the Krylov space.
440: Level: beginner
442: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
444: M*/
446: /*MC
447: KSP_DIVERGED_BREAKDOWN_BICG - A breakdown in the KSPBICG method was detected so the
448: method could not continue to enlarge the Krylov space.
451: Level: beginner
454: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
456: M*/
458: /*MC
459: KSP_DIVERGED_NONSYMMETRIC - It appears the operator or preconditioner is not
460: symmetric and this Krylov method (KSPCG, KSPMINRES, KSPCR) requires symmetry
462: Level: beginner
464: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
466: M*/
468: /*MC
469: KSP_DIVERGED_INDEFINITE_PC - It appears the preconditioner is indefinite (has both
470: positive and negative eigenvalues) and this Krylov method (KSPCG) requires it to
471: be positive definite
473: Level: beginner
475: Notes: This can happen with the PCICC preconditioner, use -pc_factor_shift_positive_definite to force
476: the PCICC preconditioner to generate a positive definite preconditioner
478: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
480: M*/
482: /*MC
483: KSP_CONVERGED_ITERATING - This flag is returned if you call KSPGetConvergedReason()
484: while the KSPSolve() is still running.
486: Level: beginner
488: .seealso: KSPSolve(), KSPGetConvergedReason(), KSPConvergedReason, KSPSetTolerances()
490: M*/
492: EXTERN PetscErrorCode KSPSetConvergenceTest(KSP,PetscErrorCode (*)(KSP,PetscInt,PetscReal,KSPConvergedReason*,void*),void *,PetscErrorCode (*)(void*));
493: EXTERN PetscErrorCode KSPGetConvergenceContext(KSP,void **);
494: EXTERN PetscErrorCode KSPDefaultConverged(KSP,PetscInt,PetscReal,KSPConvergedReason*,void *);
495: EXTERN PetscErrorCode KSPDefaultConvergedDestroy(void *);
496: EXTERN PetscErrorCode KSPDefaultConvergedCreate(void **);
497: EXTERN PetscErrorCode KSPDefaultConvergedSetUIRNorm(KSP);
498: EXTERN PetscErrorCode KSPDefaultConvergedSetUMIRNorm(KSP);
499: EXTERN PetscErrorCode KSPSkipConverged(KSP,PetscInt,PetscReal,KSPConvergedReason*,void *);
500: EXTERN PetscErrorCode KSPGetConvergedReason(KSP,KSPConvergedReason *);
502: EXTERN PetscErrorCode KSPComputeExplicitOperator(KSP,Mat *);
504: /*E
505: KSPCGType - Determines what type of CG to use
507: Level: beginner
509: .seealso: KSPCGSetType()
510: E*/
511: typedef enum {KSP_CG_SYMMETRIC=0,KSP_CG_HERMITIAN=1} KSPCGType;
514: EXTERN PetscErrorCode KSPCGSetType(KSP,KSPCGType);
516: EXTERN PetscErrorCode KSPNASHSetRadius(KSP,PetscReal);
517: EXTERN PetscErrorCode KSPNASHGetNormD(KSP,PetscReal *);
518: EXTERN PetscErrorCode KSPNASHGetObjFcn(KSP,PetscReal *);
520: EXTERN PetscErrorCode KSPSTCGSetRadius(KSP,PetscReal);
521: EXTERN PetscErrorCode KSPSTCGGetNormD(KSP,PetscReal *);
522: EXTERN PetscErrorCode KSPSTCGGetObjFcn(KSP,PetscReal *);
524: EXTERN PetscErrorCode KSPGLTRSetRadius(KSP,PetscReal);
525: EXTERN PetscErrorCode KSPGLTRGetNormD(KSP,PetscReal *);
526: EXTERN PetscErrorCode KSPGLTRGetObjFcn(KSP,PetscReal *);
527: EXTERN PetscErrorCode KSPGLTRGetMinEig(KSP,PetscReal *);
528: EXTERN PetscErrorCode KSPGLTRGetLambda(KSP,PetscReal *);
530: EXTERN PetscErrorCode PCPreSolve(PC,KSP);
531: EXTERN PetscErrorCode PCPostSolve(PC,KSP);
533: EXTERN PetscErrorCode KSPMonitorLGCreate(const char[],const char[],int,int,int,int,PetscDrawLG*);
534: EXTERN PetscErrorCode KSPMonitorLG(KSP,PetscInt,PetscReal,void*);
535: EXTERN PetscErrorCode KSPMonitorLGDestroy(PetscDrawLG);
536: EXTERN PetscErrorCode KSPMonitorLGTrueResidualNormCreate(MPI_Comm,const char[],const char[],int,int,int,int,PetscDrawLG*);
537: EXTERN PetscErrorCode KSPMonitorLGTrueResidualNorm(KSP,PetscInt,PetscReal,void*);
538: EXTERN PetscErrorCode KSPMonitorLGTrueResidualNormDestroy(PetscDrawLG);
540: EXTERN PetscErrorCode PCShellSetPreSolve(PC,PetscErrorCode (*)(void*,KSP,Vec,Vec));
541: EXTERN PetscErrorCode PCShellSetPostSolve(PC,PetscErrorCode (*)(void*,KSP,Vec,Vec));
543: EXTERN PetscErrorCode MatMFFDKSPMonitor(KSP,PetscInt,PetscReal,void *);
546: #endif