1 /* $Id: ksp.h,v 1.61 1998/01/17 17:39:51 bsmith Exp curfman $ */ 2 /* 3 Defines the interface functions for the Krylov subspace accelerators. 4 */ 5 #ifndef __KSP_PACKAGE 6 #define __KSP_PACKAGE 7 #include "petsc.h" 8 #include "vec.h" 9 #include "mat.h" 10 #include "pc.h" 11 12 #define KSP_COOKIE PETSC_COOKIE+8 13 14 typedef struct _p_KSP* KSP; 15 16 typedef enum { KSPUNKNOWN = -1,KSPRICHARDSON, KSPCHEBYCHEV, KSPCG, KSPGMRES, KSPTCQMR, KSPBCGS, 17 KSPCGS, KSPTFQMR, KSPCR, KSPLSQR, KSPPREONLY, KSPQCG, KSPNEW, KSPTRLS} KSPType; 18 19 extern int KSPCreate(MPI_Comm,KSP *); 20 extern int KSPSetType(KSP,KSPType); 21 extern int KSPSetUp(KSP); 22 extern int KSPSolve(KSP,int *); 23 extern int KSPDestroy(KSP); 24 25 extern int KSPRegisterAll(); 26 extern int KSPRegisterDestroy(); 27 extern int KSPRegister_Private(KSPType,char *,char *,int (*)(KSP),KSPType*); 28 #if defined(USE_DYNAMIC_LIBRARIES) 29 #define KSPRegister(a,b,c,d,e) KSPRegister_Private(a,b,c,0,e) 30 #else 31 #define KSPRegister(a,b,c,d,e) KSPRegister_Private(a,b,c,d,e) 32 #endif 33 34 extern int KSPGetType(KSP, KSPType *,char **); 35 extern int KSPSetPreconditionerSide(KSP,PCSide); 36 extern int KSPGetPreconditionerSide(KSP,PCSide*); 37 extern int KSPGetTolerances(KSP,double*,double*,double*,int*); 38 extern int KSPSetTolerances(KSP,double,double,double,int); 39 extern int KSPSetComputeResidual(KSP,PetscTruth); 40 extern int KSPSetUsePreconditionedResidual(KSP); 41 extern int KSPSetInitialGuessNonzero(KSP); 42 extern int KSPSetComputeEigenvalues(KSP); 43 extern int KSPSetComputeSingularValues(KSP); 44 extern int KSPSetRhs(KSP,Vec); 45 extern int KSPGetRhs(KSP,Vec *); 46 extern int KSPSetSolution(KSP,Vec); 47 extern int KSPGetSolution(KSP,Vec *); 48 extern int KSPGetResidualNorm(KSP,double*); 49 50 extern int KSPSetPC(KSP,PC); 51 extern int KSPGetPC(KSP,PC*); 52 53 extern int KSPSetMonitor(KSP,int (*)(KSP,int,double, void*), void *); 54 extern int KSPGetMonitorContext(KSP,void **); 55 extern int KSPSetResidualHistory(KSP, double *,int); 56 57 extern int KSPSetConvergenceTest(KSP,int (*)(KSP,int,double, void*), void *); 58 extern int KSPGetConvergenceContext(KSP,void **); 59 60 extern int KSPBuildSolution(KSP, Vec,Vec *); 61 extern int KSPBuildResidual(KSP, Vec, Vec,Vec *); 62 63 extern int KSPRichardsonSetScale(KSP , double); 64 extern int KSPChebychevSetEigenvalues(KSP , double, double); 65 extern int KSPComputeExtremeSingularValues(KSP, double*,double*); 66 extern int KSPComputeEigenvalues(KSP,int,double*,double*); 67 extern int KSPComputeEigenvaluesExplicitly(KSP,int,double*,double*); 68 69 extern int KSPGMRESSetRestart(KSP, int); 70 extern int KSPGMRESSetPreAllocateVectors(KSP); 71 extern int KSPGMRESSetOrthogonalization(KSP,int (*)(KSP,int)); 72 extern int KSPGMRESUnmodifiedGramSchmidtOrthogonalization(KSP,int); 73 extern int KSPGMRESModifiedGramSchmidtOrthogonalization(KSP,int); 74 extern int KSPGMRESIROrthogonalization(KSP,int); 75 extern int KSPGMRESDGKSOrthogonalization(KSP,int); 76 77 extern int KSPSetFromOptions(KSP); 78 extern int KSPAddOptionsChecker(int (*)(KSP)); 79 80 extern int KSPSingularValueMonitor(KSP,int,double, void * ); 81 extern int KSPDefaultMonitor(KSP,int,double, void *); 82 extern int KSPTrueMonitor(KSP,int,double, void *); 83 extern int KSPDefaultSMonitor(KSP,int,double, void *); 84 85 extern int KSPDefaultConverged(KSP,int,double, void *); 86 87 extern int KSPResidual(KSP,Vec,Vec,Vec,Vec,Vec,Vec); 88 extern int KSPUnwindPreconditioner(KSP,Vec,Vec); 89 extern int KSPDefaultBuildSolution(KSP,Vec,Vec*); 90 extern int KSPDefaultBuildResidual(KSP,Vec,Vec,Vec *); 91 92 extern int KSPPrintHelp(KSP); 93 94 extern int KSPSetOptionsPrefix(KSP,char*); 95 extern int KSPAppendOptionsPrefix(KSP,char*); 96 extern int KSPGetOptionsPrefix(KSP,char**); 97 98 extern int KSPView(KSP,Viewer); 99 100 extern int KSPComputeExplicitOperator(KSP,Mat *); 101 102 typedef enum {KSP_CG_SYMMETRIC=1, KSP_CG_HERMITIAN=2} KSPCGType; 103 extern int KSPCGSetType(KSP,KSPCGType); 104 105 extern int PCPreSolve(PC,KSP); 106 extern int PCPostSolve(PC,KSP); 107 108 extern int KSPLGMonitorCreate(char*,char*,int,int,int,int,DrawLG*); 109 extern int KSPLGMonitor(KSP,int,double,void*); 110 extern int KSPLGMonitorDestroy(DrawLG); 111 extern int KSPLGTrueMonitorCreate(MPI_Comm,char*,char*,int,int,int,int,DrawLG*); 112 extern int KSPLGTrueMonitor(KSP,int,double,void*); 113 extern int KSPLGTrueMonitorDestroy(DrawLG); 114 115 #endif 116 117 118