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-2012, Universitat 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 <slepc-private/epsimpl.h> /*I "slepceps.h" I*/
35: #include <slepc-private/stimpl.h> /*I "slepcst.h" I*/
36: #include <slepcblaslapack.h>
38: typedef struct _dvdFunctionList {
39: PetscErrorCode (*f)(void*);
40: void *d;
41: struct _dvdFunctionList *next;
42: } dvdFunctionList;
44: typedef PetscInt MatType_t;
45: #define DVD_MAT_HERMITIAN (1<<1)
46: #define DVD_MAT_NEG_DEF (1<<2)
47: #define DVD_MAT_POS_DEF (1<<3)
48: #define DVD_MAT_SINGULAR (1<<4)
49: #define DVD_MAT_COMPLEX (1<<5)
50: #define DVD_MAT_IMPLICIT (1<<6)
51: #define DVD_MAT_IDENTITY (1<<7)
52: #define DVD_MAT_DIAG (1<<8)
53: #define DVD_MAT_TRIANG (1<<9)
54: #define DVD_MAT_UTRIANG (1<<9)
55: #define DVD_MAT_LTRIANG (1<<10)
56: #define DVD_MAT_UNITARY (1<<11)
58: typedef PetscInt EPType_t;
59: #define DVD_EP_STD (1<<1)
60: #define DVD_EP_HERMITIAN (1<<2)
61: #define DVD_EP_INDEFINITE (1<<3)
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_KXX,
81: DVD_PROJ_KZX
82: } ProjType_t;
84: typedef enum {
85: DVD_METH_GD,
86: DVD_METH_JD,
87: DVD_METH_GD2
88: } Method_t;
90: typedef struct _dvdDashboard {
91: /**** Function steps ****/
92: /* Initialize V */
93: PetscErrorCode (*initV)(struct _dvdDashboard*);
94: void *initV_data;
95:
96: /* Find the approximate eigenpairs from V */
97: PetscErrorCode (*calcPairs)(struct _dvdDashboard*);
98: void *calcPairs_data;
100: /* Eigenpair test for convergence */
101: PetscBool (*testConv)(struct _dvdDashboard*, PetscScalar eigvr,
102: PetscScalar eigvi, PetscReal res, PetscReal *error);
103: void *testConv_data;
105: /* Number of converged eigenpairs */
106: PetscInt nconv;
108: /* Number of pairs ready to converge */
109: PetscInt npreconv;
111: /* Improve the selected eigenpairs */
112: PetscErrorCode (*improveX)(struct _dvdDashboard*, Vec *D, PetscInt max_size_D,
113: PetscInt r_s, PetscInt r_e, PetscInt *size_D);
114: void *improveX_data;
116: /* Check for restarting */
117: PetscBool (*isRestarting)(struct _dvdDashboard*);
118: void *isRestarting_data;
120: /* Perform restarting */
121: PetscErrorCode (*restartV)(struct _dvdDashboard*);
122: void *restartV_data;
124: /* Update V */
125: PetscErrorCode (*updateV)(struct _dvdDashboard*);
126: void *updateV_data;
128: /**** Problem specification ****/
129: Mat A, B; /* Problem matrices */
130: MatType_t sA, sB; /* Matrix specifications */
131: EPType_t sEP; /* Problem specifications */
132: PetscInt nev; /* number of eigenpairs */
133: EPSWhich which; /* spectrum selection */
134: PetscBool
135: withTarget; /* if there is a target */
136: PetscScalar
137: target[2]; /* target value */
138: PetscReal tol; /* tolerance */
139: PetscBool
140: correctXnorm; /* if true, norm of X are computed */
142: /**** Subspaces specification ****/
143: Vec *V, /* searching subspace */
144: *real_V, /* original V */
145: *W, /* testing subspace */
146: *real_W, /* original W */
147: *cX, /* converged right eigenvectors */
148: *cY, /* converged left eigenvectors */
149: *BcX, /* basis of B*cX */
150: *AV, /* A*V */
151: *real_AV, /* original A*V space */
152: *BV, /* B*V */
153: *real_BV, /* original B*V space */
154: *BDS; /* B * eps->DV */
155: PetscInt size_V, /* size of V */
156: size_real_V, /* original size of V */
157: size_W, /* size of W */
158: size_real_W, /* original size of W */
159: size_AV, /* size of AV */
160: size_real_AV, /* original size of AV */
161: size_BV, /* size of BV */
162: size_BDS, /* size of BDS */
163: size_real_BV, /* original size of BV */
164: size_cX, /* size of cX */
165: size_cY, /* size of cY */
166: size_D, /* active vectors */
167: size_BcX, /* size of BcX */
168: size_real_eigr, /* size of original eigr, eigi, nR, errest */
169: max_size_V, /* max size of V */
170: max_size_W, /* max size of W */
171: max_size_X, /* max new vectors in V */
172: max_size_AV, /* max size of AV */
173: max_size_BV, /* max size of BV */
174: max_size_proj, /* max size projected problem */
175: max_cX_in_proj, /* max vectors from cX in the projected problem */
176: max_cX_in_impr, /* max vectros from cX in the projector */
177: max_size_P, /* max unconverged vectors in the projector */
178: bs; /* max vectors that expands the subspace every iteration */
179: EPS eps; /* Connection to SLEPc */
180:
181: /**** Auxiliary space ****/
182: Vec *auxV; /* auxiliary vectors */
183: PetscScalar
184: *auxS; /* auxiliary scalars */
185: PetscInt
186: size_auxV, /* max size of auxV */
187: size_auxS; /* max size of auxS */
189: /**** Eigenvalues and errors ****/
190: PetscScalar
191: *ceigr, *ceigi, /* converged eigenvalues */
192: *eigr, *eigi, /* current eigenvalues */
193: *real_eigr,
194: *real_eigi; /* original eigr and eigi */
195: PetscReal
196: *nR, /* residual norm */
197: *real_nR, /* original nR */
198: *nX, /* X norm */
199: *real_nX, /* original nX */
200: *errest, /* relative error eigenpairs */
201: *real_errest, /* original errest */
202: *nBDS, /* B-norms of DS */
203: *nBV, /* B-norms of V */
204: *nBcX, /* B-norms of cX */
205: *nBpX, /* B-norms of pX */
206: *real_nBV; /* original nBDS, nBV and nBcX */
208: /**** Shared function and variables ****/
209: PetscErrorCode (*e_Vchanged)(struct _dvdDashboard*, PetscInt s_imm,
210: PetscInt e_imm, PetscInt s_new, PetscInt e_new);
211: void *e_Vchanged_data;
213: PetscErrorCode (*calcpairs_residual)(struct _dvdDashboard*, PetscInt s, PetscInt e,
214: Vec *R);
215: PetscErrorCode (*calcpairs_residual_eig)(struct _dvdDashboard*, PetscInt s, PetscInt e,
216: Vec *R);
217: PetscErrorCode (*calcpairs_selectPairs)(struct _dvdDashboard*, PetscInt n);
218: void *calcpairs_residual_data;
219: PetscErrorCode (*improvex_precond)(struct _dvdDashboard*, PetscInt i, Vec x,
220: Vec Px);
221: void *improvex_precond_data;
222: PetscErrorCode (*improvex_jd_proj_uv)(struct _dvdDashboard*, PetscInt i_s,
223: PetscInt i_e, Vec *u, Vec *v,
224: Vec *kr, Vec *auxV,
225: PetscScalar *theta,
226: PetscScalar *thetai,
227: PetscScalar *pX, PetscScalar *pY,
228: PetscInt ld);
229: PetscErrorCode (*improvex_jd_lit)(struct _dvdDashboard*, PetscInt i,
230: PetscScalar* theta, PetscScalar* thetai,
231: PetscInt *maxits, PetscReal *tol);
232: PetscErrorCode (*calcpairs_W)(struct _dvdDashboard*);
233: void *calcpairs_W_data;
234: PetscErrorCode (*calcpairs_proj_trans)(struct _dvdDashboard*);
235: PetscErrorCode (*calcpairs_eigs_trans)(struct _dvdDashboard*);
236: PetscErrorCode (*calcpairs_eig_backtrans)(struct _dvdDashboard*,PetscScalar,PetscScalar,PetscScalar*,PetscScalar*);
237: PetscErrorCode (*calcpairs_proj_res)(struct _dvdDashboard*, PetscInt r_s,
238: PetscInt r_e, Vec *R);
239: PetscErrorCode (*preTestConv)(struct _dvdDashboard*, PetscInt s, PetscInt pre,
240: PetscInt e, Vec *auxV, PetscScalar *auxS,
241: PetscInt *nConv);
243: PetscErrorCode (*e_newIteration)(struct _dvdDashboard*);
244: void *e_newIteration_data;
246: IP ipI;
247: IP ipV, /* orthogonal routine options for V subspace */
248: ipW; /* orthogonal routine options for W subspace */
250: dvdFunctionList
251: *startList, /* starting list */
252: *endList, /* ending list */
253: *destroyList; /* destructor list */
255: PetscScalar *H, /* Projected problem matrix A*/
256: *real_H, /* original H */
257: *G, /* Projected problem matrix B*/
258: *real_G, /* original G */
259: *S, /* first Schur matrix, S = pY'*H*pX */
260: *T, /* second Schur matrix, T = pY'*G*pX */
261: *cS, /* first Schur matrix of converged pairs */
262: *cT; /* second Schur matrix of converged pairs */
263: MatType_t
264: sH, /* H properties */
265: sG; /* G properties */
266: PetscInt ldH, /* leading dimension of H */
267: ldS, /* leading dimension of S */
268: ldT, /* leading dimension of T */
269: ldcS, /* leading dimension of cS */
270: ldcT, /* leading dimension of cT */
271: size_H, /* rows and columns in H */
272: size_G, /* rows and columns in G */
273: size_MT, /* rows in MT */
274: size_cS, /* dimension of cS */
275: size_cT, /* dimension of cT */
276: max_size_cS, /* max size of cS and cT */
277: cX_in_V, /* number of converged vectors in V */
278: cX_in_W, /* number of converged vectors in W */
279: cX_in_AV, /* number of converged vectors in AV */
280: cX_in_BV, /* number of converged vectors in BV */
281: cX_in_H, /* number of converged vectors in H */
282: cX_in_G; /* number of converged vectors in G */
284: PetscInt
285: V_tra_s,
286: V_tra_e, /* cX <- [cX V*MT(0:V_tra_s-1)], V <- V*MT(V_tra_s:V_tra_e) */
287: V_new_s,
288: V_new_e; /* added to V the columns V_new_s:V_new_e */
289: PetscBool
290: BV_shift, /* if true BV is shifted when vectors converge */
291: W_shift; /* if true W is shifted when vectors converge */
292: DS conv_ps, /* projected problem with the converged pairs */
293: ps; /* projected problem with the search subspace */
295: EPSOrthType orthoV_type;
297: void* prof_data; /* profiler data */
298: } dvdDashboard;
300: /* Add the function fun at the beginning of list */
301: #define DVD_FL_ADD_BEGIN(list, fun) { \
302: dvdFunctionList *fl=(list); \
304: PetscMalloc(sizeof(dvdFunctionList), &(list)); \
305: (list)->f = (PetscErrorCode(*)(void*))(fun); \
306: (list)->next = fl; }
308: /* Add the function fun at the end of list */
309: #define DVD_FL_ADD_END(list, fun) { \
310: if ((list)) {DVD_FL_ADD_END0(list, fun);} \
311: else {DVD_FL_ADD_BEGIN(list, fun);} }
313: #define DVD_FL_ADD_END0(list, fun) { \
314: dvdFunctionList *fl=(list); \
316: for(;fl->next; fl = fl->next); \
317: PetscMalloc(sizeof(dvdFunctionList), &fl->next); \
318: fl->next->f = (PetscErrorCode(*)(void*))(fun); \
319: fl->next->next = PETSC_NULL; }
321: #define DVD_FL_ADD(list, fun) DVD_FL_ADD_END(list, fun)
323: #define DVD_FL_CALL(list, arg0) { \
324: dvdFunctionList *fl; \
325: for(fl=(list); fl; fl=fl->next) \
326: if(*(dvdCallback)fl->f) (*(dvdCallback)fl->f)((arg0)); }
328: #define DVD_FL_DEL(list) { \
329: dvdFunctionList *fl=(list), *oldfl; \
331: while(fl) { \
332: oldfl = fl; fl = fl->next; PetscFree(oldfl); } \
333: (list) = PETSC_NULL;}
335: /*
336: The blackboard configuration structure: saves information about the memory
337: and other requirements.
339: The starting memory structure:
341: V W? AV BV? tKZ
342: |-----------|-----------|-----------|------------|-------|
343: nev+mpd nev+mpd scP+mpd nev?+mpd sP+scP
344: scP+mpd scP+mpd
346: The final memory structure considering W_shift and BV_shift:
348: cX V cY? W? cAV AV BcX? BV? KZ tKZ
349: |---|-------|----|------|---|-------|----|-------|---|---|
350: nev mpd nev mpd scP mpd nev mpd scP sP <- shift
351: scP scP <- !shift
352: */
353: typedef struct {
354: PetscInt max_size_V, /* max size of the searching subspace (mpd) */
355: max_size_X, /* max size of X (bs) */
356: size_V, /* real size of V (nev+size_P+mpd) */
357: max_size_oldX, /* max size of oldX */
358: max_size_auxV, /* max size of auxiliary vecs */
359: max_size_auxS, /* max size of auxiliary scalars */
360: max_nev, /* max number of converged pairs */
361: max_size_P, /* number of computed vectors for the projector */
362: max_size_cP, /* number of converged vectors in the projectors */
363: max_size_proj, /* max size projected problem */
364: max_size_cX_proj, /* max converged vectors in the projected problem */
365: own_vecs, /* number of global vecs */
366: own_scalars; /* number of local scalars */
367: Vec *free_vecs; /* free global vectors */
368: PetscScalar
369: *free_scalars; /* free scalars */
370: PetscInt state; /* method states:
371: 0: preconfiguring
372: 1: configuring
373: 2: running
374: */
375: } dvdBlackboard;
377: #define DVD_STATE_PRECONF 0
378: #define DVD_STATE_CONF 1
379: #define DVD_STATE_RUN 2
381: /* Shared types */
382: typedef PetscErrorCode (*dvdPrecond)(dvdDashboard*, PetscInt i, Vec x, Vec Px);
383: typedef PetscErrorCode (*dvdCallback)(dvdDashboard*);
384: typedef PetscErrorCode (*e_Vchanged_type)(dvdDashboard*, PetscInt s_imm,
385: PetscInt e_imm, PetscInt s_new, PetscInt e_new);
386: typedef PetscBool (*isRestarting_type)(dvdDashboard*);
387: typedef PetscErrorCode (*e_newIteration_type)(dvdDashboard*);
388: typedef PetscErrorCode (*improveX_type)(dvdDashboard*, Vec *D, PetscInt max_size_D,
389: PetscInt r_s, PetscInt r_e, PetscInt *size_D);
391: /* Structures for blas */
392: typedef PetscErrorCode (*DvdReductionPostF)(PetscScalar*,PetscInt,void*);
393: typedef struct {
394: PetscScalar
395: *out; /* final vector */
396: PetscInt
397: size_out; /* size of out */
398: DvdReductionPostF
399: f; /* function called after the reduction */
400: void *ptr;
401: } DvdReductionChunk;
403: typedef struct {
404: PetscScalar
405: *in, /* vector to sum-up with more nodes */
406: *out; /* final vector */
407: PetscInt size_in, /* size of in */
408: max_size_in; /* max size of in */
409: DvdReductionChunk
410: *ops; /* vector of reduction operations */
411: PetscInt
412: size_ops, /* size of ops */
413: max_size_ops; /* max size of ops */
414: MPI_Comm comm; /* MPI communicator */
415: } DvdReduction;
417: typedef struct {
418: PetscInt i0, i1, i2, ld, s0, e0, s1, e1;
419: PetscScalar *M;
420: } DvdMult_copy_func;
422: /* Routines for initV step */
423: PetscErrorCode dvd_initV(dvdDashboard *d, dvdBlackboard *b, PetscInt k,
424: PetscInt user, PetscBool krylov);
426: /* Routines for calcPairs step */
427: PetscErrorCode dvd_calcpairs_qz(dvdDashboard *d, dvdBlackboard *b,
428: EPSOrthType orth, IP ipI,
429: PetscInt cX_proj, PetscBool harm);
431: /* Routines for improveX step */
432: PetscErrorCode dvd_improvex_jd(dvdDashboard *d, dvdBlackboard *b, KSP ksp,
433: PetscInt max_bs, PetscInt cX_impr, PetscBool dynamic);
434: PetscErrorCode dvd_improvex_jd_proj_uv(dvdDashboard *d, dvdBlackboard *b,
435: ProjType_t p);
436: PetscErrorCode dvd_improvex_jd_lit_const(dvdDashboard *d, dvdBlackboard *b,
437: PetscInt maxits, PetscReal tol,
438: PetscReal fix);
439: PetscErrorCode dvd_improvex_gd2(dvdDashboard *d,dvdBlackboard *b,KSP ksp,PetscInt max_bs);
440: PetscErrorCode dvd_improvex_get_eigenvectors(dvdDashboard *d, PetscScalar *pX,
441: PetscScalar *pY, PetscInt ld, PetscScalar *auxS, PetscInt size_auxS);
443: /* Routines for testConv step */
444: PetscErrorCode dvd_testconv_basic(dvdDashboard *d, dvdBlackboard *b);
445: PetscErrorCode dvd_testconv_slepc(dvdDashboard *d, dvdBlackboard *b);
447: /* Routines for management of V */
448: PetscErrorCode dvd_managementV_basic(dvdDashboard *d, dvdBlackboard *b,
449: PetscInt bs, PetscInt mpd,
450: PetscInt min_size_V,
451: PetscInt plusk, PetscBool harm,
452: PetscBool allResiduals);
454: /* Some utilities */
455: PetscErrorCode dvd_static_precond_PC(dvdDashboard *d, dvdBlackboard *b, PC pc);
456: PetscErrorCode dvd_jacobi_precond(dvdDashboard *d, dvdBlackboard *b);
457: PetscErrorCode dvd_profiler(dvdDashboard *d, dvdBlackboard *b);
458: PetscErrorCode dvd_prof_init();
459: PetscErrorCode dvd_harm_conf(dvdDashboard *d, dvdBlackboard *b,
460: HarmType_t mode, PetscBool fixedTarget,
461: PetscScalar t);
463: /* Methods */
464: PetscErrorCode dvd_schm_basic_preconf(dvdDashboard *d, dvdBlackboard *b,
465: PetscInt mpd, PetscInt min_size_V, PetscInt bs,
466: PetscInt ini_size_V, PetscInt size_initV, PetscInt plusk,
467: HarmType_t harmMode, KSP ksp, InitType_t init, PetscBool allResiduals,
468: EPSOrthType orth, PetscInt cX_proj, PetscInt cX_impr, Method_t method);
469: PetscErrorCode dvd_schm_basic_conf(dvdDashboard *d, dvdBlackboard *b,
470: PetscInt mpd, PetscInt min_size_V, PetscInt bs,
471: PetscInt ini_size_V, PetscInt size_initV, PetscInt plusk,
472: IP ip, HarmType_t harmMode, PetscBool fixedTarget, PetscScalar t, KSP ksp,
473: PetscReal fix, InitType_t init, PetscBool allResiduals, EPSOrthType orth,
474: PetscInt cX_proj, PetscInt cX_impr, PetscBool dynamic, Method_t method);
476: /* BLAS routines */
477: PetscErrorCode SlepcDenseMatProd(PetscScalar *C, PetscInt _ldC, PetscScalar b,
478: PetscScalar a,
479: const PetscScalar *A, PetscInt _ldA, PetscInt rA, PetscInt cA, PetscBool At,
480: const PetscScalar *B, PetscInt _ldB, PetscInt rB, PetscInt cB, PetscBool Bt);
481: PetscErrorCode SlepcDenseMatProdTriang(
482: PetscScalar *C, MatType_t sC, PetscInt ldC,
483: const PetscScalar *A, MatType_t sA, PetscInt ldA, PetscInt rA, PetscInt cA,
484: PetscBool At,
485: const PetscScalar *B, MatType_t sB, PetscInt ldB, PetscInt rB, PetscInt cB,
486: PetscBool Bt);
487: PetscErrorCode SlepcDenseNorm(PetscScalar *A, PetscInt ldA, PetscInt _rA,
488: PetscInt cA, PetscScalar *eigi);
489: PetscErrorCode SlepcDenseCopy(PetscScalar *Y, PetscInt ldY, PetscScalar *X,
490: PetscInt ldX, PetscInt rX, PetscInt cX);
491: PetscErrorCode SlepcDenseCopyTriang(PetscScalar *Y, MatType_t sY, PetscInt ldY,
492: PetscScalar *X, MatType_t sX, PetscInt ldX,
493: PetscInt rX, PetscInt cX);
494: PetscErrorCode SlepcUpdateVectorsZ(Vec *Y, PetscScalar beta, PetscScalar alpha,
495: Vec *X, PetscInt cX, const PetscScalar *M, PetscInt ldM, PetscInt rM,
496: PetscInt cM);
497: PetscErrorCode SlepcUpdateVectorsS(Vec *Y, PetscInt dY, PetscScalar beta,
498: PetscScalar alpha, Vec *X, PetscInt cX, PetscInt dX, const PetscScalar *M,
499: PetscInt ldM, PetscInt rM, PetscInt cM);
500: PetscErrorCode SlepcUpdateVectorsD(Vec *X, PetscInt cX, PetscScalar alpha,
501: const PetscScalar *M, PetscInt ldM, PetscInt rM, PetscInt cM,
502: PetscScalar *work, PetscInt lwork);
503: PetscErrorCode VecsMult(PetscScalar *M, MatType_t sM, PetscInt ldM,
504: Vec *U, PetscInt sU, PetscInt eU,
505: Vec *V, PetscInt sV, PetscInt eV,
506: PetscScalar *workS0, PetscScalar *workS1);
507: PetscErrorCode VecsMultS(PetscScalar *M, MatType_t sM, PetscInt ldM,
508: Vec *U, PetscInt sU, PetscInt eU,
509: Vec *V, PetscInt sV, PetscInt eV, DvdReduction *r,
510: DvdMult_copy_func *sr);
511: PetscErrorCode VecsMultIb(PetscScalar *M, MatType_t sM, PetscInt ldM,
512: PetscInt rM, PetscInt cM, PetscScalar *auxS,
513: Vec V);
514: PetscErrorCode VecsMultIa(PetscScalar *M, MatType_t sM, PetscInt ldM,
515: Vec *U, PetscInt sU, PetscInt eU,
516: Vec *V, PetscInt sV, PetscInt eV);
517: PetscErrorCode SlepcAllReduceSumBegin(DvdReductionChunk *ops,
518: PetscInt max_size_ops,
519: PetscScalar *in, PetscScalar *out,
520: PetscInt max_size_in, DvdReduction *r,
521: MPI_Comm comm);
522: PetscErrorCode SlepcAllReduceSum(DvdReduction *r, PetscInt size_in,
523: DvdReductionPostF f, void *ptr,
524: PetscScalar **in);
525: PetscErrorCode SlepcAllReduceSumEnd(DvdReduction *r);
526: PetscErrorCode dvd_orthV(IP ip, Vec *DS, PetscInt size_DS, Vec *cX,
527: PetscInt size_cX, Vec *V, PetscInt V_new_s,
528: PetscInt V_new_e, PetscScalar *auxS,
529: PetscRandom rand);
530: PetscErrorCode dvd_BorthV_faster(IP ip, Vec *DS, Vec *BDS,PetscReal *BDSn, PetscInt size_DS, Vec *cX,
531: Vec *BcX, PetscReal *BcXn, PetscInt size_cX, Vec *V, Vec *BV,PetscReal *BVn,
532: PetscInt V_new_s, PetscInt V_new_e,
533: PetscScalar *auxS, PetscRandom rand);
534: PetscErrorCode dvd_BorthV_stable(IP ip,Vec *defl,PetscReal *BDSn,PetscInt size_DS,Vec *cX,PetscReal *BcXn,PetscInt size_cX,Vec *V,PetscReal *BVn,PetscInt V_new_s,PetscInt V_new_e,PetscScalar *auxS,PetscRandom rand);
536: /* SLEPc interface routines */
537: PetscErrorCode SLEPcNotImplemented();
538: PetscErrorCode EPSCreate_Davidson(EPS eps);
539: PetscErrorCode EPSReset_Davidson(EPS eps);
540: PetscErrorCode EPSSetUp_Davidson(EPS eps);
541: PetscErrorCode EPSSolve_Davidson(EPS eps);
542: PetscErrorCode EPSComputeVectors_Davidson(EPS eps);
543: PetscErrorCode EPSDavidsonSetKrylovStart_Davidson(EPS eps,PetscBool krylovstart);
544: PetscErrorCode EPSDavidsonGetKrylovStart_Davidson(EPS eps,PetscBool *krylovstart);
545: PetscErrorCode EPSDavidsonSetBlockSize_Davidson(EPS eps,PetscInt blocksize);
546: PetscErrorCode EPSDavidsonGetBlockSize_Davidson(EPS eps,PetscInt *blocksize);
547: PetscErrorCode EPSDavidsonSetRestart_Davidson(EPS eps,PetscInt minv,PetscInt plusk);
548: PetscErrorCode EPSDavidsonGetRestart_Davidson(EPS eps,PetscInt *minv,PetscInt *plusk);
549: PetscErrorCode EPSDavidsonGetInitialSize_Davidson(EPS eps,PetscInt *initialsize);
550: PetscErrorCode EPSDavidsonSetInitialSize_Davidson(EPS eps,PetscInt initialsize);
551: PetscErrorCode EPSDavidsonGetFix_Davidson(EPS eps,PetscReal *fix);
552: PetscErrorCode EPSDavidsonSetFix_Davidson(EPS eps,PetscReal fix);
553: PetscErrorCode EPSDavidsonSetBOrth_Davidson(EPS eps,EPSOrthType borth);
554: PetscErrorCode EPSDavidsonGetBOrth_Davidson(EPS eps,EPSOrthType *borth);
555: PetscErrorCode EPSDavidsonSetConstantCorrectionTolerance_Davidson(EPS eps,PetscBool constant);
556: PetscErrorCode EPSDavidsonGetConstantCorrectionTolerance_Davidson(EPS eps,PetscBool *constant);
557: PetscErrorCode EPSDavidsonSetWindowSizes_Davidson(EPS eps,PetscInt pwindow,PetscInt qwindow);
558: PetscErrorCode EPSDavidsonGetWindowSizes_Davidson(EPS eps,PetscInt *pwindow,PetscInt *qwindow);
559: PetscErrorCode EPSDavidsonSetMethod_Davidson(EPS eps,Method_t method);
560: PetscErrorCode EPSDavidsonGetMethod_Davidson(EPS eps,Method_t *method);
561: PetscErrorCode EPSGDGetDoubleExpansion_GD(EPS,PetscBool*);
562: PetscErrorCode EPSGDSetDoubleExpansion_GD(EPS,PetscBool);
564: /* Common inline function */
567: PETSC_STATIC_INLINE PetscErrorCode dvd_improvex_compute_X(dvdDashboard *d,PetscInt i_s,PetscInt i_e,Vec *u,PetscScalar *pX,PetscInt ld)
568: {
569: PetscErrorCode ierr;
570: PetscInt n = i_e - i_s, i;
573: SlepcUpdateVectorsZ(u,0.0,1.0,d->V-d->cX_in_H,d->size_V+d->cX_in_H,&pX[ld*i_s],ld,d->size_H,n);
574: /* nX(i) <- ||X(i)|| */
575: if (d->correctXnorm) {
576: for (i=0; i<n; i++) {
577: VecNormBegin(u[i],NORM_2,&d->nX[i_s+i]);
578: }
579: for (i=0; i<n; i++) {
580: VecNormEnd(u[i],NORM_2,&d->nX[i_s+i]);
581: }
582: #if !defined(PETSC_USE_COMPLEX)
583: for(i=0; i<n; i++) {
584: if(d->eigi[i_s+i] != 0.0) {
585: d->nX[i_s+i] = d->nX[i_s+i+1] = PetscSqrtScalar(d->nX[i_s+i]*d->nX[i_s+i]+d->nX[i_s+i+1]*d->nX[i_s+i+1]);
586: i++;
587: }
588: }
589: #endif
590: } else {
591: for (i=0; i<n; i++) d->nX[i_s+i] = 1.0;
592: }
593: return(0);
594: }
597: #define _Ceil(A,B) ((A)/(B)+((A)%(B)==0?0:1))
598: #define FromIntToScalar(S) ((PetscInt)_Ceil((S)*sizeof(PetscBLASInt),sizeof(PetscScalar)))
599: #define FromRealToScalar(S) ((PetscInt)_Ceil((S)*sizeof(PetscReal),sizeof(PetscScalar)))