Actual source code: epsimpl.h
2: #ifndef _EPSIMPL
3: #define _EPSIMPL
5: #include slepceps.h
10: typedef struct _EPSOps *EPSOps;
12: struct _EPSOps {
13: int (*solve)(EPS); /* one-sided solver */
14: int (*solvets)(EPS); /* two-sided solver */
15: int (*setup)(EPS);
16: int (*setfromoptions)(EPS);
17: int (*publishoptions)(EPS);
18: int (*destroy)(EPS);
19: int (*view)(EPS,PetscViewer);
20: int (*backtransform)(EPS);
21: int (*computevectors)(EPS);
22: };
24: /*
25: Maximum number of monitors you can run with a single EPS
26: */
27: #define MAXEPSMONITORS 5
29: /*
30: Defines the EPS data structure.
31: */
32: struct _p_EPS {
33: PETSCHEADER(struct _EPSOps);
34: /*------------------------- User parameters --------------------------*/
35: int max_it, /* maximum number of iterations */
36: nev, /* number of eigenvalues to compute */
37: ncv, /* number of basis vectors */
38: nv, /* number of available basis vectors (<= ncv) */
39: allocated_ncv, /* number of basis vectors allocated */
40: nds; /* number of basis vectors of deflation space */
41: PetscReal tol; /* tolerance */
42: EPSWhich which; /* which part of the spectrum to be sought */
43: PetscTruth evecsavailable; /* computed eigenvectors */
44: EPSProblemType problem_type; /* which kind of problem to be solved */
45: EPSClass solverclass; /* whether the selected solver is one- or two-sided */
47: /*------------------------- Working data --------------------------*/
48: Vec vec_initial, /* initial vector */
49: vec_initial_left, /* left initial vector for two-sided solvers */
50: *V, /* set of basis vectors */
51: *AV, /* computed eigenvectors */
52: *W, /* set of left basis vectors */
53: *AW, /* computed left eigenvectors */
54: *DS, /* deflation space */
55: *DSV; /* deflation space and basis vectors*/
56: PetscScalar *eigr, *eigi, /* real and imaginary parts of eigenvalues */
57: *T, *Tl; /* projected matrices */
58: PetscReal *errest, /* error estimates */
59: *errest_left; /* left error estimates */
60: ST OP; /* spectral transformation object */
61: void *data; /* placeholder for misc stuff associated
62: with a particular solver */
63: int nconv, /* number of converged eigenvalues */
64: its, /* number of iterations so far computed */
65: *perm; /* permutation for eigenvalue ordering */
67: /* ---------------- Default work-area and status vars -------------------- */
68: int nwork;
69: Vec *work;
71: int setupcalled;
72: PetscTruth isgeneralized,
73: ishermitian;
74: EPSConvergedReason reason;
76: int (*monitor[MAXEPSMONITORS])(EPS,int,int,PetscScalar*,PetscScalar*,PetscReal*,int,void*);
77: int (*monitordestroy[MAXEPSMONITORS])(void*);
78: void *monitorcontext[MAXEPSMONITORS];
79: int numbermonitors;
81: /* --------------- Orthogonalization --------------------- */
82: EPSOrthogonalizationType orthog_type; /* which orthogonalization to use */
83: EPSOrthogonalizationRefinementType orthog_ref; /* refinement method */
84: PetscReal orthog_eta;
85: PetscTruth ds_ortho; /* if vectors in DS have to be orthonormalized */
86: };
88: #define EPSMonitor(eps,it,nconv,eigr,eigi,errest,nest) \
89: { int _ierr,_i,_im = eps->numbermonitors; \
90: for ( _i=0; _i<_im; _i++ ) {\
91: _ierr=(*eps->monitor[_i])(eps,it,nconv,eigr,eigi,errest,nest,eps->monitorcontext[_i]);\
92: CHKERRQ(_ierr); \
93: } \
94: }
96: EXTERN PetscErrorCode EPSRegisterAll(char *);
97: EXTERN PetscErrorCode EPSRegister(const char*,const char*,const char*,int(*)(EPS));
98: #if defined(PETSC_USE_DYNAMIC_LIBRARIES)
99: #define EPSRegisterDynamic(a,b,c,d) EPSRegister(a,b,c,0)
100: #else
101: #define EPSRegisterDynamic(a,b,c,d) EPSRegister(a,b,c,d)
102: #endif
104: EXTERN PetscErrorCode EPSDestroy_Default(EPS);
105: EXTERN PetscErrorCode EPSDefaultGetWork(EPS,int);
106: EXTERN PetscErrorCode EPSDefaultFreeWork(EPS);
107: EXTERN PetscErrorCode EPSAllocateSolution(EPS);
108: EXTERN PetscErrorCode EPSFreeSolution(EPS);
109: EXTERN PetscErrorCode EPSAllocateSolutionContiguous(EPS);
110: EXTERN PetscErrorCode EPSFreeSolutionContiguous(EPS);
111: EXTERN PetscErrorCode EPSBackTransform_Default(EPS);
112: EXTERN PetscErrorCode EPSComputeVectors_Default(EPS);
113: EXTERN PetscErrorCode EPSComputeVectors_Schur(EPS);
115: /* Private functions of the solver implementations */
117: EXTERN PetscErrorCode EPSBasicArnoldi(EPS,PetscTruth,PetscScalar*,Vec*,int,int*,Vec,PetscReal*,PetscTruth*);
118: EXTERN PetscErrorCode EPSDelayedArnoldi(EPS,PetscScalar*,Vec*,int,int*,Vec,PetscReal*,PetscTruth*);
119: EXTERN PetscErrorCode EPSDelayedArnoldi1(EPS,PetscScalar*,Vec*,int,int*,Vec,PetscReal*,PetscTruth*);
120: EXTERN PetscErrorCode ArnoldiResiduals(PetscScalar*,int,PetscScalar*,PetscReal,int,int,PetscScalar*,PetscScalar*,PetscReal*,PetscScalar*);
122: #endif