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 dvd_blas_prof_init();
421: PetscErrorCode SlepcDenseMatProd(PetscScalar *C, PetscInt _ldC, PetscScalar b,
422:   PetscScalar a,
423:   const PetscScalar *A, PetscInt _ldA, PetscInt rA, PetscInt cA, PetscTruth At,
424:   const PetscScalar *B, PetscInt _ldB, PetscInt rB, PetscInt cB, PetscTruth Bt);
425: PetscErrorCode SlepcDenseMatProdTriang(
426:   PetscScalar *C, MatType_t sC, PetscInt ldC,
427:   const PetscScalar *A, MatType_t sA, PetscInt ldA, PetscInt rA, PetscInt cA,
428:   PetscTruth At,
429:   const PetscScalar *B, MatType_t sB, PetscInt ldB, PetscInt rB, PetscInt cB,
430:   PetscTruth Bt);
431: PetscErrorCode SlepcDenseNorm(PetscScalar *A, PetscInt ldA, PetscInt _rA,
432:                               PetscInt cA, PetscScalar *eigi);
433: PetscErrorCode SlepcDenseOrth(PetscScalar *A, PetscInt _ldA, PetscInt _rA,
434:                               PetscInt _cA, PetscScalar *auxS, PetscInt _lauxS,
435:                               PetscInt *ncA);
436: PetscErrorCode SlepcDenseCopy(PetscScalar *Y, PetscInt ldY, PetscScalar *X,
437:                               PetscInt ldX, PetscInt rX, PetscInt cX);
438: PetscErrorCode SlepcDenseCopyTriang(PetscScalar *Y, MatType_t sY, PetscInt ldY,
439:                                     PetscScalar *X, MatType_t sX, PetscInt ldX,
440:                                     PetscInt rX, PetscInt cX);
441: PetscErrorCode SlepcUpdateVectorsZ(Vec *Y, PetscScalar beta, PetscScalar alpha,
442:   Vec *X, PetscInt cX, const PetscScalar *M, PetscInt ldM, PetscInt rM,
443:   PetscInt cM);
444: PetscErrorCode SlepcUpdateVectorsS(Vec *Y, PetscInt dY, PetscScalar beta,
445:   PetscScalar alpha, Vec *X, PetscInt cX, PetscInt dX, const PetscScalar *M,
446:   PetscInt ldM, PetscInt rM, PetscInt cM);
447: PetscErrorCode SlepcUpdateVectorsD(Vec *X, PetscInt cX, PetscScalar alpha,
448:   const PetscScalar *M, PetscInt ldM, PetscInt rM, PetscInt cM,
449:   PetscScalar *work, PetscInt lwork);
450: PetscErrorCode VecsMult(PetscScalar *M, MatType_t sM, PetscInt ldM,
451:                         Vec *U, PetscInt sU, PetscInt eU,
452:                         Vec *V, PetscInt sV, PetscInt eV,
453:                         PetscScalar *workS0, PetscScalar *workS1);
454: PetscErrorCode VecsMultS(PetscScalar *M, MatType_t sM, PetscInt ldM,
455:                          Vec *U, PetscInt sU, PetscInt eU,
456:                          Vec *V, PetscInt sV, PetscInt eV, DvdReduction *r,
457:                          DvdMult_copy_func *sr);
458: PetscErrorCode VecsMultIb(PetscScalar *M, MatType_t sM, PetscInt ldM,
459:                           PetscInt rM, PetscInt cM, PetscScalar *auxS,
460:                           Vec V);
461: PetscErrorCode VecsMultIa(PetscScalar *M, MatType_t sM, PetscInt ldM,
462:                           Vec *U, PetscInt sU, PetscInt eU,
463:                           Vec *V, PetscInt sV, PetscInt eV);
464: PetscErrorCode SlepcAllReduceSumBegin(DvdReductionChunk *ops,
465:                                       PetscInt max_size_ops,
466:                                       PetscScalar *in, PetscScalar *out,
467:                                       PetscInt max_size_in, DvdReduction *r,
468:                                       MPI_Comm comm);
469: PetscErrorCode SlepcAllReduceSum(DvdReduction *r, PetscInt size_in,
470:                                  DvdReductionPostF f, void *ptr,
471:                                  PetscScalar **in);
472: PetscErrorCode SlepcAllReduceSumEnd(DvdReduction *r);
473: PetscErrorCode dvd_orthV(IP ip, Vec *DS, PetscInt size_DS, Vec *cX,
474:                          PetscInt size_cX, Vec *V, PetscInt V_new_s,
475:                          PetscInt V_new_e, PetscScalar *auxS, Vec auxV,
476:                          PetscRandom rand);
477: PetscErrorCode dvd_compute_eigenvectors(PetscInt n_, PetscScalar *S,
478:   PetscInt ldS_, PetscScalar *T, PetscInt ldT_, PetscScalar *pX,
479:   PetscInt ldpX_, PetscScalar *pY, PetscInt ldpY_, PetscScalar *auxS,
480:   PetscInt size_auxS, PetscTruth doProd);

482: /* SLEPc interface routines */
483: PetscErrorCode SLEPcNotImplemented();
484: PetscErrorCode EPSCreate_DAVIDSON(EPS eps);
485: PetscErrorCode EPSDestroy_DAVIDSON(EPS eps);
486: PetscErrorCode EPSSetUp_DAVIDSON(EPS eps);
487: PetscErrorCode EPSSolve_DAVIDSON(EPS eps);
488: PetscErrorCode EPSComputeVectors_QZ(EPS eps);
489: PetscErrorCode EPSDAVIDSONSetKrylovStart_DAVIDSON(EPS eps,PetscTruth krylovstart);
490: PetscErrorCode EPSDAVIDSONGetKrylovStart_DAVIDSON(EPS eps,PetscTruth *krylovstart);
491: PetscErrorCode EPSDAVIDSONSetBlockSize_DAVIDSON(EPS eps,PetscInt blocksize);
492: PetscErrorCode EPSDAVIDSONGetBlockSize_DAVIDSON(EPS eps,PetscInt *blocksize);
493: PetscErrorCode EPSDAVIDSONSetRestart_DAVIDSON(EPS eps,PetscInt minv,PetscInt plusk);
494: PetscErrorCode EPSDAVIDSONGetRestart_DAVIDSON(EPS eps,PetscInt *minv,PetscInt *plusk);
495: PetscErrorCode EPSDAVIDSONGetInitialSize_DAVIDSON(EPS eps,PetscInt *initialsize);
496: PetscErrorCode EPSDAVIDSONSetInitialSize_DAVIDSON(EPS eps,PetscInt initialsize);
497: PetscErrorCode EPSDAVIDSONGetFix_DAVIDSON(EPS eps,PetscReal *fix);
498: PetscErrorCode EPSDAVIDSONSetFix_DAVIDSON(EPS eps,PetscReal fix);

500: typedef struct {
501:   /**** Solver options ****/
502:   PetscInt blocksize,     /* block size */
503:     initialsize,          /* initial size of V */
504:     minv,                 /* size of V after restarting */
505:     plusk;                /* keep plusk eigenvectors from the last iteration */
506:   PetscTruth ipB;         /* truth if V'B*V=I */
507:   PetscInt   method;      /* method for improving the approximate solution */
508:   PetscReal  fix;         /* the fix parameter */
509:   PetscTruth krylovstart; /* true if the starting subspace is a Krylov basis */

511:   /**** Solver data ****/
512:   dvdDashboard ddb;

514:   /**** Things to destroy ****/
515:   PetscScalar *wS;
516:   Vec         *wV;
517:   PetscInt    size_wV;
518:   PC          pc;         /* pc extracted from st->ksp */
519: } EPS_DAVIDSON;