| Line | Branch | Exec | Source |
|---|---|---|---|
| 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 | 9686 | static inline PetscErrorCode EPSDavidsonFLAdd(dvdFunctionList **fl,dvdCallback f) | |
| 190 | { | ||
| 191 | 9686 | dvdFunctionList *l; | |
| 192 | |||
| 193 |
1/2✓ Branch 0 taken 12 times.
✗ Branch 1 not taken.
|
9686 | PetscFunctionBegin; |
| 194 |
4/6✓ Branch 0 taken 12 times.
✓ Branch 1 taken 48 times.
✓ Branch 2 taken 12 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 12 times.
|
9686 | PetscCall(PetscNew(&l)); |
| 195 | 9686 | l->f = f; | |
| 196 | 9686 | l->next = *fl; | |
| 197 | 9686 | *fl = l; | |
| 198 |
6/12✓ Branch 0 taken 12 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 12 times.
✓ Branch 4 taken 12 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 12 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 12 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 12 times.
|
9686 | PetscFunctionReturn(PETSC_SUCCESS); |
| 199 | } | ||
| 200 | |||
| 201 | 3654 | static inline PetscErrorCode EPSDavidsonFLCall(dvdFunctionList *fl,dvdDashboard *d) | |
| 202 | { | ||
| 203 | 3654 | dvdFunctionList *l; | |
| 204 | |||
| 205 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
3654 | PetscFunctionBegin; |
| 206 |
7/8✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 10 times.
✓ Branch 3 taken 8 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✓ Branch 6 taken 2 times.
✓ Branch 7 taken 2 times.
|
13340 | for (l=fl;l;l=l->next) PetscCall(l->f(d)); |
| 207 |
6/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
718 | PetscFunctionReturn(PETSC_SUCCESS); |
| 208 | } | ||
| 209 | |||
| 210 | 5082 | static inline PetscErrorCode EPSDavidsonFLDestroy(dvdFunctionList **fl) | |
| 211 | { | ||
| 212 | 5082 | dvdFunctionList *l,*l0; | |
| 213 | |||
| 214 |
1/2✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
|
5082 | PetscFunctionBegin; |
| 215 |
2/2✓ Branch 0 taken 10 times.
✓ Branch 1 taken 10 times.
|
14768 | for (l=*fl;l;l=l0) { |
| 216 | 9686 | l0 = l->next; | |
| 217 |
6/8✓ Branch 0 taken 2 times.
✓ Branch 1 taken 8 times.
✓ Branch 2 taken 2 times.
✓ Branch 3 taken 8 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 2 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 2 times.
|
9686 | PetscCall(PetscFree(l)); |
| 218 | } | ||
| 219 | 5082 | *fl = NULL; | |
| 220 |
6/12✓ Branch 0 taken 2 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 2 times.
|
5082 | 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*); | ||
| 295 |