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 <slepcnep.h> | ||
14 | #include <slepc/private/slepcimpl.h> | ||
15 | |||
16 | /* SUBMANSEC = NEP */ | ||
17 | |||
18 | SLEPC_EXTERN PetscBool NEPRegisterAllCalled; | ||
19 | SLEPC_EXTERN PetscBool NEPMonitorRegisterAllCalled; | ||
20 | SLEPC_EXTERN PetscErrorCode NEPRegisterAll(void); | ||
21 | SLEPC_EXTERN PetscErrorCode NEPMonitorRegisterAll(void); | ||
22 | SLEPC_EXTERN PetscLogEvent NEP_SetUp,NEP_Solve,NEP_Refine,NEP_FunctionEval,NEP_JacobianEval,NEP_Resolvent,NEP_CISS_SVD; | ||
23 | |||
24 | typedef struct _NEPOps *NEPOps; | ||
25 | |||
26 | struct _NEPOps { | ||
27 | PetscErrorCode (*solve)(NEP); | ||
28 | PetscErrorCode (*setup)(NEP); | ||
29 | PetscErrorCode (*setfromoptions)(NEP,PetscOptionItems); | ||
30 | PetscErrorCode (*publishoptions)(NEP); | ||
31 | PetscErrorCode (*destroy)(NEP); | ||
32 | PetscErrorCode (*reset)(NEP); | ||
33 | PetscErrorCode (*view)(NEP,PetscViewer); | ||
34 | PetscErrorCode (*computevectors)(NEP); | ||
35 | PetscErrorCode (*setdstype)(NEP); | ||
36 | }; | ||
37 | |||
38 | /* | ||
39 | Maximum number of monitors you can run with a single NEP | ||
40 | */ | ||
41 | #define MAXNEPMONITORS 5 | ||
42 | |||
43 | typedef enum { NEP_STATE_INITIAL, | ||
44 | NEP_STATE_SETUP, | ||
45 | NEP_STATE_SOLVED, | ||
46 | NEP_STATE_EIGENVECTORS } NEPStateType; | ||
47 | |||
48 | /* | ||
49 | How the problem function T(lambda) has been defined by the user | ||
50 | - Callback: one callback to build the function matrix, another one for the Jacobian | ||
51 | - Split: in split form sum_j(A_j*f_j(lambda)) | ||
52 | */ | ||
53 | typedef enum { NEP_USER_INTERFACE_CALLBACK=1, | ||
54 | NEP_USER_INTERFACE_SPLIT } NEPUserInterface; | ||
55 | |||
56 | /* | ||
57 | To check for unsupported features at NEPSetUp_XXX() | ||
58 | */ | ||
59 | typedef enum { NEP_FEATURE_CALLBACK=1, /* callback user interface */ | ||
60 | NEP_FEATURE_REGION=4, /* nontrivial region for filtering */ | ||
61 | NEP_FEATURE_CONVERGENCE=16, /* convergence test selected by user */ | ||
62 | NEP_FEATURE_STOPPING=32, /* stopping test */ | ||
63 | NEP_FEATURE_TWOSIDED=64 /* two-sided variant */ | ||
64 | } NEPFeatureType; | ||
65 | |||
66 | /* | ||
67 | Defines the NEP data structure. | ||
68 | */ | ||
69 | struct _p_NEP { | ||
70 | PETSCHEADER(struct _NEPOps); | ||
71 | /*------------------------- User parameters ---------------------------*/ | ||
72 | PetscInt max_it; /* maximum number of iterations */ | ||
73 | PetscInt nev; /* number of eigenvalues to compute */ | ||
74 | PetscInt ncv; /* number of basis vectors */ | ||
75 | PetscInt mpd; /* maximum dimension of projected problem */ | ||
76 | PetscInt nini; /* number of initial vectors (negative means not copied yet) */ | ||
77 | PetscScalar target; /* target value */ | ||
78 | PetscReal tol; /* tolerance */ | ||
79 | NEPConv conv; /* convergence test */ | ||
80 | NEPStop stop; /* stopping test */ | ||
81 | NEPWhich which; /* which part of the spectrum to be sought */ | ||
82 | NEPProblemType problem_type; /* which kind of problem to be solved */ | ||
83 | NEPRefine refine; /* type of refinement to be applied after solve */ | ||
84 | PetscInt npart; /* number of partitions of the communicator */ | ||
85 | PetscReal rtol; /* tolerance for refinement */ | ||
86 | PetscInt rits; /* number of iterations of the refinement method */ | ||
87 | NEPRefineScheme scheme; /* scheme for solving linear systems within refinement */ | ||
88 | PetscBool trackall; /* whether all the residuals must be computed */ | ||
89 | PetscBool twosided; /* whether to compute left eigenvectors (two-sided solver) */ | ||
90 | |||
91 | /*-------------- User-provided functions and contexts -----------------*/ | ||
92 | NEPFunctionFn *computefunction; | ||
93 | NEPJacobianFn *computejacobian; | ||
94 | void *functionctx; | ||
95 | void *jacobianctx; | ||
96 | NEPConvergenceTestFn *converged; | ||
97 | NEPConvergenceTestFn *convergeduser; | ||
98 | PetscCtxDestroyFn *convergeddestroy; | ||
99 | NEPStoppingTestFn *stopping; | ||
100 | NEPStoppingTestFn *stoppinguser; | ||
101 | PetscCtxDestroyFn *stoppingdestroy; | ||
102 | void *convergedctx; | ||
103 | void *stoppingctx; | ||
104 | NEPMonitorFn *monitor[MAXNEPMONITORS]; | ||
105 | PetscCtxDestroyFn *monitordestroy[MAXNEPMONITORS]; | ||
106 | void *monitorcontext[MAXNEPMONITORS]; | ||
107 | PetscInt numbermonitors; | ||
108 | |||
109 | /*----------------- Child objects and working data -------------------*/ | ||
110 | DS ds; /* direct solver object */ | ||
111 | BV V; /* set of basis vectors and computed eigenvectors */ | ||
112 | BV W; /* left basis vectors (if left eigenvectors requested) */ | ||
113 | RG rg; /* optional region for filtering */ | ||
114 | SlepcSC sc; /* sorting criterion data */ | ||
115 | Mat function; /* function matrix */ | ||
116 | Mat function_pre; /* function matrix (preconditioner) */ | ||
117 | Mat jacobian; /* Jacobian matrix */ | ||
118 | Mat *A; /* matrix coefficients of split form */ | ||
119 | FN *f; /* matrix functions of split form */ | ||
120 | PetscInt nt; /* number of terms in split form */ | ||
121 | MatStructure mstr; /* pattern of split matrices */ | ||
122 | Mat *P; /* matrix coefficients of split form (preconditioner) */ | ||
123 | MatStructure mstrp; /* pattern of split matrices (preconditioner) */ | ||
124 | Vec *IS; /* references to user-provided initial space */ | ||
125 | PetscScalar *eigr,*eigi; /* real and imaginary parts of eigenvalues */ | ||
126 | PetscReal *errest; /* error estimates */ | ||
127 | PetscInt *perm; /* permutation for eigenvalue ordering */ | ||
128 | PetscInt nwork; /* number of work vectors */ | ||
129 | Vec *work; /* work vectors */ | ||
130 | KSP refineksp; /* ksp used in refinement */ | ||
131 | PetscSubcomm refinesubc; /* context for sub-communicators */ | ||
132 | void *data; /* placeholder for solver-specific stuff */ | ||
133 | |||
134 | /* ----------------------- Status variables --------------------------*/ | ||
135 | NEPStateType state; /* initial -> setup -> solved -> eigenvectors */ | ||
136 | PetscInt nconv; /* number of converged eigenvalues */ | ||
137 | PetscInt its; /* number of iterations so far computed */ | ||
138 | PetscInt n,nloc; /* problem dimensions (global, local) */ | ||
139 | PetscReal *nrma; /* computed matrix norms */ | ||
140 | NEPUserInterface fui; /* how the user has defined the nonlinear operator */ | ||
141 | PetscBool useds; /* whether the solver uses the DS object or not */ | ||
142 | Mat resolvent; /* shell matrix to be used in NEPApplyResolvent */ | ||
143 | NEPConvergedReason reason; | ||
144 | }; | ||
145 | |||
146 | /* | ||
147 | Macros to test valid NEP arguments | ||
148 | */ | ||
149 | #if !defined(PETSC_USE_DEBUG) | ||
150 | |||
151 | #define NEPCheckProblem(h,arg) do {(void)(h);} while (0) | ||
152 | #define NEPCheckCallback(h,arg) do {(void)(h);} while (0) | ||
153 | #define NEPCheckSplit(h,arg) do {(void)(h);} while (0) | ||
154 | #define NEPCheckDerivatives(h,arg) do {(void)(h);} while (0) | ||
155 | #define NEPCheckSolved(h,arg) do {(void)(h);} while (0) | ||
156 | |||
157 | #else | ||
158 | |||
159 | #define NEPCheckProblem(h,arg) \ | ||
160 | do { \ | ||
161 | PetscCheck(((h)->fui),PetscObjectComm((PetscObject)(h)),PETSC_ERR_ARG_WRONGSTATE,"The nonlinear eigenproblem has not been specified yet. Parameter #%d",arg); \ | ||
162 | } while (0) | ||
163 | |||
164 | #define NEPCheckCallback(h,arg) \ | ||
165 | do { \ | ||
166 | PetscCheck((h)->fui==NEP_USER_INTERFACE_CALLBACK,PetscObjectComm((PetscObject)(h)),PETSC_ERR_ARG_WRONGSTATE,"This operation requires the nonlinear eigenproblem specified with callbacks. Parameter #%d",arg); \ | ||
167 | } while (0) | ||
168 | |||
169 | #define NEPCheckSplit(h,arg) \ | ||
170 | do { \ | ||
171 | PetscCheck((h)->fui==NEP_USER_INTERFACE_SPLIT,PetscObjectComm((PetscObject)(h)),PETSC_ERR_ARG_WRONGSTATE,"This operation requires the nonlinear eigenproblem in split form. Parameter #%d",arg); \ | ||
172 | } while (0) | ||
173 | |||
174 | #define NEPCheckSolved(h,arg) \ | ||
175 | do { \ | ||
176 | PetscCheck((h)->state>=NEP_STATE_SOLVED,PetscObjectComm((PetscObject)(h)),PETSC_ERR_ARG_WRONGSTATE,"Must call NEPSolve() first: Parameter #%d",arg); \ | ||
177 | } while (0) | ||
178 | |||
179 | #endif | ||
180 | |||
181 | /* Check for unsupported features */ | ||
182 | #define NEPCheckUnsupportedCondition(nep,mask,condition,msg) \ | ||
183 | do { \ | ||
184 | if (condition) { \ | ||
185 | PetscCheck(!((mask) & NEP_FEATURE_CALLBACK) || (nep)->fui!=NEP_USER_INTERFACE_CALLBACK,PetscObjectComm((PetscObject)(nep)),PETSC_ERR_SUP,"The solver '%s'%s cannot be used with callback functions (use the split operator)",((PetscObject)(nep))->type_name,(msg)); \ | ||
186 | if ((mask) & NEP_FEATURE_REGION) { \ | ||
187 | PetscBool __istrivial; \ | ||
188 | PetscCall(RGIsTrivial((nep)->rg,&__istrivial)); \ | ||
189 | PetscCheck(__istrivial,PetscObjectComm((PetscObject)(nep)),PETSC_ERR_SUP,"The solver '%s'%s does not support region filtering",((PetscObject)(nep))->type_name,(msg)); \ | ||
190 | } \ | ||
191 | PetscCheck(!((mask) & NEP_FEATURE_CONVERGENCE) || (nep)->converged==NEPConvergedRelative,PetscObjectComm((PetscObject)(nep)),PETSC_ERR_SUP,"The solver '%s'%s only supports the default convergence test",((PetscObject)(nep))->type_name,(msg)); \ | ||
192 | PetscCheck(!((mask) & NEP_FEATURE_STOPPING) || (nep)->stopping==NEPStoppingBasic,PetscObjectComm((PetscObject)(nep)),PETSC_ERR_SUP,"The solver '%s'%s only supports the default stopping test",((PetscObject)(nep))->type_name,(msg)); \ | ||
193 | PetscCheck(!((mask) & NEP_FEATURE_TWOSIDED) || !(nep)->twosided,PetscObjectComm((PetscObject)(nep)),PETSC_ERR_SUP,"The solver '%s'%s cannot compute left eigenvectors (no two-sided variant)",((PetscObject)(nep))->type_name,(msg)); \ | ||
194 | } \ | ||
195 | } while (0) | ||
196 | #define NEPCheckUnsupported(nep,mask) NEPCheckUnsupportedCondition(nep,mask,PETSC_TRUE,"") | ||
197 | |||
198 | /* Check for ignored features */ | ||
199 | #define NEPCheckIgnoredCondition(nep,mask,condition,msg) \ | ||
200 | do { \ | ||
201 | if (condition) { \ | ||
202 | if (((mask) & NEP_FEATURE_CALLBACK) && (nep)->fui==NEP_USER_INTERFACE_CALLBACK) PetscCall(PetscInfo((nep),"The solver '%s'%s ignores the user interface settings\n",((PetscObject)(nep))->type_name,(msg))); \ | ||
203 | if ((mask) & NEP_FEATURE_REGION) { \ | ||
204 | PetscBool __istrivial; \ | ||
205 | PetscCall(RGIsTrivial((nep)->rg,&__istrivial)); \ | ||
206 | if (!__istrivial) PetscCall(PetscInfo((nep),"The solver '%s'%s ignores the specified region\n",((PetscObject)(nep))->type_name,(msg))); \ | ||
207 | } \ | ||
208 | if (((mask) & NEP_FEATURE_CONVERGENCE) && (nep)->converged!=NEPConvergedRelative) PetscCall(PetscInfo((nep),"The solver '%s'%s ignores the convergence test settings\n",((PetscObject)(nep))->type_name,(msg))); \ | ||
209 | if (((mask) & NEP_FEATURE_STOPPING) && (nep)->stopping!=NEPStoppingBasic) PetscCall(PetscInfo((nep),"The solver '%s'%s ignores the stopping test settings\n",((PetscObject)(nep))->type_name,(msg))); \ | ||
210 | if (((mask) & NEP_FEATURE_TWOSIDED) && (nep)->twosided) PetscCall(PetscInfo((nep),"The solver '%s'%s ignores the two-sided flag\n",((PetscObject)(nep))->type_name,(msg))); \ | ||
211 | } \ | ||
212 | } while (0) | ||
213 | #define NEPCheckIgnored(nep,mask) NEPCheckIgnoredCondition(nep,mask,PETSC_TRUE,"") | ||
214 | |||
215 | /* | ||
216 | NEP_KSPSetOperators - Sets the KSP matrices | ||
217 | */ | ||
218 | 7344 | static inline PetscErrorCode NEP_KSPSetOperators(KSP ksp,Mat A,Mat B) | |
219 | { | ||
220 | 7344 | const char *prefix; | |
221 | |||
222 |
1/2✓ Branch 0 taken 9 times.
✗ Branch 1 not taken.
|
7344 | PetscFunctionBegin; |
223 |
4/6✓ Branch 0 taken 9 times.
✓ Branch 1 taken 36 times.
✓ Branch 2 taken 9 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 9 times.
|
7344 | PetscCall(KSPSetOperators(ksp,A,B)); |
224 |
4/6✓ Branch 0 taken 9 times.
✓ Branch 1 taken 36 times.
✓ Branch 2 taken 9 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 9 times.
|
7344 | PetscCall(MatGetOptionsPrefix(B,&prefix)); |
225 |
2/2✓ Branch 0 taken 45 times.
✓ Branch 1 taken 30 times.
|
7344 | if (!prefix) { |
226 | /* set Mat prefix to be the same as KSP to enable setting command-line options (e.g. MUMPS) | ||
227 | only applies if the Mat has no user-defined prefix */ | ||
228 |
4/6✓ Branch 0 taken 9 times.
✓ Branch 1 taken 36 times.
✓ Branch 2 taken 9 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 9 times.
|
5168 | PetscCall(KSPGetOptionsPrefix(ksp,&prefix)); |
229 |
4/6✓ Branch 0 taken 9 times.
✓ Branch 1 taken 36 times.
✓ Branch 2 taken 9 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 9 times.
|
5168 | PetscCall(MatSetOptionsPrefix(B,prefix)); |
230 | } | ||
231 |
6/12✓ Branch 0 taken 9 times.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✓ Branch 3 taken 9 times.
✓ Branch 4 taken 9 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 9 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 9 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✓ Branch 11 taken 9 times.
|
1465 | PetscFunctionReturn(PETSC_SUCCESS); |
232 | } | ||
233 | |||
234 | SLEPC_INTERN PetscErrorCode NEPSetDimensions_Default(NEP,PetscInt,PetscInt*,PetscInt*); | ||
235 | SLEPC_INTERN PetscErrorCode NEPComputeVectors(NEP); | ||
236 | SLEPC_INTERN PetscErrorCode NEPReset_Problem(NEP); | ||
237 | SLEPC_INTERN PetscErrorCode NEPGetDefaultShift(NEP,PetscScalar*); | ||
238 | SLEPC_INTERN PetscErrorCode NEPComputeVectors_Schur(NEP); | ||
239 | SLEPC_INTERN PetscErrorCode NEPComputeResidualNorm_Private(NEP,PetscBool,PetscScalar,Vec,Vec*,PetscReal*); | ||
240 | SLEPC_INTERN PetscErrorCode NEPNewtonRefinementSimple(NEP,PetscInt*,PetscReal,PetscInt); | ||
241 |