Actual source code: epsimpl.h
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2007, Universidad Politecnica de Valencia, Spain
6: This file is part of SLEPc. See the README file for conditions of use
7: and additional information.
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
11: #ifndef _EPSIMPL
12: #define _EPSIMPL
14: #include slepceps.h
19: typedef struct _EPSOps *EPSOps;
21: struct _EPSOps {
22: int (*solve)(EPS); /* one-sided solver */
23: int (*solvets)(EPS); /* two-sided solver */
24: int (*setup)(EPS);
25: int (*setfromoptions)(EPS);
26: int (*publishoptions)(EPS);
27: int (*destroy)(EPS);
28: int (*view)(EPS,PetscViewer);
29: int (*backtransform)(EPS);
30: int (*computevectors)(EPS);
31: };
33: /*
34: Maximum number of monitors you can run with a single EPS
35: */
36: #define MAXEPSMONITORS 5
38: /*
39: Defines the EPS data structure.
40: */
41: struct _p_EPS {
42: PETSCHEADER(struct _EPSOps);
43: /*------------------------- User parameters --------------------------*/
44: int max_it, /* maximum number of iterations */
45: nev, /* number of eigenvalues to compute */
46: ncv, /* number of basis vectors */
47: nv, /* number of available basis vectors (<= ncv) */
48: allocated_ncv, /* number of basis vectors allocated */
49: nds; /* number of basis vectors of deflation space */
50: PetscReal tol; /* tolerance */
51: EPSWhich which; /* which part of the spectrum to be sought */
52: PetscTruth evecsavailable; /* computed eigenvectors */
53: EPSProblemType problem_type; /* which kind of problem to be solved */
54: EPSClass solverclass; /* whether the selected solver is one- or two-sided */
56: /*------------------------- Working data --------------------------*/
57: Vec vec_initial, /* initial vector */
58: vec_initial_left, /* left initial vector for two-sided solvers */
59: *V, /* set of basis vectors */
60: *AV, /* computed eigenvectors */
61: *W, /* set of left basis vectors */
62: *AW, /* computed left eigenvectors */
63: *DS, /* deflation space */
64: *DSV; /* deflation space and basis vectors*/
65: PetscScalar *eigr, *eigi, /* real and imaginary parts of eigenvalues */
66: *T, *Tl; /* projected matrices */
67: PetscReal *errest, /* error estimates */
68: *errest_left; /* left error estimates */
69: ST OP; /* spectral transformation object */
70: IP ip; /* innerproduct object */
71: void *data; /* placeholder for misc stuff associated
72: with a particular solver */
73: int nconv, /* number of converged eigenvalues */
74: its, /* number of iterations so far computed */
75: *perm; /* permutation for eigenvalue ordering */
77: /* ---------------- Default work-area and status vars -------------------- */
78: int nwork;
79: Vec *work;
81: int setupcalled;
82: PetscTruth isgeneralized,
83: ispositive,
84: ishermitian;
85: EPSConvergedReason reason;
87: int (*monitor[MAXEPSMONITORS])(EPS,int,int,PetscScalar*,PetscScalar*,PetscReal*,int,void*);
88: int (*monitordestroy[MAXEPSMONITORS])(void*);
89: void *monitorcontext[MAXEPSMONITORS];
90: int numbermonitors;
92: PetscTruth ds_ortho; /* if vectors in DS have to be orthonormalized */
93: };
95: #define EPSMonitor(eps,it,nconv,eigr,eigi,errest,nest) \
96: { int _ierr,_i,_im = eps->numbermonitors; \
97: for ( _i=0; _i<_im; _i++ ) {\
98: _ierr=(*eps->monitor[_i])(eps,it,nconv,eigr,eigi,errest,nest,eps->monitorcontext[_i]);\
99: CHKERRQ(_ierr); \
100: } \
101: }
103: EXTERN PetscErrorCode EPSRegisterAll(char *);
105: EXTERN PetscErrorCode EPSDestroy_Default(EPS);
106: EXTERN PetscErrorCode EPSDefaultGetWork(EPS,int);
107: EXTERN PetscErrorCode EPSDefaultFreeWork(EPS);
108: EXTERN PetscErrorCode EPSAllocateSolution(EPS);
109: EXTERN PetscErrorCode EPSFreeSolution(EPS);
110: EXTERN PetscErrorCode EPSAllocateSolutionContiguous(EPS);
111: EXTERN PetscErrorCode EPSFreeSolutionContiguous(EPS);
112: EXTERN PetscErrorCode EPSBackTransform_Default(EPS);
113: EXTERN PetscErrorCode EPSComputeVectors_Default(EPS);
114: EXTERN PetscErrorCode EPSComputeVectors_Hermitian(EPS);
115: EXTERN PetscErrorCode EPSComputeVectors_Schur(EPS);
117: /* Private functions of the solver implementations */
119: EXTERN PetscErrorCode EPSBasicArnoldi(EPS,PetscTruth,PetscScalar*,Vec*,int,int*,Vec,PetscReal*,PetscTruth*);
120: EXTERN PetscErrorCode EPSDelayedArnoldi(EPS,PetscScalar*,Vec*,int,int*,Vec,PetscReal*,PetscTruth*);
121: EXTERN PetscErrorCode EPSDelayedArnoldi1(EPS,PetscScalar*,Vec*,int,int*,Vec,PetscReal*,PetscTruth*);
122: EXTERN PetscErrorCode ArnoldiResiduals(PetscScalar*,int,PetscScalar*,PetscReal,int,int,PetscScalar*,PetscScalar*,PetscReal*,PetscScalar*);
124: #endif