Actual source code: epsimpl.h
slepc-3.13.4 2020-09-02
1: /*
2: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3: SLEPc - Scalable Library for Eigenvalue Problem Computations
4: Copyright (c) 2002-2020, Universitat Politecnica de Valencia, Spain
6: This file is part of SLEPc.
7: SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: */
11: #if !defined(SLEPCEPSIMPL_H)
12: #define SLEPCEPSIMPL_H
14: #include <slepceps.h>
15: #include <slepc/private/slepcimpl.h>
17: SLEPC_EXTERN PetscBool EPSRegisterAllCalled;
18: SLEPC_EXTERN PetscErrorCode EPSRegisterAll(void);
19: SLEPC_EXTERN PetscLogEvent EPS_SetUp,EPS_Solve;
21: typedef struct _EPSOps *EPSOps;
23: struct _EPSOps {
24: PetscErrorCode (*solve)(EPS);
25: PetscErrorCode (*setup)(EPS);
26: PetscErrorCode (*setfromoptions)(PetscOptionItems*,EPS);
27: PetscErrorCode (*publishoptions)(EPS);
28: PetscErrorCode (*destroy)(EPS);
29: PetscErrorCode (*reset)(EPS);
30: PetscErrorCode (*view)(EPS,PetscViewer);
31: PetscErrorCode (*backtransform)(EPS);
32: PetscErrorCode (*computevectors)(EPS);
33: PetscErrorCode (*setdefaultst)(EPS);
34: };
36: /*
37: Maximum number of monitors you can run with a single EPS
38: */
39: #define MAXEPSMONITORS 5
41: /*
42: The solution process goes through several states
43: */
44: typedef enum { EPS_STATE_INITIAL,
45: EPS_STATE_SETUP,
46: EPS_STATE_SOLVED,
47: EPS_STATE_EIGENVECTORS } EPSStateType;
49: /*
50: To classify the different solvers into categories
51: */
52: typedef enum { EPS_CATEGORY_KRYLOV, /* Krylov solver: relies on STApply and STBackTransform (same as OTHER) */
53: EPS_CATEGORY_PRECOND, /* Preconditioned solver: uses ST only to manage preconditioner */
54: EPS_CATEGORY_CONTOUR, /* Contour integral: ST used to solve linear systems at integration points */
55: EPS_CATEGORY_OTHER } EPSSolverType;
57: /*
58: Defines the EPS data structure
59: */
60: struct _p_EPS {
61: PETSCHEADER(struct _EPSOps);
62: /*------------------------- User parameters ---------------------------*/
63: PetscInt max_it; /* maximum number of iterations */
64: PetscInt nev; /* number of eigenvalues to compute */
65: PetscInt ncv; /* number of basis vectors */
66: PetscInt mpd; /* maximum dimension of projected problem */
67: PetscInt nini,ninil; /* number of initial vectors (negative means not copied yet) */
68: PetscInt nds; /* number of basis vectors of deflation space */
69: PetscScalar target; /* target value */
70: PetscReal tol; /* tolerance */
71: EPSConv conv; /* convergence test */
72: EPSStop stop; /* stopping test */
73: EPSWhich which; /* which part of the spectrum to be sought */
74: PetscReal inta,intb; /* interval [a,b] for spectrum slicing */
75: EPSProblemType problem_type; /* which kind of problem to be solved */
76: EPSExtraction extraction; /* which kind of extraction to be applied */
77: EPSBalance balance; /* the balancing method */
78: PetscInt balance_its; /* number of iterations of the balancing method */
79: PetscReal balance_cutoff; /* cutoff value for balancing */
80: PetscBool trueres; /* whether the true residual norm must be computed */
81: PetscBool trackall; /* whether all the residuals must be computed */
82: PetscBool purify; /* whether eigenvectors need to be purified */
83: PetscBool twosided; /* whether to compute left eigenvectors (two-sided solver) */
85: /*-------------- User-provided functions and contexts -----------------*/
86: PetscErrorCode (*converged)(EPS,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*);
87: PetscErrorCode (*convergeduser)(EPS,PetscScalar,PetscScalar,PetscReal,PetscReal*,void*);
88: PetscErrorCode (*convergeddestroy)(void*);
89: PetscErrorCode (*stopping)(EPS,PetscInt,PetscInt,PetscInt,PetscInt,EPSConvergedReason*,void*);
90: PetscErrorCode (*stoppinguser)(EPS,PetscInt,PetscInt,PetscInt,PetscInt,EPSConvergedReason*,void*);
91: PetscErrorCode (*stoppingdestroy)(void*);
92: PetscErrorCode (*arbitrary)(PetscScalar,PetscScalar,Vec,Vec,PetscScalar*,PetscScalar*,void*);
93: void *convergedctx;
94: void *stoppingctx;
95: void *arbitraryctx;
96: PetscErrorCode (*monitor[MAXEPSMONITORS])(EPS,PetscInt,PetscInt,PetscScalar*,PetscScalar*,PetscReal*,PetscInt,void*);
97: PetscErrorCode (*monitordestroy[MAXEPSMONITORS])(void**);
98: void *monitorcontext[MAXEPSMONITORS];
99: PetscInt numbermonitors;
101: /*----------------- Child objects and working data -------------------*/
102: ST st; /* spectral transformation object */
103: DS ds; /* direct solver object */
104: DS dsts; /* auxiliary direct solver object used in two-sided case */
105: BV V; /* set of basis vectors and computed eigenvectors */
106: BV W; /* left basis vectors (if left eigenvectors requested) */
107: RG rg; /* optional region for filtering */
108: SlepcSC sc; /* sorting criterion data */
109: Vec D; /* diagonal matrix for balancing */
110: Vec *IS,*ISL; /* references to user-provided initial spaces */
111: Vec *defl; /* references to user-provided deflation space */
112: PetscScalar *eigr,*eigi; /* real and imaginary parts of eigenvalues */
113: PetscReal *errest; /* error estimates */
114: PetscScalar *rr,*ri; /* values computed by user's arbitrary selection function */
115: PetscInt *perm; /* permutation for eigenvalue ordering */
116: PetscInt nwork; /* number of work vectors */
117: Vec *work; /* work vectors */
118: void *data; /* placeholder for solver-specific stuff */
120: /* ----------------------- Status variables --------------------------*/
121: EPSStateType state; /* initial -> setup -> solved -> eigenvectors */
122: EPSSolverType categ; /* solver category */
123: PetscInt nconv; /* number of converged eigenvalues */
124: PetscInt its; /* number of iterations so far computed */
125: PetscInt n,nloc; /* problem dimensions (global, local) */
126: PetscReal nrma,nrmb; /* computed matrix norms */
127: PetscBool useds; /* whether the solver uses the DS object or not */
128: PetscBool hasts; /* whether the solver has two-sided variant */
129: PetscBool isgeneralized;
130: PetscBool ispositive;
131: PetscBool ishermitian;
132: EPSConvergedReason reason;
133: };
135: /*
136: Macros to test valid EPS arguments
137: */
138: #if !defined(PETSC_USE_DEBUG)
140: #define EPSCheckSolved(h,arg) do {} while (0)
142: #else
144: #define EPSCheckSolved(h,arg) \
145: do { \
146: if ((h)->state<EPS_STATE_SOLVED) SETERRQ1(PetscObjectComm((PetscObject)(h)),PETSC_ERR_ARG_WRONGSTATE,"Must call EPSSolve() first: Parameter #%d",arg); \
147: } while (0)
149: #endif
151: /*
152: EPS_SetInnerProduct - set B matrix for inner product if appropriate.
153: */
154: PETSC_STATIC_INLINE PetscErrorCode EPS_SetInnerProduct(EPS eps)
155: {
157: Mat B;
160: if (!eps->V) { EPSGetBV(eps,&eps->V); }
161: if (eps->ispositive || (eps->isgeneralized && eps->ishermitian)) {
162: STGetBilinearForm(eps->st,&B);
163: BVSetMatrix(eps->V,B,PetscNot(eps->ispositive));
164: MatDestroy(&B);
165: } else {
166: BVSetMatrix(eps->V,NULL,PETSC_FALSE);
167: }
168: return(0);
169: }
171: /*
172: EPS_Purify - purify the first k vectors in the V basis
173: */
174: PETSC_STATIC_INLINE PetscErrorCode EPS_Purify(EPS eps,PetscInt k)
175: {
177: PetscInt i;
178: Vec v,z;
181: BVCreateVec(eps->V,&v);
182: for (i=0;i<k;i++) {
183: BVCopyVec(eps->V,i,v);
184: BVGetColumn(eps->V,i,&z);
185: STApply(eps->st,v,z);
186: BVRestoreColumn(eps->V,i,&z);
187: }
188: VecDestroy(&v);
189: return(0);
190: }
192: SLEPC_INTERN PetscErrorCode EPSSetWhichEigenpairs_Default(EPS);
193: SLEPC_INTERN PetscErrorCode EPSSetDimensions_Default(EPS,PetscInt,PetscInt*,PetscInt*);
194: SLEPC_INTERN PetscErrorCode EPSBackTransform_Default(EPS);
195: SLEPC_INTERN PetscErrorCode EPSComputeVectors(EPS);
196: SLEPC_INTERN PetscErrorCode EPSComputeVectors_Hermitian(EPS);
197: SLEPC_INTERN PetscErrorCode EPSComputeVectors_Schur(EPS);
198: SLEPC_INTERN PetscErrorCode EPSComputeVectors_Indefinite(EPS);
199: SLEPC_INTERN PetscErrorCode EPSComputeVectors_Twosided(EPS);
200: SLEPC_INTERN PetscErrorCode EPSComputeVectors_Slice(EPS);
201: SLEPC_INTERN PetscErrorCode EPSComputeResidualNorm_Private(EPS,PetscBool,PetscScalar,PetscScalar,Vec,Vec,Vec*,PetscReal*);
202: SLEPC_INTERN PetscErrorCode EPSComputeRitzVector(EPS,PetscScalar*,PetscScalar*,BV,Vec,Vec);
203: SLEPC_INTERN PetscErrorCode EPSGetStartVector(EPS,PetscInt,PetscBool*);
204: SLEPC_INTERN PetscErrorCode EPSGetLeftStartVector(EPS,PetscInt,PetscBool*);
206: /* Private functions of the solver implementations */
208: SLEPC_INTERN PetscErrorCode EPSDelayedArnoldi(EPS,PetscScalar*,PetscInt,PetscInt,PetscInt*,PetscReal*,PetscBool*);
209: SLEPC_INTERN PetscErrorCode EPSDelayedArnoldi1(EPS,PetscScalar*,PetscInt,PetscInt,PetscInt*,PetscReal*,PetscBool*);
210: SLEPC_INTERN PetscErrorCode EPSKrylovConvergence(EPS,PetscBool,PetscInt,PetscInt,PetscReal,PetscReal,PetscReal,PetscInt*);
211: SLEPC_INTERN PetscErrorCode EPSPseudoLanczos(EPS,PetscReal*,PetscReal*,PetscReal*,PetscInt,PetscInt*,PetscBool*,PetscBool*,PetscReal*,Vec);
212: SLEPC_INTERN PetscErrorCode EPSBuildBalance_Krylov(EPS);
213: SLEPC_INTERN PetscErrorCode EPSSetDefaultST(EPS);
214: SLEPC_INTERN PetscErrorCode EPSSetDefaultST_Precond(EPS);
215: SLEPC_INTERN PetscErrorCode EPSSetDefaultST_GMRES(EPS);
217: #endif