Line data Source code
1 : /*
2 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
3 : SLEPc - Scalable Library for Eigenvalue Problem Computations
4 : Copyright (c) 2002-, Universitat Politecnica de Valencia, Spain
5 :
6 : This file is part of SLEPc.
7 : SLEPc is distributed under a 2-clause BSD license (see LICENSE).
8 : - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9 : */
10 : /*
11 : Method: General Davidson Method (includes GD and JD)
12 :
13 : References:
14 : - Ernest R. Davidson. Super-matrix methods. Computer Physics Communications,
15 : 53:49-60, May 1989.
16 : */
17 :
18 : #pragma once
19 :
20 : #include <slepc/private/epsimpl.h>
21 : #include <slepc/private/vecimplslepc.h>
22 :
23 : typedef PetscInt MatType_t;
24 : #define DVD_MAT_HERMITIAN (1<<1)
25 : #define DVD_MAT_NEG_DEF (1<<2)
26 : #define DVD_MAT_POS_DEF (1<<3)
27 : #define DVD_MAT_SINGULAR (1<<4)
28 : #define DVD_MAT_COMPLEX (1<<5)
29 : #define DVD_MAT_IMPLICIT (1<<6)
30 : #define DVD_MAT_IDENTITY (1<<7)
31 : #define DVD_MAT_DIAG (1<<8)
32 : #define DVD_MAT_TRIANG (1<<9)
33 : #define DVD_MAT_UTRIANG (1<<9)
34 : #define DVD_MAT_LTRIANG (1<<10)
35 : #define DVD_MAT_UNITARY (1<<11)
36 :
37 : typedef PetscInt EPType_t;
38 : #define DVD_EP_STD (1<<1)
39 : #define DVD_EP_HERMITIAN (1<<2)
40 : #define DVD_EP_INDEFINITE (1<<3)
41 :
42 : #define DVD_IS(T,P) ((T) & (P))
43 : #define DVD_ISNOT(T,P) (((T) & (P)) ^ (P))
44 :
45 : struct _dvdDashboard;
46 : typedef PetscErrorCode (*dvdCallback)(struct _dvdDashboard*);
47 : typedef struct _dvdFunctionList {
48 : dvdCallback f;
49 : struct _dvdFunctionList *next;
50 : } dvdFunctionList;
51 :
52 : typedef enum {
53 : DVD_HARM_NONE,
54 : DVD_HARM_RR,
55 : DVD_HARM_RRR,
56 : DVD_HARM_REIGS,
57 : DVD_HARM_LEIGS
58 : } HarmType_t;
59 :
60 : typedef enum {
61 : DVD_INITV_CLASSIC,
62 : DVD_INITV_KRYLOV
63 : } InitType_t;
64 :
65 : /*
66 : Dashboard struct: contains the methods that will be employed and the tuning
67 : options.
68 : */
69 :
70 : typedef struct _dvdDashboard {
71 : /**** Function steps ****/
72 : /* Initialize V */
73 : PetscErrorCode (*initV)(struct _dvdDashboard*);
74 : void *initV_data;
75 :
76 : /* Find the approximate eigenpairs from V */
77 : PetscErrorCode (*calcPairs)(struct _dvdDashboard*);
78 : void *calcPairs_data;
79 :
80 : /* Eigenpair test for convergence */
81 : PetscBool (*testConv)(struct _dvdDashboard*,PetscScalar,PetscScalar,PetscReal,PetscReal*);
82 : void *testConv_data;
83 :
84 : /* Improve the selected eigenpairs */
85 : PetscErrorCode (*improveX)(struct _dvdDashboard*,PetscInt,PetscInt,PetscInt*);
86 : void *improveX_data;
87 :
88 : /* Check for restarting */
89 : PetscErrorCode (*isRestarting)(struct _dvdDashboard*,PetscBool*);
90 : void *isRestarting_data;
91 :
92 : /* Perform restarting */
93 : PetscErrorCode (*restartV)(struct _dvdDashboard*);
94 : void *restartV_data;
95 :
96 : /* Update V */
97 : PetscErrorCode (*updateV)(struct _dvdDashboard*);
98 : void *updateV_data;
99 :
100 : /**** Problem specification ****/
101 : Mat A,B; /* problem matrices */
102 : MatType_t sA,sB; /* matrix specifications */
103 : EPType_t sEP; /* problem specifications */
104 : PetscInt nev; /* number of eigenpairs */
105 : EPSWhich which; /* spectrum selection */
106 : PetscBool withTarget; /* if there is a target */
107 : PetscScalar target[2]; /* target value */
108 : PetscReal tol; /* tolerance */
109 : PetscBool correctXnorm; /* if true, norm of X are computed */
110 :
111 : /**** Subspaces specification ****/
112 : PetscInt nconv; /* number of converged eigenpairs */
113 : PetscInt npreconv; /* number of pairs ready to converge */
114 :
115 : BV W; /* left basis for harmonic case */
116 : BV AX; /* A*V */
117 : BV BX; /* B*V */
118 : PetscInt size_D; /* active vectors */
119 : PetscInt max_size_proj; /* max size projected problem */
120 : PetscInt max_size_P; /* max unconverged vectors in the projector */
121 : PetscInt bs; /* max vectors that expands the subspace every iteration */
122 : EPS eps; /* connection to SLEPc */
123 :
124 : /**** Auxiliary space ****/
125 : VecPool auxV; /* auxiliary vectors */
126 : BV auxBV; /* auxiliary vectors */
127 :
128 : /**** Eigenvalues and errors ****/
129 : PetscScalar *ceigr,*ceigi; /* converged eigenvalues */
130 : PetscScalar *eigr,*eigi; /* current eigenvalues */
131 : PetscReal *nR; /* residual norm */
132 : PetscReal *real_nR; /* original nR */
133 : PetscReal *nX; /* X norm */
134 : PetscReal *real_nX; /* original nX */
135 : PetscReal *errest; /* relative error eigenpairs */
136 : PetscReal *nBds; /* B-norms of projected problem */
137 :
138 : /**** Shared function and variables ****/
139 : PetscErrorCode (*e_Vchanged)(struct _dvdDashboard*,PetscInt,PetscInt,PetscInt,PetscInt);
140 : void *e_Vchanged_data;
141 : PetscErrorCode (*calcpairs_residual)(struct _dvdDashboard*,PetscInt,PetscInt);
142 : PetscErrorCode (*calcpairs_selectPairs)(struct _dvdDashboard*,PetscInt);
143 : void *calcpairs_residual_data;
144 : PetscErrorCode (*improvex_precond)(struct _dvdDashboard*,PetscInt,Vec,Vec);
145 : void *improvex_precond_data;
146 : PetscErrorCode (*improvex_jd_proj_uv)(struct _dvdDashboard*,PetscInt,PetscInt,Vec*,Vec*,Vec*,PetscScalar*,PetscScalar*,PetscScalar*,PetscScalar*,PetscInt);
147 : PetscErrorCode (*improvex_jd_lit)(struct _dvdDashboard*,PetscInt,PetscScalar*,PetscScalar*,PetscInt*,PetscReal*);
148 : PetscErrorCode (*calcpairs_W)(struct _dvdDashboard*);
149 : void *calcpairs_W_data;
150 : PetscErrorCode (*calcpairs_proj_trans)(struct _dvdDashboard*);
151 : PetscErrorCode (*calcpairs_eigs_trans)(struct _dvdDashboard*);
152 : PetscErrorCode (*calcpairs_eig_backtrans)(struct _dvdDashboard*,PetscScalar,PetscScalar,PetscScalar*,PetscScalar*);
153 : PetscErrorCode (*calcpairs_proj_res)(struct _dvdDashboard*,PetscInt,PetscInt,Vec*);
154 : PetscErrorCode (*preTestConv)(struct _dvdDashboard*,PetscInt,PetscInt,PetscInt,PetscInt*);
155 : PetscErrorCode (*e_newIteration)(struct _dvdDashboard*);
156 : void *e_newIteration_data;
157 :
158 : dvdFunctionList *startList; /* starting list */
159 : dvdFunctionList *endList; /* ending list */
160 : dvdFunctionList *destroyList;/* destructor list */
161 :
162 : Mat H,G; /* projected problem matrices */
163 : Mat auxM; /* auxiliary dense matrix */
164 : PetscInt size_MT; /* rows in MT */
165 :
166 : PetscInt V_tra_s;
167 : PetscInt V_tra_e; /* cX <- [cX V*MT(0:V_tra_s-1)], V <- V*MT(V_tra_s:V_tra_e) */
168 : PetscInt V_new_s;
169 : PetscInt V_new_e; /* added to V the columns V_new_s:V_new_e */
170 : PetscBool W_shift; /* if true W is shifted when vectors converge */
171 : } dvdDashboard;
172 :
173 : typedef struct {
174 : /*------------------------- User parameters ---------------------------*/
175 : PetscInt blocksize; /* block size */
176 : PetscInt initialsize; /* initial size of V */
177 : PetscInt minv; /* size of V after restarting */
178 : PetscInt plusk; /* keep plusk eigenvectors from the last iteration */
179 : PetscBool ipB; /* true if B-ortho is used */
180 : PetscReal fix; /* the fix parameter */
181 : PetscBool krylovstart; /* true if the starting subspace is a Krylov basis */
182 : PetscBool dynamic; /* true if dynamic stopping criterion is used */
183 : PetscBool doubleexp; /* double expansion in GD (GD2) */
184 :
185 : /*----------------- Child objects and working data -------------------*/
186 : dvdDashboard ddb;
187 : } EPS_DAVIDSON;
188 :
189 939 : static inline PetscErrorCode EPSDavidsonFLAdd(dvdFunctionList **fl,dvdCallback f)
190 : {
191 939 : dvdFunctionList *l;
192 :
193 939 : PetscFunctionBegin;
194 939 : PetscCall(PetscNew(&l));
195 939 : l->f = f;
196 939 : l->next = *fl;
197 939 : *fl = l;
198 939 : PetscFunctionReturn(PETSC_SUCCESS);
199 : }
200 :
201 355 : static inline PetscErrorCode EPSDavidsonFLCall(dvdFunctionList *fl,dvdDashboard *d)
202 : {
203 355 : dvdFunctionList *l;
204 :
205 355 : PetscFunctionBegin;
206 1294 : for (l=fl;l;l=l->next) PetscCall(l->f(d));
207 355 : PetscFunctionReturn(PETSC_SUCCESS);
208 : }
209 :
210 495 : static inline PetscErrorCode EPSDavidsonFLDestroy(dvdFunctionList **fl)
211 : {
212 495 : dvdFunctionList *l,*l0;
213 :
214 495 : PetscFunctionBegin;
215 1434 : for (l=*fl;l;l=l0) {
216 939 : l0 = l->next;
217 939 : PetscCall(PetscFree(l));
218 : }
219 495 : *fl = NULL;
220 495 : PetscFunctionReturn(PETSC_SUCCESS);
221 : }
222 :
223 : /*
224 : The blackboard configuration structure: saves information about the memory
225 : and other requirements.
226 :
227 : The starting memory structure:
228 :
229 : V W? AV BV? tKZ
230 : |-----------|-----------|-----------|------------|-------|
231 : nev+mpd nev+mpd scP+mpd nev?+mpd sP+scP
232 : scP+mpd scP+mpd
233 :
234 : The final memory structure considering W_shift:
235 :
236 : cX V cY? W? cAV AV BcX? BV? KZ tKZ
237 : |---|-------|----|------|---|-------|----|-------|---|---|
238 : nev mpd nev mpd scP mpd nev mpd scP sP <- shift
239 : scP scP <- !shift
240 : */
241 : typedef struct {
242 : PetscInt max_size_V; /* max size of the searching subspace (mpd) */
243 : PetscInt max_size_X; /* max size of X (bs) */
244 : PetscInt size_V; /* real size of V (nev+size_P+mpd) */
245 : PetscInt max_size_oldX; /* max size of oldX */
246 : PetscInt max_nev; /* max number of converged pairs */
247 : PetscInt max_size_P; /* number of computed vectors for the projector */
248 : PetscInt max_size_cP; /* number of converged vectors in the projectors */
249 : PetscInt max_size_proj; /* max size projected problem */
250 : PetscInt max_size_cX_proj; /* max converged vectors in the projected problem */
251 : PetscInt state; /* method states:
252 : 0: preconfiguring
253 : 1: configuring
254 : 2: running */
255 : } dvdBlackboard;
256 :
257 : #define DVD_STATE_PRECONF 0
258 : #define DVD_STATE_CONF 1
259 : #define DVD_STATE_RUN 2
260 :
261 : /* Prototypes of non-static auxiliary functions */
262 : SLEPC_INTERN PetscErrorCode dvd_calcpairs_qz(dvdDashboard*,dvdBlackboard*,PetscBool,PetscBool);
263 : SLEPC_INTERN PetscErrorCode dvd_improvex_gd2(dvdDashboard*,dvdBlackboard*,KSP,PetscInt);
264 : SLEPC_INTERN PetscErrorCode dvd_improvex_jd(dvdDashboard*,dvdBlackboard*,KSP,PetscInt,PetscBool);
265 : SLEPC_INTERN PetscErrorCode dvd_improvex_jd_proj_uv(dvdDashboard*,dvdBlackboard*);
266 : SLEPC_INTERN PetscErrorCode dvd_improvex_jd_lit_const(dvdDashboard*,dvdBlackboard*,PetscInt,PetscReal,PetscReal);
267 : SLEPC_INTERN PetscErrorCode dvd_improvex_compute_X(dvdDashboard*,PetscInt,PetscInt,Vec*,PetscScalar*,PetscInt);
268 : SLEPC_INTERN PetscErrorCode dvd_initV(dvdDashboard*,dvdBlackboard*,PetscInt,PetscInt,PetscBool);
269 : SLEPC_INTERN PetscErrorCode dvd_orthV(BV,PetscInt,PetscInt);
270 : SLEPC_INTERN PetscErrorCode dvd_schm_basic_preconf(dvdDashboard*,dvdBlackboard*,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,HarmType_t,KSP,InitType_t,PetscBool,PetscBool,PetscBool);
271 : SLEPC_INTERN PetscErrorCode dvd_schm_basic_conf(dvdDashboard*,dvdBlackboard*,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,PetscInt,HarmType_t,PetscBool,PetscScalar,KSP,PetscReal,InitType_t,PetscBool,PetscBool,PetscBool,PetscBool);
272 : SLEPC_INTERN PetscErrorCode dvd_testconv_slepc(dvdDashboard*,dvdBlackboard*);
273 : SLEPC_INTERN PetscErrorCode dvd_managementV_basic(dvdDashboard*,dvdBlackboard*,PetscInt,PetscInt,PetscInt,PetscInt,PetscBool,PetscBool);
274 : SLEPC_INTERN PetscErrorCode dvd_static_precond_PC(dvdDashboard*,dvdBlackboard*,PC);
275 : SLEPC_INTERN PetscErrorCode dvd_harm_updateproj(dvdDashboard*);
276 : SLEPC_INTERN PetscErrorCode dvd_harm_conf(dvdDashboard*,dvdBlackboard*,HarmType_t,PetscBool,PetscScalar);
277 :
278 : /* Internal interface routines */
279 : SLEPC_INTERN PetscErrorCode EPSReset_XD(EPS);
280 : SLEPC_INTERN PetscErrorCode EPSSetUp_XD(EPS);
281 : SLEPC_INTERN PetscErrorCode EPSSolve_XD(EPS);
282 : SLEPC_INTERN PetscErrorCode EPSComputeVectors_XD(EPS);
283 : SLEPC_INTERN PetscErrorCode EPSXDSetKrylovStart_XD(EPS,PetscBool);
284 : SLEPC_INTERN PetscErrorCode EPSXDGetKrylovStart_XD(EPS,PetscBool*);
285 : SLEPC_INTERN PetscErrorCode EPSXDSetBlockSize_XD(EPS,PetscInt);
286 : SLEPC_INTERN PetscErrorCode EPSXDGetBlockSize_XD(EPS,PetscInt*);
287 : SLEPC_INTERN PetscErrorCode EPSXDSetRestart_XD(EPS,PetscInt,PetscInt);
288 : SLEPC_INTERN PetscErrorCode EPSXDGetRestart_XD(EPS,PetscInt*,PetscInt*);
289 : SLEPC_INTERN PetscErrorCode EPSXDGetInitialSize_XD(EPS,PetscInt*);
290 : SLEPC_INTERN PetscErrorCode EPSXDSetInitialSize_XD(EPS,PetscInt);
291 : SLEPC_INTERN PetscErrorCode EPSXDSetBOrth_XD(EPS,PetscBool);
292 : SLEPC_INTERN PetscErrorCode EPSXDGetBOrth_XD(EPS,PetscBool*);
293 : SLEPC_INTERN PetscErrorCode EPSJDGetFix_JD(EPS,PetscReal*);
294 : SLEPC_INTERN PetscErrorCode EPSJDGetConstCorrectionTol_JD(EPS,PetscBool*);
|