GCC Code Coverage Report

Directory: ./
File: include/slepc/private/pepimpl.h
Date: 2025-10-04 04:19:13
Exec Total Coverage
Lines: 9 9 100.0%
Functions: 1 1 100.0%
Branches: 25 40 62.5%
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