1 /* $Id: petscksp.h,v 1.96 2000/09/23 01:05:25 buschelm Exp balay $ */ 2 /* 3 Defines the interface functions for the Krylov subspace accelerators. 4 */ 5 #ifndef __PETSCKSP_H 6 #define __PETSCKSP_H 7 #include "petscpc.h" 8 9 #define KSP_COOKIE PETSC_COOKIE+8 10 11 typedef struct _p_KSP* KSP; 12 13 #define KSPRICHARDSON "richardson" 14 #define KSPCHEBYCHEV "chebychev" 15 #define KSPCG "cg" 16 #define KSPGMRES "gmres" 17 #define KSPTCQMR "tcqmr" 18 #define KSPBCGS "bcgs" 19 #define KSPCGS "cgs" 20 #define KSPTFQMR "tfqmr" 21 #define KSPCR "cr" 22 #define KSPLSQR "lsqr" 23 #define KSPPREONLY "preonly" 24 #define KSPQCG "qcg" 25 #define KSPBICG "bicg" 26 #define KSPFGMRES "fgmres" 27 #define KSPMINRES "minres" 28 #define KSPSYMMLQ "symmlq" 29 typedef char * KSPType; 30 31 EXTERN int KSPCreate(MPI_Comm,KSP *); 32 EXTERN int KSPSetType(KSP,KSPType); 33 EXTERN int KSPSetUp(KSP); 34 EXTERN int KSPSolve(KSP,int *); 35 EXTERN int KSPSolveTranspose(KSP,int *); 36 EXTERN int KSPDestroy(KSP); 37 38 extern FList KSPList; 39 EXTERN int KSPRegisterAll(char *); 40 EXTERN int KSPRegisterDestroy(void); 41 42 EXTERN int KSPRegister(char*,char*,char*,int(*)(KSP)); 43 #if defined(PETSC_USE_DYNAMIC_LIBRARIES) 44 #define KSPRegisterDynamic(a,b,c,d) KSPRegister(a,b,c,0) 45 #else 46 #define KSPRegisterDynamic(a,b,c,d) KSPRegister(a,b,c,d) 47 #endif 48 49 EXTERN int KSPGetType(KSP,KSPType *); 50 EXTERN int KSPSetPreconditionerSide(KSP,PCSide); 51 EXTERN int KSPGetPreconditionerSide(KSP,PCSide*); 52 EXTERN int KSPGetTolerances(KSP,double*,double*,double*,int*); 53 EXTERN int KSPSetTolerances(KSP,double,double,double,int); 54 EXTERN int KSPSetComputeResidual(KSP,PetscTruth); 55 EXTERN int KSPSetUsePreconditionedResidual(KSP); 56 EXTERN int KSPSetInitialGuessNonzero(KSP); 57 EXTERN int KSPGetInitialGuessNonzero(KSP,PetscTruth *); 58 EXTERN int KSPSetComputeEigenvalues(KSP); 59 EXTERN int KSPSetComputeSingularValues(KSP); 60 EXTERN int KSPSetRhs(KSP,Vec); 61 EXTERN int KSPGetRhs(KSP,Vec *); 62 EXTERN int KSPSetSolution(KSP,Vec); 63 EXTERN int KSPGetSolution(KSP,Vec *); 64 EXTERN int KSPGetResidualNorm(KSP,double*); 65 EXTERN int KSPGetIterationNumber(KSP,int*); 66 67 EXTERN int KSPSetPC(KSP,PC); 68 EXTERN int KSPGetPC(KSP,PC*); 69 70 EXTERN int KSPSetAvoidNorms(KSP); 71 72 EXTERN int KSPSetMonitor(KSP,int (*)(KSP,int,double,void*),void *,int (*)(void*)); 73 EXTERN int KSPClearMonitor(KSP); 74 EXTERN int KSPGetMonitorContext(KSP,void **); 75 EXTERN int KSPGetResidualHistory(KSP,double **,int *); 76 EXTERN int KSPSetResidualHistory(KSP,double *,int,PetscTruth); 77 78 79 EXTERN int KSPBuildSolution(KSP,Vec,Vec *); 80 EXTERN int KSPBuildResidual(KSP,Vec,Vec,Vec *); 81 82 EXTERN int KSPRichardsonSetScale(KSP,double); 83 EXTERN int KSPChebychevSetEigenvalues(KSP,double,double); 84 EXTERN int KSPComputeExtremeSingularValues(KSP,double*,double*); 85 EXTERN int KSPComputeEigenvalues(KSP,int,double*,double*,int *); 86 EXTERN int KSPComputeEigenvaluesExplicitly(KSP,int,double*,double*); 87 88 EXTERN int KSPGMRESSetRestart(KSP,int); 89 EXTERN int KSPGMRESSetHapTol(KSP,double); 90 EXTERN int KSPGMRESSetPreAllocateVectors(KSP); 91 EXTERN int KSPGMRESSetOrthogonalization(KSP,int (*)(KSP,int)); 92 EXTERN int KSPGMRESUnmodifiedGramSchmidtOrthogonalization(KSP,int); 93 EXTERN int KSPGMRESModifiedGramSchmidtOrthogonalization(KSP,int); 94 EXTERN int KSPGMRESIROrthogonalization(KSP,int); 95 96 EXTERN int KSPFGMRESModifyPCNoChange(KSP,int,int,double,void*); 97 EXTERN int KSPFGMRESModifyPCSLES(KSP,int,int,double,void*); 98 EXTERN int KSPFGMRESSetModifyPC(KSP,int (*)(KSP,int,int,double,void*),void*,int(*)(void*)); 99 100 EXTERN int KSPSetFromOptions(KSP); 101 EXTERN int KSPSetTypeFromOptions(KSP); 102 EXTERN int KSPAddOptionsChecker(int (*)(KSP)); 103 104 EXTERN int KSPSingularValueMonitor(KSP,int,double,void *); 105 EXTERN int KSPDefaultMonitor(KSP,int,double,void *); 106 EXTERN int KSPTrueMonitor(KSP,int,double,void *); 107 EXTERN int KSPDefaultSMonitor(KSP,int,double,void *); 108 EXTERN int KSPVecViewMonitor(KSP,int,double,void *); 109 EXTERN int KSPGMRESKrylovMonitor(KSP,int,double,void *); 110 111 112 EXTERN int KSPResidual(KSP,Vec,Vec,Vec,Vec,Vec,Vec); 113 EXTERN int KSPUnwindPreconditioner(KSP,Vec,Vec); 114 EXTERN int KSPDefaultBuildSolution(KSP,Vec,Vec*); 115 EXTERN int KSPDefaultBuildResidual(KSP,Vec,Vec,Vec *); 116 117 EXTERN int KSPSetOptionsPrefix(KSP,char*); 118 EXTERN int KSPAppendOptionsPrefix(KSP,char*); 119 EXTERN int KSPGetOptionsPrefix(KSP,char**); 120 121 EXTERN int KSPView(KSP,Viewer); 122 123 /* this table must match finclude/petscksp.h */ 124 typedef enum {/* converged */ 125 KSP_CONVERGED_RTOL = 2, 126 KSP_CONVERGED_ATOL = 3, 127 KSP_CONVERGED_ITS = 4, 128 KSP_CONVERGED_QCG_NEG_CURVE = 5, 129 KSP_CONVERGED_QCG_CONSTRAINED = 6, 130 KSP_CONVERGED_STEP_LENGTH = 7, 131 /* diverged */ 132 KSP_DIVERGED_ITS = -3, 133 KSP_DIVERGED_DTOL = -4, 134 KSP_DIVERGED_BREAKDOWN = -5, 135 KSP_DIVERGED_BREAKDOWN_BICG = -6, 136 KSP_DIVERGED_NONSYMMETRIC = -7, 137 KSP_DIVERGED_INDEFINITE_PC = -8, 138 139 KSP_CONVERGED_ITERATING = 0} KSPConvergedReason; 140 141 EXTERN int KSPSetConvergenceTest(KSP,int (*)(KSP,int,double,KSPConvergedReason*,void*),void *); 142 EXTERN int KSPGetConvergenceContext(KSP,void **); 143 EXTERN int KSPDefaultConverged(KSP,int,double,KSPConvergedReason*,void *); 144 EXTERN int KSPSkipConverged(KSP,int,double,KSPConvergedReason*,void *); 145 EXTERN int KSPGetConvergedReason(KSP,KSPConvergedReason *); 146 147 EXTERN int KSPComputeExplicitOperator(KSP,Mat *); 148 149 typedef enum {KSP_CG_SYMMETRIC=1,KSP_CG_HERMITIAN=2} KSPCGType; 150 EXTERN int KSPCGSetType(KSP,KSPCGType); 151 152 EXTERN int PCPreSolve(PC,KSP); 153 EXTERN int PCPostSolve(PC,KSP); 154 155 EXTERN int KSPLGMonitorCreate(char*,char*,int,int,int,int,DrawLG*); 156 EXTERN int KSPLGMonitor(KSP,int,double,void*); 157 EXTERN int KSPLGMonitorDestroy(DrawLG); 158 EXTERN int KSPLGTrueMonitorCreate(MPI_Comm,char*,char*,int,int,int,int,DrawLG*); 159 EXTERN int KSPLGTrueMonitor(KSP,int,double,void*); 160 EXTERN int KSPLGTrueMonitorDestroy(DrawLG); 161 162 #endif 163 164 165