Actual source code: epsimpl.h
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2013, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
8: SLEPc is free software: you can redistribute it and/or modify it under the
9: terms of version 3 of the GNU Lesser General Public License as published by
10: the Free Software Foundation.
12: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
13: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
14: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
15: more details.
17: You should have received a copy of the GNU Lesser General Public License
18: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
19: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
20: */
22: #if !defined(_EPSIMPL)
23: #define _EPSIMPL
25: #include <slepceps.h>
26: #include <slepc-private/slepcimpl.h>
28: PETSC_EXTERN PetscLogEvent EPS_SetUp,EPS_Solve;
30: typedef struct _EPSOps *EPSOps;
32: struct _EPSOps {
33: PetscErrorCode (*solve)(EPS);
34: PetscErrorCode (*setup)(EPS);
35: PetscErrorCode (*setfromoptions)(EPS);
36: PetscErrorCode (*publishoptions)(EPS);
37: PetscErrorCode (*destroy)(EPS);
38: PetscErrorCode (*reset)(EPS);
39: PetscErrorCode (*view)(EPS,PetscViewer);
40: PetscErrorCode (*backtransform)(EPS);
41: PetscErrorCode (*computevectors)(EPS);
42: };
44: /*
45: Maximum number of monitors you can run with a single EPS
46: */
47: #define MAXEPSMONITORS 5
49: /*
50: Defines the EPS data structure.
51: */
52: struct _p_EPS {
53: PETSCHEADER(struct _EPSOps);
54: /*------------------------- User parameters --------------------------*/
55: PetscInt max_it; /* maximum number of iterations */
56: PetscInt nev; /* number of eigenvalues to compute */
57: PetscInt ncv; /* number of basis vectors */
58: PetscInt mpd; /* maximum dimension of projected problem */
59: PetscInt nini,ninil; /* number of initial vectors (negative means not copied yet) */
60: PetscInt nds; /* number of basis vectors of deflation space */
61: PetscScalar target; /* target value */
62: PetscReal tol; /* tolerance */
63: EPSConv conv; /* convergence test */
64: EPSWhich which; /* which part of the spectrum to be sought */
65: PetscBool leftvecs; /* if left eigenvectors are requested */
66: PetscReal inta,intb; /* interval [a,b] for spectrum slicing */
67: EPSProblemType problem_type; /* which kind of problem to be solved */
68: EPSExtraction extraction; /* which kind of extraction to be applied */
69: EPSBalance balance; /* the balancing method */
70: PetscInt balance_its; /* number of iterations of the balancing method */
71: PetscReal balance_cutoff; /* cutoff value for balancing */
72: PetscReal nrma,nrmb; /* matrix norms */
73: PetscBool adaptive; /* whether matrix norms are adaptively improved */
74: PetscBool trueres; /* whether the true residual norm must be computed */
75: PetscBool trackall; /* whether all the residuals must be computed */
77: /*-------------- User-provided functions and contexts -----------------*/
78: PetscErrorCode (*comparison)(PetscScalar,PetscScalar,PetscScalar,PetscScalar,PetscInt*,void*);
79: PetscErrorCode (*converged)(EPS,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*);
80: PetscErrorCode (*arbitrary)(PetscScalar,PetscScalar,Vec,Vec,PetscScalar*,PetscScalar*,void*);
81: void *comparisonctx;
82: void *convergedctx;
83: void *arbitraryctx;
85: /*------------------------- Working data --------------------------*/
86: Vec D; /* diagonal matrix for balancing */
87: Vec *V; /* set of basis vectors and computed eigenvectors */
88: Vec *W; /* set of left basis vectors and computed left eigenvectors */
89: Vec *IS,*ISL; /* placeholder for references to user-provided initial space */
90: Vec *defl; /* deflation space */
91: PetscScalar *eigr,*eigi; /* real and imaginary parts of eigenvalues */
92: PetscReal *errest; /* error estimates */
93: PetscReal *errest_left; /* left error estimates */
94: PetscScalar *rr,*ri; /* values computed by user's arbitrary selection function */
95: ST st; /* spectral transformation object */
96: IP ip; /* innerproduct object */
97: DS ds; /* direct solver object */
98: void *data; /* placeholder for misc stuff associated
99: with a particular solver */
100: PetscInt nconv; /* number of converged eigenvalues */
101: PetscInt its; /* number of iterations so far computed */
102: PetscInt *perm; /* permutation for eigenvalue ordering */
103: PetscInt nv; /* size of current Schur decomposition */
104: PetscInt n,nloc; /* problem dimensions (global, local) */
105: PetscInt allocated_ncv; /* number of basis vectors allocated */
106: PetscBool evecsavailable; /* computed eigenvectors */
107: PetscRandom rand; /* random number generator */
108: Vec t; /* template vector */
110: /* ---------------- Default work-area and status vars -------------------- */
111: PetscInt nwork;
112: Vec *work;
114: PetscBool ds_ortho; /* if defl vectors have been stored & orthonormalized */
115: PetscInt setupcalled;
116: PetscBool isgeneralized;
117: PetscBool ispositive;
118: PetscBool ishermitian;
119: EPSConvergedReason reason;
121: PetscErrorCode (*monitor[MAXEPSMONITORS])(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
122: PetscErrorCode (*monitordestroy[MAXEPSMONITORS])(void**);
123: void *monitorcontext[MAXEPSMONITORS];
124: PetscInt numbermonitors;
125: };
127: PETSC_INTERN PetscErrorCode EPSReset_Default(EPS);
128: PETSC_INTERN PetscErrorCode EPSSetWhichEigenpairs_Default(EPS);
129: PETSC_INTERN PetscErrorCode EPSAllocateSolution(EPS);
130: PETSC_INTERN PetscErrorCode EPSFreeSolution(EPS);
131: PETSC_INTERN PetscErrorCode EPSBackTransform_Default(EPS);
132: PETSC_INTERN PetscErrorCode EPSComputeVectors_Default(EPS);
133: PETSC_INTERN PetscErrorCode EPSComputeVectors_Hermitian(EPS);
134: PETSC_INTERN PetscErrorCode EPSComputeVectors_Schur(EPS);
135: PETSC_INTERN PetscErrorCode EPSComputeResidualNorm_Private(EPS,PetscScalar,PetscScalar,Vec,Vec,PetscReal*);
136: PETSC_INTERN PetscErrorCode EPSComputeRelativeError_Private(EPS,PetscScalar,PetscScalar,Vec,Vec,PetscReal*);
137: PETSC_INTERN PetscErrorCode EPSComputeRitzVector(EPS,PetscScalar*,PetscScalar*,Vec*,PetscInt,Vec,Vec);
139: /* Private functions of the solver implementations */
141: PETSC_INTERN PetscErrorCode EPSBasicArnoldi(EPS,PetscBool,PetscScalar*,PetscInt,Vec*,PetscInt,PetscInt*,Vec,PetscReal*,PetscBool*);
142: PETSC_INTERN PetscErrorCode EPSDelayedArnoldi(EPS,PetscScalar*,PetscInt,Vec*,PetscInt,PetscInt*,Vec,PetscReal*,PetscBool*);
143: PETSC_INTERN PetscErrorCode EPSDelayedArnoldi1(EPS,PetscScalar*,PetscInt,Vec*,PetscInt,PetscInt*,Vec,PetscReal*,PetscBool*);
144: PETSC_INTERN PetscErrorCode EPSKrylovConvergence(EPS,PetscBool,PetscInt,PetscInt,Vec*,PetscInt,PetscReal,PetscReal,PetscInt*);
145: PETSC_INTERN PetscErrorCode EPSFullLanczos(EPS,PetscReal*,PetscReal*,Vec*,PetscInt,PetscInt*,Vec,PetscBool*);
146: PETSC_INTERN PetscErrorCode EPSBuildBalance_Krylov(EPS);
148: #endif