| 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 | #pragma once | ||
| 12 | |||
| 13 | #include <slepcpep.h> | ||
| 14 | #include <slepc/private/slepcimpl.h> | ||
| 15 | |||
| 16 | /* SUBMANSEC = PEP */ | ||
| 17 | |||
| 18 | SLEPC_EXTERN PetscBool PEPRegisterAllCalled; | ||
| 19 | SLEPC_EXTERN PetscBool PEPMonitorRegisterAllCalled; | ||
| 20 | SLEPC_EXTERN PetscErrorCode PEPRegisterAll(void); | ||
| 21 | SLEPC_EXTERN PetscErrorCode PEPMonitorRegisterAll(void); | ||
| 22 | SLEPC_EXTERN PetscLogEvent PEP_SetUp,PEP_Solve,PEP_Refine,PEP_CISS_SVD; | ||
| 23 | |||
| 24 | typedef struct _PEPOps *PEPOps; | ||
| 25 | |||
| 26 | struct _PEPOps { | ||
| 27 | PetscErrorCode (*solve)(PEP); | ||
| 28 | PetscErrorCode (*setup)(PEP); | ||
| 29 | PetscErrorCode (*setfromoptions)(PEP,PetscOptionItems); | ||
| 30 | PetscErrorCode (*publishoptions)(PEP); | ||
| 31 | PetscErrorCode (*destroy)(PEP); | ||
| 32 | PetscErrorCode (*reset)(PEP); | ||
| 33 | PetscErrorCode (*view)(PEP,PetscViewer); | ||
| 34 | PetscErrorCode (*backtransform)(PEP); | ||
| 35 | PetscErrorCode (*computevectors)(PEP); | ||
| 36 | PetscErrorCode (*extractvectors)(PEP); | ||
| 37 | PetscErrorCode (*setdefaultst)(PEP); | ||
| 38 | PetscErrorCode (*setdstype)(PEP); | ||
| 39 | }; | ||
| 40 | |||
| 41 | /* | ||
| 42 | Maximum number of monitors you can run with a single PEP | ||
| 43 | */ | ||
| 44 | #define MAXPEPMONITORS 5 | ||
| 45 | |||
| 46 | typedef enum { PEP_STATE_INITIAL, | ||
| 47 | PEP_STATE_SETUP, | ||
| 48 | PEP_STATE_SOLVED, | ||
| 49 | PEP_STATE_EIGENVECTORS } PEPStateType; | ||
| 50 | |||
| 51 | /* | ||
| 52 | To check for unsupported features at PEPSetUp_XXX() | ||
| 53 | */ | ||
| 54 | typedef enum { PEP_FEATURE_NONMONOMIAL=1, /* non-monomial bases */ | ||
| 55 | PEP_FEATURE_REGION=4, /* nontrivial region for filtering */ | ||
| 56 | PEP_FEATURE_EXTRACT=8, /* eigenvector extraction */ | ||
| 57 | PEP_FEATURE_CONVERGENCE=16, /* convergence test selected by user */ | ||
| 58 | PEP_FEATURE_STOPPING=32, /* stopping test */ | ||
| 59 | PEP_FEATURE_SCALE=64 /* scaling */ | ||
| 60 | } PEPFeatureType; | ||
| 61 | |||
| 62 | /* | ||
| 63 | Defines the PEP data structure. | ||
| 64 | */ | ||
| 65 | struct _p_PEP { | ||
| 66 | PETSCHEADER(struct _PEPOps); | ||
| 67 | /*------------------------- User parameters ---------------------------*/ | ||
| 68 | PetscInt max_it; /* maximum number of iterations */ | ||
| 69 | PetscInt nev; /* number of eigenvalues to compute */ | ||
| 70 | PetscInt ncv; /* number of basis vectors */ | ||
| 71 | PetscInt mpd; /* maximum dimension of projected problem */ | ||
| 72 | PetscInt nini; /* number of initial vectors (negative means not copied yet) */ | ||
| 73 | PetscScalar target; /* target value */ | ||
| 74 | PetscReal tol; /* tolerance */ | ||
| 75 | PEPConv conv; /* convergence test */ | ||
| 76 | PEPStop stop; /* stopping test */ | ||
| 77 | PEPWhich which; /* which part of the spectrum to be sought */ | ||
| 78 | PetscReal inta,intb; /* interval [a,b] for spectrum slicing */ | ||
| 79 | PEPBasis basis; /* polynomial basis used to represent the problem */ | ||
| 80 | PEPProblemType problem_type; /* which kind of problem to be solved */ | ||
| 81 | PEPScale scale; /* scaling strategy to be used */ | ||
| 82 | PetscReal sfactor,dsfactor; /* scaling factors */ | ||
| 83 | PetscInt sits; /* number of iterations of the scaling method */ | ||
| 84 | PetscReal slambda; /* norm eigenvalue approximation for scaling */ | ||
| 85 | PEPRefine refine; /* type of refinement to be applied after solve */ | ||
| 86 | PetscInt npart; /* number of partitions of the communicator */ | ||
| 87 | PetscReal rtol; /* tolerance for refinement */ | ||
| 88 | PetscInt rits; /* number of iterations of the refinement method */ | ||
| 89 | PEPRefineScheme scheme; /* scheme for solving linear systems within refinement */ | ||
| 90 | PEPExtract extract; /* type of extraction used */ | ||
| 91 | PetscBool trackall; /* whether all the residuals must be computed */ | ||
| 92 | |||
| 93 | /*-------------- User-provided functions and contexts -----------------*/ | ||
| 94 | PEPConvergenceTestFn *converged; | ||
| 95 | PEPConvergenceTestFn *convergeduser; | ||
| 96 | PetscCtxDestroyFn *convergeddestroy; | ||
| 97 | PEPStoppingTestFn *stopping; | ||
| 98 | PEPStoppingTestFn *stoppinguser; | ||
| 99 | PetscCtxDestroyFn *stoppingdestroy; | ||
| 100 | void *convergedctx; | ||
| 101 | void *stoppingctx; | ||
| 102 | PEPMonitorFn *monitor[MAXPEPMONITORS]; | ||
| 103 | PetscCtxDestroyFn *monitordestroy[MAXPEPMONITORS]; | ||
| 104 | void *monitorcontext[MAXPEPMONITORS]; | ||
| 105 | PetscInt numbermonitors; | ||
| 106 | |||
| 107 | /*----------------- Child objects and working data -------------------*/ | ||
| 108 | ST st; /* spectral transformation object */ | ||
| 109 | DS ds; /* direct solver object */ | ||
| 110 | BV V; /* set of basis vectors and computed eigenvectors */ | ||
| 111 | RG rg; /* optional region for filtering */ | ||
| 112 | SlepcSC sc; /* sorting criterion data */ | ||
| 113 | Mat *A; /* coefficient matrices of the polynomial */ | ||
| 114 | PetscInt nmat; /* number of matrices */ | ||
| 115 | Vec Dl,Dr; /* diagonal matrices for balancing */ | ||
| 116 | Vec *IS; /* references to user-provided initial space */ | ||
| 117 | PetscScalar *eigr,*eigi; /* real and imaginary parts of eigenvalues */ | ||
| 118 | PetscReal *errest; /* error estimates */ | ||
| 119 | PetscInt *perm; /* permutation for eigenvalue ordering */ | ||
| 120 | PetscReal *pbc; /* coefficients defining the polynomial basis */ | ||
| 121 | PetscScalar *solvematcoeffs; /* coefficients to compute the matrix to be inverted */ | ||
| 122 | PetscInt nwork; /* number of work vectors */ | ||
| 123 | Vec *work; /* work vectors */ | ||
| 124 | KSP refineksp; /* ksp used in refinement */ | ||
| 125 | PetscSubcomm refinesubc; /* context for sub-communicators */ | ||
| 126 | void *data; /* placeholder for solver-specific stuff */ | ||
| 127 | |||
| 128 | /* ----------------------- Status variables --------------------------*/ | ||
| 129 | PEPStateType state; /* initial -> setup -> solved -> eigenvectors */ | ||
| 130 | PetscInt nconv; /* number of converged eigenvalues */ | ||
| 131 | PetscInt its; /* number of iterations so far computed */ | ||
| 132 | PetscInt n,nloc; /* problem dimensions (global, local) */ | ||
| 133 | PetscReal *nrma; /* computed matrix norms */ | ||
| 134 | PetscReal nrml[2]; /* computed matrix norms for the linearization */ | ||
| 135 | PetscBool sfactor_set; /* flag to indicate the user gave sfactor */ | ||
| 136 | PetscBool lineariz; /* current solver is based on linearization */ | ||
| 137 | PEPConvergedReason reason; | ||
| 138 | }; | ||
| 139 | |||
| 140 | /* | ||
| 141 | Macros to test valid PEP arguments | ||
| 142 | */ | ||
| 143 | #if !defined(PETSC_USE_DEBUG) | ||
| 144 | |||
| 145 | #define PEPCheckSolved(h,arg) do {(void)(h);} while (0) | ||
| 146 | |||
| 147 | #else | ||
| 148 | |||
| 149 | #define PEPCheckSolved(h,arg) \ | ||
| 150 | do { \ | ||
| 151 | PetscCheck((h)->state>=PEP_STATE_SOLVED,PetscObjectComm((PetscObject)(h)),PETSC_ERR_ARG_WRONGSTATE,"Must call PEPSolve() first: Parameter #%d",arg); \ | ||
| 152 | } while (0) | ||
| 153 | |||
| 154 | #endif | ||
| 155 | |||
| 156 | /* | ||
| 157 | Macros to check settings at PEPSetUp() | ||
| 158 | */ | ||
| 159 | |||
| 160 | /* PEPCheckHermitian: the problem is Hermitian or Hyperbolic */ | ||
| 161 | #define PEPCheckHermitianCondition(pep,condition,msg) \ | ||
| 162 | do { \ | ||
| 163 | if (condition) { \ | ||
| 164 | PetscCheck((pep)->problem_type==PEP_HERMITIAN || (pep)->problem_type==PEP_HYPERBOLIC,PetscObjectComm((PetscObject)(pep)),PETSC_ERR_SUP,"The solver '%s'%s can only be used for Hermitian (or hyperbolic) problems",((PetscObject)(pep))->type_name,(msg)); \ | ||
| 165 | } \ | ||
| 166 | } while (0) | ||
| 167 | #define PEPCheckHermitian(pep) PEPCheckHermitianCondition(pep,PETSC_TRUE,"") | ||
| 168 | |||
| 169 | /* PEPCheckQuadratic: the polynomial has degree 2 */ | ||
| 170 | #define PEPCheckQuadraticCondition(pep,condition,msg) \ | ||
| 171 | do { \ | ||
| 172 | if (condition) { \ | ||
| 173 | PetscCheck((pep)->nmat==3,PetscObjectComm((PetscObject)(pep)),PETSC_ERR_SUP,"The solver '%s'%s is only available for quadratic problems",((PetscObject)(pep))->type_name,(msg)); \ | ||
| 174 | } \ | ||
| 175 | } while (0) | ||
| 176 | #define PEPCheckQuadratic(pep) PEPCheckQuadraticCondition(pep,PETSC_TRUE,"") | ||
| 177 | |||
| 178 | /* PEPCheckShiftSinvert: shift or shift-and-invert ST */ | ||
| 179 | #define PEPCheckShiftSinvertCondition(pep,condition,msg) \ | ||
| 180 | do { \ | ||
| 181 | if (condition) { \ | ||
| 182 | PetscBool __flg; \ | ||
| 183 | PetscCall(PetscObjectTypeCompareAny((PetscObject)(pep)->st,&__flg,STSINVERT,STSHIFT,"")); \ | ||
| 184 | PetscCheck(__flg,PetscObjectComm((PetscObject)(pep)),PETSC_ERR_SUP,"The solver '%s'%s requires shift or shift-and-invert spectral transform",((PetscObject)(pep))->type_name,(msg)); \ | ||
| 185 | } \ | ||
| 186 | } while (0) | ||
| 187 | #define PEPCheckShiftSinvert(pep) PEPCheckShiftSinvertCondition(pep,PETSC_TRUE,"") | ||
| 188 | |||
| 189 | /* PEPCheckSinvertCayley: shift-and-invert or Cayley ST */ | ||
| 190 | #define PEPCheckSinvertCayleyCondition(pep,condition,msg) \ | ||
| 191 | do { \ | ||
| 192 | if (condition) { \ | ||
| 193 | PetscBool __flg; \ | ||
| 194 | PetscCall(PetscObjectTypeCompareAny((PetscObject)(pep)->st,&__flg,STSINVERT,STCAYLEY,"")); \ | ||
| 195 | PetscCheck(__flg,PetscObjectComm((PetscObject)(pep)),PETSC_ERR_SUP,"The solver '%s'%s requires shift-and-invert or Cayley transform",((PetscObject)(pep))->type_name,(msg)); \ | ||
| 196 | } \ | ||
| 197 | } while (0) | ||
| 198 | #define PEPCheckSinvertCayley(pep) PEPCheckSinvertCayleyCondition(pep,PETSC_TRUE,"") | ||
| 199 | |||
| 200 | /* Check for unsupported features */ | ||
| 201 | #define PEPCheckUnsupportedCondition(pep,mask,condition,msg) \ | ||
| 202 | do { \ | ||
| 203 | if (condition) { \ | ||
| 204 | PetscCheck(!((mask) & PEP_FEATURE_NONMONOMIAL) || (pep)->basis==PEP_BASIS_MONOMIAL,PetscObjectComm((PetscObject)(pep)),PETSC_ERR_SUP,"The solver '%s'%s is not implemented for non-monomial bases",((PetscObject)(pep))->type_name,(msg)); \ | ||
| 205 | if ((mask) & PEP_FEATURE_REGION) { \ | ||
| 206 | PetscBool __istrivial; \ | ||
| 207 | PetscCall(RGIsTrivial((pep)->rg,&__istrivial)); \ | ||
| 208 | PetscCheck(__istrivial,PetscObjectComm((PetscObject)(pep)),PETSC_ERR_SUP,"The solver '%s'%s does not support region filtering",((PetscObject)(pep))->type_name,(msg)); \ | ||
| 209 | } \ | ||
| 210 | PetscCheck(!((mask) & PEP_FEATURE_EXTRACT) || !(pep)->extract || (pep)->extract==PEP_EXTRACT_NONE,PetscObjectComm((PetscObject)(pep)),PETSC_ERR_SUP,"The solver '%s'%s does not support extraction variants",((PetscObject)(pep))->type_name,(msg)); \ | ||
| 211 | PetscCheck(!((mask) & PEP_FEATURE_CONVERGENCE) || (pep)->converged==PEPConvergedRelative,PetscObjectComm((PetscObject)(pep)),PETSC_ERR_SUP,"The solver '%s'%s only supports the default convergence test",((PetscObject)(pep))->type_name,(msg)); \ | ||
| 212 | PetscCheck(!((mask) & PEP_FEATURE_STOPPING) || (pep)->stopping==PEPStoppingBasic,PetscObjectComm((PetscObject)(pep)),PETSC_ERR_SUP,"The solver '%s'%s only supports the default stopping test",((PetscObject)(pep))->type_name,(msg)); \ | ||
| 213 | } \ | ||
| 214 | } while (0) | ||
| 215 | #define PEPCheckUnsupported(pep,mask) PEPCheckUnsupportedCondition(pep,mask,PETSC_TRUE,"") | ||
| 216 | |||
| 217 | /* Check for ignored features */ | ||
| 218 | #define PEPCheckIgnoredCondition(pep,mask,condition,msg) \ | ||
| 219 | do { \ | ||
| 220 | if (condition) { \ | ||
| 221 | if (((mask) & PEP_FEATURE_NONMONOMIAL) && (pep)->basis!=PEP_BASIS_MONOMIAL) PetscCall(PetscInfo((pep),"The solver '%s'%s ignores the basis settings\n",((PetscObject)(pep))->type_name,(msg))); \ | ||
| 222 | if ((mask) & PEP_FEATURE_REGION) { \ | ||
| 223 | PetscBool __istrivial; \ | ||
| 224 | PetscCall(RGIsTrivial((pep)->rg,&__istrivial)); \ | ||
| 225 | if (!__istrivial) PetscCall(PetscInfo((pep),"The solver '%s'%s ignores the specified region\n",((PetscObject)(pep))->type_name,(msg))); \ | ||
| 226 | } \ | ||
| 227 | if (((mask) & PEP_FEATURE_EXTRACT) && (pep)->extract && (pep)->extract!=PEP_EXTRACT_NONE) PetscCall(PetscInfo((pep),"The solver '%s'%s ignores the extract settings\n",((PetscObject)(pep))->type_name,(msg))); \ | ||
| 228 | if (((mask) & PEP_FEATURE_CONVERGENCE) && (pep)->converged!=PEPConvergedRelative) PetscCall(PetscInfo((pep),"The solver '%s'%s ignores the convergence test settings\n",((PetscObject)(pep))->type_name,(msg))); \ | ||
| 229 | if (((mask) & PEP_FEATURE_STOPPING) && (pep)->stopping!=PEPStoppingBasic) PetscCall(PetscInfo((pep),"The solver '%s'%s ignores the stopping test settings\n",((PetscObject)(pep))->type_name,(msg))); \ | ||
| 230 | if (((mask) & PEP_FEATURE_SCALE) && (pep)->scale!=PEP_SCALE_NONE) PetscCall(PetscInfo((pep),"The solver '%s'%s ignores the scaling settings\n",((PetscObject)(pep))->type_name,(msg))); \ | ||
| 231 | } \ | ||
| 232 | } while (0) | ||
| 233 | #define PEPCheckIgnored(pep,mask) PEPCheckIgnoredCondition(pep,mask,PETSC_TRUE,"") | ||
| 234 | |||
| 235 | /* | ||
| 236 | PEP_KSPSetOperators - Sets the KSP matrices | ||
| 237 | */ | ||
| 238 | 2920 | static inline PetscErrorCode PEP_KSPSetOperators(KSP ksp,Mat A,Mat B) | |
| 239 | { | ||
| 240 | 2920 | const char *prefix; | |
| 241 | |||
| 242 |
1/2✓ Branch 0 taken 7 times.
✗ Branch 1 not taken.
|
2920 | PetscFunctionBegin; |
| 243 |
4/6✓ Branch 0 taken 7 times.
✓ Branch 1 taken 28 times.
✓ Branch 2 taken 7 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 7 times.
|
2920 | PetscCall(KSPSetOperators(ksp,A,B)); |
| 244 |
4/6✓ Branch 0 taken 7 times.
✓ Branch 1 taken 28 times.
✓ Branch 2 taken 7 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 7 times.
|
2920 | PetscCall(MatGetOptionsPrefix(B,&prefix)); |
| 245 |
2/2✓ Branch 0 taken 35 times.
✓ Branch 1 taken 20 times.
|
2920 | if (!prefix) { |
| 246 | /* set Mat prefix to be the same as KSP to enable setting command-line options (e.g. MUMPS) | ||
| 247 | only applies if the Mat has no user-defined prefix */ | ||
| 248 |
4/6✓ Branch 0 taken 7 times.
✓ Branch 1 taken 28 times.
✓ Branch 2 taken 7 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 7 times.
|
2261 | PetscCall(KSPGetOptionsPrefix(ksp,&prefix)); |
| 249 |
4/6✓ Branch 0 taken 7 times.
✓ Branch 1 taken 28 times.
✓ Branch 2 taken 7 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 7 times.
|
2261 | PetscCall(MatSetOptionsPrefix(B,prefix)); |
| 250 | } | ||
| 251 |
6/12✓ Branch 0 taken 7 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 7 times.
✓ Branch 4 taken 7 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 7 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 7 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 7 times.
|
580 | PetscFunctionReturn(PETSC_SUCCESS); |
| 252 | } | ||
| 253 | |||
| 254 | SLEPC_INTERN PetscErrorCode PEPSetWhichEigenpairs_Default(PEP); | ||
| 255 | SLEPC_INTERN PetscErrorCode PEPSetDimensions_Default(PEP,PetscInt,PetscInt*,PetscInt*); | ||
| 256 | SLEPC_INTERN PetscErrorCode PEPExtractVectors(PEP); | ||
| 257 | SLEPC_INTERN PetscErrorCode PEPBackTransform_Default(PEP); | ||
| 258 | SLEPC_INTERN PetscErrorCode PEPComputeVectors(PEP); | ||
| 259 | SLEPC_INTERN PetscErrorCode PEPComputeVectors_Default(PEP); | ||
| 260 | SLEPC_INTERN PetscErrorCode PEPComputeVectors_Indefinite(PEP); | ||
| 261 | SLEPC_INTERN PetscErrorCode PEPComputeResidualNorm_Private(PEP,PetscScalar,PetscScalar,Vec,Vec,Vec*,PetscReal*); | ||
| 262 | SLEPC_INTERN PetscErrorCode PEPKrylovConvergence(PEP,PetscBool,PetscInt,PetscInt,PetscReal,PetscInt*); | ||
| 263 | SLEPC_INTERN PetscErrorCode PEPComputeScaleFactor(PEP); | ||
| 264 | SLEPC_INTERN PetscErrorCode PEPBuildDiagonalScaling(PEP); | ||
| 265 | SLEPC_INTERN PetscErrorCode PEPBasisCoefficients(PEP,PetscReal*); | ||
| 266 | SLEPC_INTERN PetscErrorCode PEPEvaluateBasis(PEP,PetscScalar,PetscScalar,PetscScalar*,PetscScalar*); | ||
| 267 | SLEPC_INTERN PetscErrorCode PEPEvaluateBasisDerivative(PEP,PetscScalar,PetscScalar,PetscScalar*,PetscScalar*); | ||
| 268 | SLEPC_INTERN PetscErrorCode PEPEvaluateBasisMat(PEP,PetscInt,PetscScalar*,PetscInt,PetscInt,PetscScalar*,PetscInt,PetscScalar*,PetscInt,PetscScalar*,PetscInt); | ||
| 269 | SLEPC_INTERN PetscErrorCode PEPNewtonRefinement_TOAR(PEP,PetscScalar,PetscInt*,PetscReal*,PetscInt,PetscScalar*,PetscInt); | ||
| 270 | SLEPC_INTERN PetscErrorCode PEPNewtonRefinementSimple(PEP,PetscInt*,PetscReal,PetscInt); | ||
| 271 | SLEPC_INTERN PetscErrorCode PEPSetDefaultST(PEP); | ||
| 272 | SLEPC_INTERN PetscErrorCode PEPSetDefaultST_Transform(PEP); | ||
| 273 |