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-2013, Universitat Politecnica de Valencia, Spain
12: This file is part of SLEPc.
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: */
28: /*
29: Dashboard struct: contains the methods that will be employed and the tunning
30: options.
31: */
33: #include <slepc-private/epsimpl.h> /*I "slepceps.h" I*/
34: #include <slepc-private/stimpl.h> /*I "slepcst.h" I*/
35: #include <slepcblaslapack.h>
37: typedef struct _dvdFunctionList {
38: PetscErrorCode (*f)(void*);
39: void *d;
40: struct _dvdFunctionList *next;
41: } dvdFunctionList;
43: typedef PetscInt MatType_t;
44: #define DVD_MAT_HERMITIAN (1<<1)
45: #define DVD_MAT_NEG_DEF (1<<2)
46: #define DVD_MAT_POS_DEF (1<<3)
47: #define DVD_MAT_SINGULAR (1<<4)
48: #define DVD_MAT_COMPLEX (1<<5)
49: #define DVD_MAT_IMPLICIT (1<<6)
50: #define DVD_MAT_IDENTITY (1<<7)
51: #define DVD_MAT_DIAG (1<<8)
52: #define DVD_MAT_TRIANG (1<<9)
53: #define DVD_MAT_UTRIANG (1<<9)
54: #define DVD_MAT_LTRIANG (1<<10)
55: #define DVD_MAT_UNITARY (1<<11)
57: typedef PetscInt EPType_t;
58: #define DVD_EP_STD (1<<1)
59: #define DVD_EP_HERMITIAN (1<<2)
60: #define DVD_EP_INDEFINITE (1<<3)
62: #define DVD_IS(T,P) ((T) & (P))
63: #define DVD_ISNOT(T,P) (((T) & (P)) ^ (P))
65: typedef enum {
66: DVD_HARM_NONE,
67: DVD_HARM_RR,
68: DVD_HARM_RRR,
69: DVD_HARM_REIGS,
70: DVD_HARM_LEIGS
71: } HarmType_t;
73: typedef enum {
74: DVD_INITV_CLASSIC,
75: DVD_INITV_KRYLOV
76: } InitType_t;
78: typedef enum {
79: DVD_PROJ_KXX,
80: DVD_PROJ_KZX
81: } ProjType_t;
83: typedef enum {
84: DVD_METH_GD,
85: DVD_METH_JD,
86: DVD_METH_GD2
87: } Method_t;
89: typedef struct _dvdDashboard {
90: /**** Function steps ****/
91: /* Initialize V */
92: PetscErrorCode (*initV)(struct _dvdDashboard*);
93: void *initV_data;
95: /* Find the approximate eigenpairs from V */
96: PetscErrorCode (*calcPairs)(struct _dvdDashboard*);
97: void *calcPairs_data;
99: /* Eigenpair test for convergence */
100: PetscBool (*testConv)(struct _dvdDashboard*,PetscScalar eigvr,PetscScalar eigvi,PetscReal res,PetscReal *error);
101: void *testConv_data;
103: /* Number of converged eigenpairs */
104: PetscInt nconv;
106: /* Number of pairs ready to converge */
107: PetscInt npreconv;
109: /* Improve the selected eigenpairs */
110: PetscErrorCode (*improveX)(struct _dvdDashboard*,Vec *D,PetscInt max_size_D,PetscInt r_s,PetscInt r_e,PetscInt *size_D);
111: void *improveX_data;
113: /* Check for restarting */
114: PetscBool (*isRestarting)(struct _dvdDashboard*);
115: void *isRestarting_data;
117: /* Perform restarting */
118: PetscErrorCode (*restartV)(struct _dvdDashboard*);
119: void *restartV_data;
121: /* Update V */
122: PetscErrorCode (*updateV)(struct _dvdDashboard*);
123: void *updateV_data;
125: /**** Problem specification ****/
126: Mat A, B; /* Problem matrices */
127: MatType_t sA, sB; /* Matrix specifications */
128: EPType_t sEP; /* Problem specifications */
129: PetscInt nev; /* number of eigenpairs */
130: EPSWhich which; /* spectrum selection */
131: PetscBool
132: withTarget; /* if there is a target */
133: PetscScalar
134: target[2]; /* target value */
135: PetscReal tol; /* tolerance */
136: PetscBool
137: correctXnorm; /* if true, norm of X are computed */
139: /**** Subspaces specification ****/
140: Vec *V, /* searching subspace */
141: *real_V, /* original V */
142: *W, /* testing subspace */
143: *real_W, /* original W */
144: *cX, /* converged right eigenvectors */
145: *cY, /* converged left eigenvectors */
146: *BcX, /* basis of B*cX */
147: *AV, /* A*V */
148: *real_AV, /* original A*V space */
149: *BV, /* B*V */
150: *real_BV, /* original B*V space */
151: *BDS; /* B * eps->DV */
152: PetscInt size_V, /* size of V */
153: size_real_V, /* original size of V */
154: size_W, /* size of W */
155: size_real_W, /* original size of W */
156: size_AV, /* size of AV */
157: size_real_AV, /* original size of AV */
158: size_BV, /* size of BV */
159: size_BDS, /* size of BDS */
160: size_real_BV, /* original size of BV */
161: size_cX, /* size of cX */
162: size_cY, /* size of cY */
163: size_D, /* active vectors */
164: size_BcX, /* size of BcX */
165: size_real_eigr, /* size of original eigr, eigi, nR, errest */
166: max_size_V, /* max size of V */
167: max_size_W, /* max size of W */
168: max_size_X, /* max new vectors in V */
169: max_size_AV, /* max size of AV */
170: max_size_BV, /* max size of BV */
171: max_size_proj, /* max size projected problem */
172: max_cX_in_proj, /* max vectors from cX in the projected problem */
173: max_cX_in_impr, /* max vectros from cX in the projector */
174: max_size_P, /* max unconverged vectors in the projector */
175: bs; /* max vectors that expands the subspace every iteration */
176: EPS eps; /* Connection to SLEPc */
178: /**** Auxiliary space ****/
179: Vec *auxV; /* auxiliary vectors */
180: PetscScalar
181: *auxS; /* auxiliary scalars */
182: PetscInt
183: size_auxV, /* max size of auxV */
184: size_auxS; /* max size of auxS */
186: /**** Eigenvalues and errors ****/
187: PetscScalar
188: *ceigr, *ceigi, /* converged eigenvalues */
189: *eigr, *eigi, /* current eigenvalues */
190: *real_eigr,
191: *real_eigi; /* original eigr and eigi */
192: PetscReal
193: *nR, /* residual norm */
194: *real_nR, /* original nR */
195: *nX, /* X norm */
196: *real_nX, /* original nX */
197: *errest, /* relative error eigenpairs */
198: *real_errest, /* original errest */
199: *nBDS, /* B-norms of DS */
200: *nBV, /* B-norms of V */
201: *nBcX, /* B-norms of cX */
202: *nBpX, /* B-norms of pX */
203: *real_nBV; /* original nBDS, nBV and nBcX */
205: /**** Shared function and variables ****/
206: PetscErrorCode (*e_Vchanged)(struct _dvdDashboard*,PetscInt s_imm,PetscInt e_imm,PetscInt s_new,PetscInt e_new);
207: void *e_Vchanged_data;
209: PetscErrorCode (*calcpairs_residual)(struct _dvdDashboard*,PetscInt s,PetscInt e,Vec *R);
210: PetscErrorCode (*calcpairs_residual_eig)(struct _dvdDashboard*,PetscInt s,PetscInt e,Vec *R);
211: PetscErrorCode (*calcpairs_selectPairs)(struct _dvdDashboard*,PetscInt n);
212: void *calcpairs_residual_data;
213: PetscErrorCode (*improvex_precond)(struct _dvdDashboard*,PetscInt i,Vec x,Vec Px);
214: void *improvex_precond_data;
215: PetscErrorCode (*improvex_jd_proj_uv)(struct _dvdDashboard*,PetscInt i_s,PetscInt i_e,Vec *u,Vec *v,Vec *kr,Vec *auxV,PetscScalar *theta,PetscScalar *thetai,PetscScalar *pX,PetscScalar *pY,PetscInt ld);
216: PetscErrorCode (*improvex_jd_lit)(struct _dvdDashboard*,PetscInt i,PetscScalar* theta,PetscScalar* thetai,PetscInt *maxits,PetscReal *tol);
217: PetscErrorCode (*calcpairs_W)(struct _dvdDashboard*);
218: void *calcpairs_W_data;
219: PetscErrorCode (*calcpairs_proj_trans)(struct _dvdDashboard*);
220: PetscErrorCode (*calcpairs_eigs_trans)(struct _dvdDashboard*);
221: PetscErrorCode (*calcpairs_eig_backtrans)(struct _dvdDashboard*,PetscScalar,PetscScalar,PetscScalar*,PetscScalar*);
222: PetscErrorCode (*calcpairs_proj_res)(struct _dvdDashboard*,PetscInt r_s,PetscInt r_e,Vec *R);
223: PetscErrorCode (*preTestConv)(struct _dvdDashboard*,PetscInt s,PetscInt pre,PetscInt e,Vec *auxV,PetscScalar *auxS,PetscInt *nConv);
225: PetscErrorCode (*e_newIteration)(struct _dvdDashboard*);
226: void *e_newIteration_data;
228: IP ipI;
229: IP ipV, /* orthogonal routine options for V subspace */
230: ipW; /* orthogonal routine options for W subspace */
232: dvdFunctionList
233: *startList, /* starting list */
234: *endList, /* ending list */
235: *destroyList; /* destructor list */
237: PetscScalar *H, /* Projected problem matrix A*/
238: *real_H, /* original H */
239: *G, /* Projected problem matrix B*/
240: *real_G, /* original G */
241: *S, /* first Schur matrix, S = pY'*H*pX */
242: *T, /* second Schur matrix, T = pY'*G*pX */
243: *cS, /* first Schur matrix of converged pairs */
244: *cT; /* second Schur matrix of converged pairs */
245: MatType_t
246: sH, /* H properties */
247: sG; /* G properties */
248: PetscInt ldH, /* leading dimension of H */
249: ldS, /* leading dimension of S */
250: ldT, /* leading dimension of T */
251: ldcS, /* leading dimension of cS */
252: ldcT, /* leading dimension of cT */
253: size_H, /* rows and columns in H */
254: size_G, /* rows and columns in G */
255: size_MT, /* rows in MT */
256: size_cS, /* dimension of cS */
257: size_cT, /* dimension of cT */
258: max_size_cS, /* max size of cS and cT */
259: cX_in_V, /* number of converged vectors in V */
260: cX_in_W, /* number of converged vectors in W */
261: cX_in_AV, /* number of converged vectors in AV */
262: cX_in_BV, /* number of converged vectors in BV */
263: cX_in_H, /* number of converged vectors in H */
264: cX_in_G; /* number of converged vectors in G */
266: PetscInt
267: V_tra_s,
268: V_tra_e, /* cX <- [cX V*MT(0:V_tra_s-1)], V <- V*MT(V_tra_s:V_tra_e) */
269: V_new_s,
270: V_new_e; /* added to V the columns V_new_s:V_new_e */
271: PetscBool
272: BV_shift, /* if true BV is shifted when vectors converge */
273: W_shift; /* if true W is shifted when vectors converge */
274: DS conv_ps, /* projected problem with the converged pairs */
275: ps; /* projected problem with the search subspace */
277: EPSOrthType orthoV_type;
279: void* prof_data; /* profiler data */
280: } dvdDashboard;
282: /* Add the function fun at the beginning of list */
283: #define DVD_FL_ADD_BEGIN(list, fun) { \
284: dvdFunctionList *fl=(list); \
286: PetscMalloc(sizeof(dvdFunctionList), &(list)); \
287: (list)->f = (PetscErrorCode(*)(void*))(fun); \
288: (list)->next = fl; \
289: }
291: /* Add the function fun at the end of list */
292: #define DVD_FL_ADD_END(list, fun) { \
293: if ((list)) {DVD_FL_ADD_END0(list, fun);} \
294: else {DVD_FL_ADD_BEGIN(list, fun);} }
296: #define DVD_FL_ADD_END0(list, fun) { \
297: dvdFunctionList *fl=(list); \
299: for (;fl->next; fl = fl->next); \
300: PetscMalloc(sizeof(dvdFunctionList), &fl->next); \
301: fl->next->f = (PetscErrorCode(*)(void*))(fun); \
302: fl->next->next = NULL; \
303: }
305: #define DVD_FL_ADD(list, fun) DVD_FL_ADD_END(list, fun)
307: #define DVD_FL_CALL(list, arg0) { \
308: dvdFunctionList *fl; \
309: for (fl=(list); fl; fl=fl->next) \
310: if (*(dvdCallback)fl->f) (*(dvdCallback)fl->f)((arg0)); \
311: }
313: #define DVD_FL_DEL(list) { \
314: dvdFunctionList *fl=(list), *oldfl; \
316: while (fl) { \
317: oldfl = fl; fl = fl->next; PetscFree(oldfl);\
318: } \
319: (list) = NULL; \
320: }
322: /*
323: The blackboard configuration structure: saves information about the memory
324: and other requirements.
326: The starting memory structure:
328: V W? AV BV? tKZ
329: |-----------|-----------|-----------|------------|-------|
330: nev+mpd nev+mpd scP+mpd nev?+mpd sP+scP
331: scP+mpd scP+mpd
333: The final memory structure considering W_shift and BV_shift:
335: cX V cY? W? cAV AV BcX? BV? KZ tKZ
336: |---|-------|----|------|---|-------|----|-------|---|---|
337: nev mpd nev mpd scP mpd nev mpd scP sP <- shift
338: scP scP <- !shift
339: */
340: typedef struct {
341: PetscInt max_size_V, /* max size of the searching subspace (mpd) */
342: max_size_X, /* max size of X (bs) */
343: size_V, /* real size of V (nev+size_P+mpd) */
344: max_size_oldX, /* max size of oldX */
345: max_size_auxV, /* max size of auxiliary vecs */
346: max_size_auxS, /* max size of auxiliary scalars */
347: max_nev, /* max number of converged pairs */
348: max_size_P, /* number of computed vectors for the projector */
349: max_size_cP, /* number of converged vectors in the projectors */
350: max_size_proj, /* max size projected problem */
351: max_size_cX_proj, /* max converged vectors in the projected problem */
352: own_vecs, /* number of global vecs */
353: own_scalars; /* number of local scalars */
354: Vec *free_vecs; /* free global vectors */
355: PetscScalar
356: *free_scalars; /* free scalars */
357: PetscInt state; /* method states:
358: 0: preconfiguring
359: 1: configuring
360: 2: running
361: */
362: } dvdBlackboard;
364: #define DVD_STATE_PRECONF 0
365: #define DVD_STATE_CONF 1
366: #define DVD_STATE_RUN 2
368: /* Shared types */
369: typedef PetscErrorCode (*dvdPrecond)(dvdDashboard*,PetscInt i,Vec x,Vec Px);
370: typedef PetscErrorCode (*dvdCallback)(dvdDashboard*);
371: typedef PetscErrorCode (*e_Vchanged_type)(dvdDashboard*,PetscInt s_imm,PetscInt e_imm,PetscInt s_new,PetscInt e_new);
372: typedef PetscBool (*isRestarting_type)(dvdDashboard*);
373: typedef PetscErrorCode (*e_newIteration_type)(dvdDashboard*);
374: typedef PetscErrorCode (*improveX_type)(dvdDashboard*,Vec *D,PetscInt max_size_D,PetscInt r_s,PetscInt r_e,PetscInt *size_D);
376: /* Structures for blas */
377: typedef PetscErrorCode (*DvdReductionPostF)(PetscScalar*,PetscInt,void*);
378: typedef struct {
379: PetscScalar
380: *out; /* final vector */
381: PetscInt
382: size_out; /* size of out */
383: DvdReductionPostF
384: f; /* function called after the reduction */
385: void *ptr;
386: } DvdReductionChunk;
388: typedef struct {
389: PetscScalar
390: *in, /* vector to sum-up with more nodes */
391: *out; /* final vector */
392: PetscInt size_in, /* size of in */
393: max_size_in; /* max size of in */
394: DvdReductionChunk
395: *ops; /* vector of reduction operations */
396: PetscInt
397: size_ops, /* size of ops */
398: max_size_ops; /* max size of ops */
399: MPI_Comm comm; /* MPI communicator */
400: } DvdReduction;
402: typedef struct {
403: PetscInt i0, i1, i2, ld, s0, e0, s1, e1;
404: PetscScalar *M;
405: } DvdMult_copy_func;
407: /* Routines for initV step */
408: PETSC_INTERN PetscErrorCode dvd_initV(dvdDashboard *d,dvdBlackboard *b,PetscInt k,PetscInt user,PetscBool krylov);
410: /* Routines for calcPairs step */
411: PETSC_INTERN PetscErrorCode dvd_calcpairs_qz(dvdDashboard *d,dvdBlackboard *b,EPSOrthType orth,IP ipI,PetscInt cX_proj,PetscBool harm);
413: /* Routines for improveX step */
414: PETSC_INTERN PetscErrorCode dvd_improvex_jd(dvdDashboard *d,dvdBlackboard *b,KSP ksp,PetscInt max_bs,PetscInt cX_impr,PetscBool dynamic);
415: PETSC_INTERN PetscErrorCode dvd_improvex_jd_proj_uv(dvdDashboard *d,dvdBlackboard *b,ProjType_t p);
416: PETSC_INTERN PetscErrorCode dvd_improvex_jd_lit_const(dvdDashboard *d,dvdBlackboard *b,PetscInt maxits,PetscReal tol,PetscReal fix);
417: PETSC_INTERN PetscErrorCode dvd_improvex_gd2(dvdDashboard *d,dvdBlackboard *b,KSP ksp,PetscInt max_bs);
418: PETSC_INTERN PetscErrorCode dvd_improvex_get_eigenvectors(dvdDashboard *d,PetscScalar *pX,
419: PetscScalar *pY,PetscInt ld,PetscScalar *auxS,PetscInt size_auxS);
421: /* Routines for testConv step */
422: PETSC_INTERN PetscErrorCode dvd_testconv_basic(dvdDashboard *d,dvdBlackboard *b);
423: PETSC_INTERN PetscErrorCode dvd_testconv_slepc(dvdDashboard *d,dvdBlackboard *b);
425: /* Routines for management of V */
426: PETSC_INTERN PetscErrorCode dvd_managementV_basic(dvdDashboard *d,dvdBlackboard *b,PetscInt bs,PetscInt mpd,PetscInt min_size_V,PetscInt plusk,PetscBool harm,PetscBool allResiduals);
428: /* Some utilities */
429: PETSC_INTERN PetscErrorCode dvd_static_precond_PC(dvdDashboard *d,dvdBlackboard *b,PC pc);
430: PETSC_INTERN PetscErrorCode dvd_jacobi_precond(dvdDashboard *d,dvdBlackboard *b);
431: PETSC_INTERN PetscErrorCode dvd_profiler(dvdDashboard *d,dvdBlackboard *b);
432: PETSC_INTERN PetscErrorCode dvd_prof_init();
433: PETSC_INTERN PetscErrorCode dvd_harm_conf(dvdDashboard *d,dvdBlackboard *b,HarmType_t mode,PetscBool fixedTarget,PetscScalar t);
435: /* Methods */
436: PETSC_INTERN PetscErrorCode dvd_schm_basic_preconf(dvdDashboard *d,dvdBlackboard *b,PetscInt mpd,PetscInt min_size_V,PetscInt bs,PetscInt ini_size_V,PetscInt size_initV,PetscInt plusk,HarmType_t harmMode,KSP ksp,InitType_t init,PetscBool allResiduals,EPSOrthType orth,PetscInt cX_proj,PetscInt cX_impr,Method_t method);
437: PETSC_INTERN PetscErrorCode dvd_schm_basic_conf(dvdDashboard *d,dvdBlackboard *b,PetscInt mpd,PetscInt min_size_V,PetscInt bs,PetscInt ini_size_V,PetscInt size_initV,PetscInt plusk,IP ip,HarmType_t harmMode,PetscBool fixedTarget,PetscScalar t,KSP ksp,PetscReal fix,InitType_t init,PetscBool allResiduals,EPSOrthType orth,PetscInt cX_proj,PetscInt cX_impr,PetscBool dynamic,Method_t method);
439: /* BLAS routines */
440: PETSC_INTERN PetscErrorCode SlepcDenseMatProdTriang(PetscScalar *C,MatType_t sC,PetscInt ldC,const PetscScalar *A,MatType_t sA,PetscInt ldA,PetscInt rA,PetscInt cA,PetscBool At,const PetscScalar *B,MatType_t sB,PetscInt ldB,PetscInt rB,PetscInt cB,PetscBool Bt);
441: PETSC_INTERN PetscErrorCode SlepcDenseNorm(PetscScalar *A,PetscInt ldA,PetscInt _rA,PetscInt cA,PetscScalar *eigi);
442: PETSC_INTERN PetscErrorCode SlepcDenseCopy(PetscScalar *Y,PetscInt ldY,PetscScalar *X,PetscInt ldX,PetscInt rX,PetscInt cX);
443: PETSC_INTERN PetscErrorCode SlepcDenseCopyTriang(PetscScalar *Y,MatType_t sY,PetscInt ldY,PetscScalar *X,MatType_t sX,PetscInt ldX,PetscInt rX,PetscInt cX);
444: PETSC_INTERN PetscErrorCode SlepcUpdateVectorsZ(Vec *Y,PetscScalar beta,PetscScalar alpha,Vec *X,PetscInt cX,const PetscScalar *M,PetscInt ldM,PetscInt rM,PetscInt cM);
445: PETSC_INTERN PetscErrorCode SlepcUpdateVectorsS(Vec *Y,PetscInt dY,PetscScalar beta,PetscScalar alpha,Vec *X,PetscInt cX,PetscInt dX,const PetscScalar *M,PetscInt ldM,PetscInt rM,PetscInt cM);
446: PETSC_INTERN PetscErrorCode SlepcUpdateVectorsD(Vec *X,PetscInt cX,PetscScalar alpha,const PetscScalar *M,PetscInt ldM,PetscInt rM,PetscInt cM,PetscScalar *work,PetscInt lwork);
447: PETSC_INTERN PetscErrorCode VecsMult(PetscScalar *M,MatType_t sM,PetscInt ldM,Vec *U,PetscInt sU,PetscInt eU,Vec *V,PetscInt sV,PetscInt eV,PetscScalar *workS0,PetscScalar *workS1);
448: PETSC_INTERN PetscErrorCode VecsMultS(PetscScalar *M,MatType_t sM,PetscInt ldM,Vec *U,PetscInt sU,PetscInt eU,Vec *V,PetscInt sV,PetscInt eV,DvdReduction *r,DvdMult_copy_func *sr);
449: PETSC_INTERN PetscErrorCode VecsMultIb(PetscScalar *M,MatType_t sM,PetscInt ldM,PetscInt rM,PetscInt cM,PetscScalar *auxS,Vec V);
450: PETSC_INTERN PetscErrorCode VecsMultIa(PetscScalar *M,MatType_t sM,PetscInt ldM,Vec *U,PetscInt sU,PetscInt eU,Vec *V,PetscInt sV,PetscInt eV);
451: PETSC_INTERN PetscErrorCode SlepcAllReduceSumBegin(DvdReductionChunk *ops,PetscInt max_size_ops,PetscScalar *in,PetscScalar *out,PetscInt max_size_in,DvdReduction *r,MPI_Comm comm);
452: PETSC_INTERN PetscErrorCode SlepcAllReduceSum(DvdReduction *r,PetscInt size_in,DvdReductionPostF f,void *ptr,PetscScalar **in);
453: PETSC_INTERN PetscErrorCode SlepcAllReduceSumEnd(DvdReduction *r);
454: PETSC_INTERN PetscErrorCode dvd_orthV(IP ip,Vec *DS,PetscInt size_DS,Vec *cX,PetscInt size_cX,Vec *V,PetscInt V_new_s,PetscInt V_new_e,PetscScalar *auxS,PetscRandom rand);
455: PETSC_INTERN PetscErrorCode dvd_BorthV_faster(IP ip,Vec *DS,Vec *BDS,PetscReal *BDSn,PetscInt size_DS,Vec *cX,Vec *BcX,PetscReal *BcXn,PetscInt size_cX,Vec *V,Vec *BV,PetscReal *BVn,PetscInt V_new_s,PetscInt V_new_e,PetscScalar *auxS,PetscRandom rand);
456: PETSC_INTERN 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);
458: /* SLEPc interface routines */
459: PETSC_INTERN PetscErrorCode SLEPcNotImplemented();
460: PETSC_INTERN PetscErrorCode EPSCreate_XD(EPS eps);
461: PETSC_INTERN PetscErrorCode EPSReset_XD(EPS eps);
462: PETSC_INTERN PetscErrorCode EPSSetUp_XD(EPS eps);
463: PETSC_INTERN PetscErrorCode EPSSolve_XD(EPS eps);
464: PETSC_INTERN PetscErrorCode EPSComputeVectors_XD(EPS eps);
465: PETSC_INTERN PetscErrorCode EPSXDSetKrylovStart_XD(EPS eps,PetscBool krylovstart);
466: PETSC_INTERN PetscErrorCode EPSXDGetKrylovStart_XD(EPS eps,PetscBool *krylovstart);
467: PETSC_INTERN PetscErrorCode EPSXDSetBlockSize_XD(EPS eps,PetscInt blocksize);
468: PETSC_INTERN PetscErrorCode EPSXDGetBlockSize_XD(EPS eps,PetscInt *blocksize);
469: PETSC_INTERN PetscErrorCode EPSXDSetRestart_XD(EPS eps,PetscInt minv,PetscInt plusk);
470: PETSC_INTERN PetscErrorCode EPSXDGetRestart_XD(EPS eps,PetscInt *minv,PetscInt *plusk);
471: PETSC_INTERN PetscErrorCode EPSXDGetInitialSize_XD(EPS eps,PetscInt *initialsize);
472: PETSC_INTERN PetscErrorCode EPSXDSetInitialSize_XD(EPS eps,PetscInt initialsize);
473: PETSC_INTERN PetscErrorCode EPSXDGetFix_XD(EPS eps,PetscReal *fix);
474: PETSC_INTERN PetscErrorCode EPSJDSetFix_JD(EPS eps,PetscReal fix);
475: PETSC_INTERN PetscErrorCode EPSXDSetBOrth_XD(EPS eps,EPSOrthType borth);
476: PETSC_INTERN PetscErrorCode EPSXDGetBOrth_XD(EPS eps,EPSOrthType *borth);
477: PETSC_INTERN PetscErrorCode EPSJDSetConstCorrectionTol_JD(EPS eps,PetscBool constant);
478: PETSC_INTERN PetscErrorCode EPSJDGetConstCorrectionTol_JD(EPS eps,PetscBool *constant);
479: PETSC_INTERN PetscErrorCode EPSXDSetWindowSizes_XD(EPS eps,PetscInt pwindow,PetscInt qwindow);
480: PETSC_INTERN PetscErrorCode EPSXDGetWindowSizes_XD(EPS eps,PetscInt *pwindow,PetscInt *qwindow);
481: PETSC_INTERN PetscErrorCode EPSXDSetMethod(EPS eps,Method_t method);
482: PETSC_INTERN PetscErrorCode EPSXDGetMethod_XD(EPS eps,Method_t *method);
484: /* Common inline function */
487: PETSC_STATIC_INLINE PetscErrorCode dvd_improvex_compute_X(dvdDashboard *d,PetscInt i_s,PetscInt i_e,Vec *u,PetscScalar *pX,PetscInt ld)
488: {
489: PetscErrorCode ierr;
490: PetscInt n = i_e - i_s,i;
493: 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);
494: /* nX(i) <- ||X(i)|| */
495: if (d->correctXnorm) {
496: for (i=0; i<n; i++) {
497: VecNormBegin(u[i],NORM_2,&d->nX[i_s+i]);
498: }
499: for (i=0; i<n; i++) {
500: VecNormEnd(u[i],NORM_2,&d->nX[i_s+i]);
501: }
502: #if !defined(PETSC_USE_COMPLEX)
503: for (i=0;i<n;i++) {
504: if (d->eigi[i_s+i] != 0.0) {
505: 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]);
506: i++;
507: }
508: }
509: #endif
510: } else {
511: for (i=0; i<n; i++) d->nX[i_s+i] = 1.0;
512: }
513: return(0);
514: }
517: #define _Ceil(A,B) ((A)/(B)+((A)%(B)==0?0:1))
518: #define FromIntToScalar(S) ((PetscInt)_Ceil((S)*sizeof(PetscBLASInt),sizeof(PetscScalar)))
519: #define FromRealToScalar(S) ((PetscInt)_Ceil((S)*sizeof(PetscReal),sizeof(PetscScalar)))