Actual source code: davidson.h
1: /*
2: Method: General Davidson Method (includes GD and JD)
4: References:
5: - Ernest R. Davidson. Super-matrix methods. Computer Physics Communications,
6: 53:49–60, May 1989.
8: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
9: SLEPc - Scalable Library for Eigenvalue Problem Computations
10: Copyright (c) 2002-2010, Universidad Politecnica de Valencia, Spain
12: This file is part of SLEPc.
13:
14: SLEPc is free software: you can redistribute it and/or modify it under the
15: terms of version 3 of the GNU Lesser General Public License as published by
16: the Free Software Foundation.
18: SLEPc is distributed in the hope that it will be useful, but WITHOUT ANY
19: WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
20: FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for
21: more details.
23: You should have received a copy of the GNU Lesser General Public License
24: along with SLEPc. If not, see <http://www.gnu.org/licenses/>.
25: - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - - -
26: */
29: /*
30: Dashboard struct: contains the methods that will be employed and the tunning
31: options.
32: */
34: #include "petsc.h"
35: #include private/epsimpl.h
36: #include private/stimpl.h
39: typedef struct _dvdFunctionList {
40: PetscErrorCode (*f)(void*);
41: void *d;
42: struct _dvdFunctionList *next;
43: } dvdFunctionList;
45: typedef PetscInt MatType_t;
46: #define DVD_MAT_HERMITIAN (1<<1)
47: #define DVD_MAT_NEG_DEF (1<<2)
48: #define DVD_MAT_POS_DEF (1<<3)
49: #define DVD_MAT_SINGULAR (1<<4)
50: #define DVD_MAT_COMPLEX (1<<5)
51: #define DVD_MAT_IMPLICIT (1<<6)
52: #define DVD_MAT_IDENTITY (1<<7)
53: #define DVD_MAT_DIAG (1<<8)
54: #define DVD_MAT_TRIANG (1<<9)
55: #define DVD_MAT_UTRIANG (1<<9)
56: #define DVD_MAT_LTRIANG (1<<10)
57: #define DVD_MAT_UNITARY (1<<11)
59: typedef PetscInt EPType_t;
60: #define DVD_EP_STD (1<<1)
61: #define DVD_EP_HERMITIAN (1<<2)
63: #define DVD_IS(T,P) ((T) & (P))
64: #define DVD_ISNOT(T,P) (((T) & (P)) ^ (P))
66: typedef enum {
67: DVD_HARM_NONE,
68: DVD_HARM_RR,
69: DVD_HARM_RRR,
70: DVD_HARM_REIGS,
71: DVD_HARM_LEIGS
72: } HarmType_t;
74: typedef enum {
75: DVD_INITV_CLASSIC,
76: DVD_INITV_KRYLOV
77: } InitType_t;
79: typedef enum {
80: DVD_PROJ_KBXX,
81: DVD_PROJ_KBXY,
82: DVD_PROJ_KBXZ,
83: DVD_PROJ_KBXZY
84: } ProjType_t;
86: typedef enum {
87: DVD_MT_IDENTITY,/* without transformation */
88: DVD_MT_pX, /* using the projected problem eigenvectors */
89: DVD_MT_ORTHO /* using an orthonormal transformation */
90: } MT_type_t;
92: typedef enum {
93: DVD_ORTHOV_NONE,/* V isn't orthonormalized */
94: DVD_ORTHOV_I, /* V is orthonormalized */
95: DVD_ORTHOV_B /* V is B-orthonormalized */
96: } orthoV_type_t;
99: typedef struct _dvdDashboard {
100: /**** Function steps ****/
101: /* Initialize V */
102: PetscErrorCode (*initV)(struct _dvdDashboard*);
103: void *initV_data;
104:
105: /* Find the approximate eigenpairs from V */
106: PetscErrorCode (*calcPairs)(struct _dvdDashboard*);
107: void *calcPairs_data;
109: /* Eigenpair test for convergence */
110: PetscTruth (*testConv)(struct _dvdDashboard*, PetscScalar eigvr,
111: PetscScalar eigvi, PetscReal res, PetscReal *error);
112: void *testConv_data;
114: /* Number of converged eigenpairs */
115: PetscInt nconv;
117: /* Number of pairs ready to converge */
118: PetscInt npreconv;
120: /* Improve the selected eigenpairs */
121: PetscErrorCode (*improveX)(struct _dvdDashboard*, Vec *D, PetscInt max_size_D,
122: PetscInt r_s, PetscInt r_e, PetscInt *size_D);
123: void *improveX_data;
125: /* Check for restarting */
126: PetscTruth (*isRestarting)(struct _dvdDashboard*);
127: void *isRestarting_data;
129: /* Perform restarting */
130: PetscErrorCode (*restartV)(struct _dvdDashboard*);
131: void *restartV_data;
133: /* Update V */
134: PetscErrorCode (*updateV)(struct _dvdDashboard*);
135: void *updateV_data;
137: /**** Problem specification ****/
138: Mat A, B; /* Problem matrices */
139: MatType_t sA, sB; /* Matrix specifications */
140: EPType_t sEP; /* Problem specifications */
141: PetscInt nev; /* number of eigenpairs */
142: EPSWhich which; /* spectrum selection */
143: PetscTruth
144: withTarget; /* if there is a target */
145: PetscScalar
146: target[2]; /* target value */
147: PetscReal tol; /* tolerance */
148: PetscTruth
149: correctXnorm; /* if true, tol < |r|/|x| */
151: /**** Subspaces specification ****/
152: Vec *V, /* searching subspace */
153: *W, /* testing subspace */
154: *cX, /* converged right eigenvectors */
155: *cY, /* converged left eigenvectors */
156: *BcX, /* basis of B*cX */
157: *AV, /* A*V */
158: *real_AV, /* original A*V space */
159: *BV, /* B*V */
160: *real_BV; /* original B*V space */
161: PetscInt size_V, /* size of V */
162: size_AV, /* size of AV */
163: size_BV, /* size of BV */
164: size_cX, /* size of cX */
165: size_cY, /* size of cY */
166: size_D, /* active vectors */
167: max_size_V, /* max size of V */
168: max_size_X, /* max size of X */
169: max_size_AV, /* max size of AV */
170: max_size_BV; /* max size of BV */
171: EPS eps; /* Connection to SLEPc */
172:
173: /**** Auxiliary space ****/
174: Vec *auxV; /* auxiliary vectors */
175: PetscScalar
176: *auxS; /* auxiliary scalars */
177: PetscInt
178: size_auxV, /* max size of auxV */
179: size_auxS; /* max size of auxS */
181: /**** Eigenvalues and errors ****/
182: PetscScalar
183: *ceigr, *ceigi, /* converged eigenvalues */
184: *eigr, *eigi; /* current eigenvalues */
185: PetscReal
186: *nR, /* residual norm */
187: *nX, /* X norm */
188: *errest; /* relative error eigenpairs */
190: /**** Shared function and variables ****/
191: PetscErrorCode (*e_Vchanged)(struct _dvdDashboard*, PetscInt s_imm,
192: PetscInt e_imm, PetscInt s_new, PetscInt e_new);
193: void *e_Vchanged_data;
195: PetscErrorCode (*calcpairs_residual)(struct _dvdDashboard*, PetscInt s, PetscInt e,
196: Vec *R, PetscScalar *auxS, Vec auxV);
197: PetscErrorCode (*calcpairs_selectPairs)(struct _dvdDashboard*, PetscInt n);
198: void *calcpairs_residual_data;
199: PetscErrorCode (*calcpairs_X)(struct _dvdDashboard*, PetscInt s, PetscInt e,
200: Vec *X);
201: PetscErrorCode (*calcpairs_Y)(struct _dvdDashboard*, PetscInt s, PetscInt e,
202: Vec *Y);
203: PetscErrorCode (*improvex_precond)(struct _dvdDashboard*, PetscInt i, Vec x,
204: Vec Px);
205: void *improvex_precond_data;
206: PetscErrorCode (*improvex_jd_proj_uv)(struct _dvdDashboard*, PetscInt i_s,
207: PetscInt i_e, Vec **u, Vec **v,
208: Vec **kr, Vec **auxV_,
209: PetscScalar *theta,
210: PetscScalar *thetai,
211: PetscScalar *pX, PetscScalar *pY,
212: PetscInt ld);
213: PetscErrorCode (*improvex_jd_lit)(struct _dvdDashboard*, PetscInt i,
214: PetscScalar* theta, PetscScalar* thetai,
215: PetscInt *maxits, PetscReal *tol);
216: PetscErrorCode (*calcpairs_W)(struct _dvdDashboard*);
217: void *calcpairs_W_data;
218: PetscErrorCode (*calcpairs_proj_trans)(struct _dvdDashboard*);
219: PetscErrorCode (*calcpairs_eigs_trans)(struct _dvdDashboard*);
220: PetscErrorCode (*calcpairs_proj_res)(struct _dvdDashboard*, PetscInt r_s,
221: PetscInt r_e, Vec *R);
222: PetscErrorCode (*preTestConv)(struct _dvdDashboard*, PetscInt s, PetscInt pre,
223: PetscInt e, Vec *auxV, PetscScalar *auxS,
224: PetscInt *nConv);
226: PetscErrorCode (*e_newIteration)(struct _dvdDashboard*);
227: void *e_newIteration_data;
229: IP ipI;
230: IP ipV, /* orthogonal routine options for V subspace */
231: ipW; /* orthogonal routine options for W subspace */
233: dvdFunctionList
234: *startList, /* starting list */
235: *endList, /* ending list */
236: *destroyList; /* destructor list */
238: PetscScalar *H, /* Projected problem matrix A*/
239: *real_H, /* original H */
240: *G, /* Projected problem matrix B*/
241: *real_G, /* original G */
242: *pX, /* projected problem right eigenvectors */
243: *pY, /* projected problem left eigenvectors */
244: *MTX, /* right transformation matrix */
245: *MTY, /* left transformation matrix */
246: *S, /* first Schur matrix, S = pY'*H*pX */
247: *T, /* second Schur matrix, T = pY'*G*pX */
248: *cS, /* first Schur matrix of converged pairs */
249: *cT; /* second Schur matrix of converged pairs */
250: MatType_t
251: pX_type, /* pX properties */
252: pY_type, /* pY properties */
253: sH, /* H properties */
254: sG; /* G properties */
255: PetscInt ldH, /* leading dimension of H */
256: ldpX, /* leading dimension of pX */
257: ldpY, /* leading dimension of pY */
258: ldMTX, /* leading dimension of MTX */
259: ldMTY, /* leading dimension of MTY */
260: ldS, /* leading dimension of S */
261: ldT, /* leading dimension of T */
262: ldcS, /* leading dimension of cS */
263: ldcT, /* leading dimension of cT */
264: size_H, /* rows and columns in H */
265: size_G, /* rows and columns in G */
266: size_MT; /* rows in MT */
268: PetscInt V_imm_s,
269: V_imm_e, /* unchanged V columns, V_imm_s:V_imm_e-1 */
270: V_tra_s,
271: V_tra_e, /* V(V_tra_e:) = V*MT(V_tra_s:V_tra_e-1) */
272: V_new_s,
273: V_new_e; /* added to V the columns V_new_s:V_new_e */
275: MT_type_t MT_type;
276: orthoV_type_t orthoV_type;
278: PetscRandom rand; /* random seed */
279: void* prof_data; /* profiler data */
280: } dvdDashboard;
282: #define DVD_FL_ADD(list, fun) { \
283: dvdFunctionList *fl=(list); \
285: PetscMalloc(sizeof(dvdFunctionList), &(list)); \
286: (list)->f = (PetscErrorCode(*)(void*))(fun); \
287: (list)->next = fl; }
289: #define DVD_FL_CALL(list, arg0) { \
290: dvdFunctionList *fl; \
291: for(fl=(list); fl; fl=fl->next) (*(dvdCallback)fl->f)((arg0)); }
293: #define DVD_FL_DEL(list) { \
294: dvdFunctionList *fl=(list), *oldfl; \
296: while(fl) { \
297: oldfl = fl; fl = fl->next; PetscFree(oldfl); }}
299: /*
300: The blackboard configuration structure: saves information about the memory
301: and other requirements
302: */
303: typedef struct {
304: PetscInt max_size_V, /* max size of the searching subspace */
305: max_size_X, /* max size of X */
306: size_V, /* real size of V */
307: max_size_oldX, /* max size of oldX */
308: max_size_auxV, /* max size of auxiliary vecs */
309: max_size_auxS, /* max size of auxiliary scalars */
310: max_nev, /* max number of converged pairs */
311: own_vecs, /* number of global vecs */
312: own_scalars; /* number of local scalars */
313: Vec *free_vecs; /* free global vectors */
314: PetscScalar
315: *free_scalars; /* free scalars */
316: PetscInt state; /* method states:
317: 0: preconfiguring
318: 1: configuring
319: 2: running
320: */
321: } dvdBlackboard;
323: #define DVD_STATE_PRECONF 0
324: #define DVD_STATE_CONF 1
325: #define DVD_STATE_RUN 2
327: /* Shared types */
328: typedef void* dvdPrecondData; // DEPRECATED!!
329: typedef PetscErrorCode (*dvdPrecond)(dvdDashboard*, PetscInt i, Vec x, Vec Px);
330: typedef PetscErrorCode (*dvdCallback)(dvdDashboard*);
331: typedef PetscErrorCode (*e_Vchanged_type)(dvdDashboard*, PetscInt s_imm,
332: PetscInt e_imm, PetscInt s_new, PetscInt e_new);
333: typedef PetscTruth (*isRestarting_type)(dvdDashboard*);
334: typedef PetscErrorCode (*e_newIteration_type)(dvdDashboard*);
335: typedef PetscErrorCode (*improveX_type)(dvdDashboard*, Vec *D, PetscInt max_size_D,
336: PetscInt r_s, PetscInt r_e, PetscInt *size_D);
338: /* Structures for blas */
339: typedef PetscErrorCode (*DvdReductionPostF)(PetscScalar*,PetscInt,void*);
340: typedef struct {
341: PetscScalar
342: *out; /* final vector */
343: PetscInt
344: size_out; /* size of out */
345: DvdReductionPostF
346: f; /* function called after the reduction */
347: void *ptr;
348: } DvdReductionChunk;
350: typedef struct {
351: PetscScalar
352: *in, /* vector to sum-up with more nodes */
353: *out; /* final vector */
354: PetscInt size_in, /* size of in */
355: max_size_in; /* max size of in */
356: DvdReductionChunk
357: *ops; /* vector of reduction operations */
358: PetscInt
359: size_ops, /* size of ops */
360: max_size_ops; /* max size of ops */
361: MPI_Comm comm; /* MPI communicator */
362: } DvdReduction;
364: typedef struct {
365: PetscInt i0, i1, i2, ld, s0, e0, s1, e1;
366: PetscScalar *M;
367: } DvdMult_copy_func;
369: /* Routines for initV step */
370: PetscErrorCode dvd_initV_classic(dvdDashboard *d, dvdBlackboard *b, PetscInt k);
371: PetscErrorCode dvd_initV_user(dvdDashboard *d, dvdBlackboard *b,
372: PetscInt size_userV, PetscInt k);
373: PetscErrorCode dvd_initV_krylov(dvdDashboard *d, dvdBlackboard *b, PetscInt k);
375: /* Routines for calcPairs step */
376: PetscErrorCode dvd_calcpairs_rr(dvdDashboard *d, dvdBlackboard *b);
377: PetscErrorCode dvd_calcpairs_qz(dvdDashboard *d, dvdBlackboard *b, IP ip);
379: /* Routines for improveX step */
380: PetscErrorCode dvd_improvex_jd(dvdDashboard *d, dvdBlackboard *b, KSP ksp,
381: PetscInt max_bs);
382: PetscErrorCode dvd_improvex_jd_proj_uv(dvdDashboard *d, dvdBlackboard *b,
383: ProjType_t p);
384: PetscErrorCode dvd_improvex_jd_lit_const(dvdDashboard *d, dvdBlackboard *b,
385: PetscInt maxits, PetscReal tol,
386: PetscReal fix);
388: /* Routines for testConv step */
389: PetscErrorCode dvd_testconv_basic(dvdDashboard *d, dvdBlackboard *b);
390: PetscErrorCode dvd_testconv_slepc(dvdDashboard *d, dvdBlackboard *b);
392: /* Routines for management of V */
393: PetscErrorCode dvd_managementV_basic(dvdDashboard *d, dvdBlackboard *b,
394: PetscInt bs, PetscInt max_size_V,
395: PetscInt mpd, PetscInt min_size_V,
396: PetscInt plusk, PetscTruth harm,
397: PetscTruth allResiduals);
399: /* Some utilities */
400: PetscErrorCode dvd_static_precond_PC(dvdDashboard *d, dvdBlackboard *b, PC pc);
401: PetscErrorCode dvd_jacobi_precond(dvdDashboard *d, dvdBlackboard *b);
402: PetscErrorCode dvd_profiler(dvdDashboard *d, dvdBlackboard *b);
403: PetscErrorCode dvd_prof_init();
404: PetscErrorCode dvd_harm_conf(dvdDashboard *d, dvdBlackboard *b,
405: HarmType_t mode, PetscTruth fixedTarget,
406: PetscScalar t);
408: /* Methods */
409: PetscErrorCode dvd_schm_basic_preconf(dvdDashboard *d, dvdBlackboard *b,
410: PetscInt max_size_V, PetscInt mpd, PetscInt min_size_V, PetscInt bs,
411: PetscInt ini_size_V, PetscInt size_initV, PetscInt plusk, PC pc,
412: HarmType_t harmMode, KSP ksp, InitType_t init, PetscTruth allResiduals);
413: PetscErrorCode dvd_schm_basic_conf(dvdDashboard *d, dvdBlackboard *b,
414: PetscInt max_size_V, PetscInt mpd, PetscInt min_size_V, PetscInt bs,
415: PetscInt ini_size_V, PetscInt size_initV, PetscInt plusk, PC pc,
416: IP ip, HarmType_t harmMode, PetscTruth fixedTarget, PetscScalar t, KSP ksp,
417: PetscReal fix, InitType_t init, PetscTruth allResiduals);
419: /* BLAS routines */
420: PetscErrorCode SlepcDenseMatProd(PetscScalar *C, PetscInt _ldC, PetscScalar b,
421: PetscScalar a,
422: const PetscScalar *A, PetscInt _ldA, PetscInt rA, PetscInt cA, PetscTruth At,
423: const PetscScalar *B, PetscInt _ldB, PetscInt rB, PetscInt cB, PetscTruth Bt);
424: PetscErrorCode SlepcDenseMatProdTriang(
425: PetscScalar *C, MatType_t sC, PetscInt ldC,
426: const PetscScalar *A, MatType_t sA, PetscInt ldA, PetscInt rA, PetscInt cA,
427: PetscTruth At,
428: const PetscScalar *B, MatType_t sB, PetscInt ldB, PetscInt rB, PetscInt cB,
429: PetscTruth Bt);
430: PetscErrorCode SlepcDenseNorm(PetscScalar *A, PetscInt ldA, PetscInt _rA,
431: PetscInt cA, PetscScalar *eigi);
432: PetscErrorCode SlepcDenseOrth(PetscScalar *A, PetscInt _ldA, PetscInt _rA,
433: PetscInt _cA, PetscScalar *auxS, PetscInt _lauxS,
434: PetscInt *ncA);
435: PetscErrorCode SlepcDenseCopy(PetscScalar *Y, PetscInt ldY, PetscScalar *X,
436: PetscInt ldX, PetscInt rX, PetscInt cX);
437: PetscErrorCode SlepcDenseCopyTriang(PetscScalar *Y, MatType_t sY, PetscInt ldY,
438: PetscScalar *X, MatType_t sX, PetscInt ldX,
439: PetscInt rX, PetscInt cX);
440: PetscErrorCode SlepcUpdateVectorsZ(Vec *Y, PetscScalar beta, PetscScalar alpha,
441: Vec *X, PetscInt cX, const PetscScalar *M, PetscInt ldM, PetscInt rM,
442: PetscInt cM);
443: PetscErrorCode SlepcUpdateVectorsS(Vec *Y, PetscInt dY, PetscScalar beta,
444: PetscScalar alpha, Vec *X, PetscInt cX, PetscInt dX, const PetscScalar *M,
445: PetscInt ldM, PetscInt rM, PetscInt cM);
446: PetscErrorCode SlepcUpdateVectorsD(Vec *X, PetscInt cX, PetscScalar alpha,
447: const PetscScalar *M, PetscInt ldM, PetscInt rM, PetscInt cM,
448: PetscScalar *work, PetscInt lwork);
449: PetscErrorCode VecsMult(PetscScalar *M, MatType_t sM, PetscInt ldM,
450: Vec *U, PetscInt sU, PetscInt eU,
451: Vec *V, PetscInt sV, PetscInt eV,
452: PetscScalar *workS0, PetscScalar *workS1);
453: PetscErrorCode VecsMultS(PetscScalar *M, MatType_t sM, PetscInt ldM,
454: Vec *U, PetscInt sU, PetscInt eU,
455: Vec *V, PetscInt sV, PetscInt eV, DvdReduction *r,
456: DvdMult_copy_func *sr);
457: PetscErrorCode VecsMultIb(PetscScalar *M, MatType_t sM, PetscInt ldM,
458: PetscInt rM, PetscInt cM, PetscScalar *auxS,
459: Vec V);
460: PetscErrorCode VecsMultIa(PetscScalar *M, MatType_t sM, PetscInt ldM,
461: Vec *U, PetscInt sU, PetscInt eU,
462: Vec *V, PetscInt sV, PetscInt eV);
463: PetscErrorCode SlepcAllReduceSumBegin(DvdReductionChunk *ops,
464: PetscInt max_size_ops,
465: PetscScalar *in, PetscScalar *out,
466: PetscInt max_size_in, DvdReduction *r,
467: MPI_Comm comm);
468: PetscErrorCode SlepcAllReduceSum(DvdReduction *r, PetscInt size_in,
469: DvdReductionPostF f, void *ptr,
470: PetscScalar **in);
471: PetscErrorCode SlepcAllReduceSumEnd(DvdReduction *r);
472: PetscErrorCode dvd_orthV(IP ip, Vec *DS, PetscInt size_DS, Vec *cX,
473: PetscInt size_cX, Vec *V, PetscInt V_new_s,
474: PetscInt V_new_e, PetscScalar *auxS, Vec auxV,
475: PetscRandom rand);
476: PetscErrorCode dvd_compute_eigenvectors(PetscInt n_, PetscScalar *S,
477: PetscInt ldS_, PetscScalar *T, PetscInt ldT_, PetscScalar *pX,
478: PetscInt ldpX_, PetscScalar *pY, PetscInt ldpY_, PetscScalar *auxS,
479: PetscInt size_auxS, PetscTruth doProd);
481: /* SLEPc interface routines */
482: PetscErrorCode SLEPcNotImplemented();
483: PetscErrorCode EPSCreate_DAVIDSON(EPS eps);
484: PetscErrorCode EPSDestroy_DAVIDSON(EPS eps);
485: PetscErrorCode EPSSetUp_DAVIDSON(EPS eps);
486: PetscErrorCode EPSSolve_DAVIDSON(EPS eps);
487: PetscErrorCode EPSComputeVectors_QZ(EPS eps);
488: PetscErrorCode EPSDAVIDSONSetKrylovStart_DAVIDSON(EPS eps,PetscTruth krylovstart);
489: PetscErrorCode EPSDAVIDSONGetKrylovStart_DAVIDSON(EPS eps,PetscTruth *krylovstart);
490: PetscErrorCode EPSDAVIDSONSetBlockSize_DAVIDSON(EPS eps,PetscInt blocksize);
491: PetscErrorCode EPSDAVIDSONGetBlockSize_DAVIDSON(EPS eps,PetscInt *blocksize);
492: PetscErrorCode EPSDAVIDSONSetRestart_DAVIDSON(EPS eps,PetscInt minv,PetscInt plusk);
493: PetscErrorCode EPSDAVIDSONGetRestart_DAVIDSON(EPS eps,PetscInt *minv,PetscInt *plusk);
494: PetscErrorCode EPSDAVIDSONGetInitialSize_DAVIDSON(EPS eps,PetscInt *initialsize);
495: PetscErrorCode EPSDAVIDSONSetInitialSize_DAVIDSON(EPS eps,PetscInt initialsize);
496: PetscErrorCode EPSDAVIDSONGetFix_DAVIDSON(EPS eps,PetscReal *fix);
497: PetscErrorCode EPSDAVIDSONSetFix_DAVIDSON(EPS eps,PetscReal fix);
499: typedef struct {
500: /**** Solver options ****/
501: PetscInt blocksize, /* block size */
502: initialsize, /* initial size of V */
503: minv, /* size of V after restarting */
504: plusk; /* keep plusk eigenvectors from the last iteration */
505: PetscTruth ipB; /* truth if V'B*V=I */
506: PetscInt method; /* method for improving the approximate solution */
507: PetscReal fix; /* the fix parameter */
508: PetscTruth krylovstart; /* true if the starting subspace is a Krylov basis */
510: /**** Solver data ****/
511: dvdDashboard ddb;
513: /**** Things to destroy ****/
514: PetscScalar *wS;
515: Vec *wV;
516: PetscInt size_wV;
517: PC pc; /* pc extracted from st->ksp */
518: } EPS_DAVIDSON;