eigen_solver.h
Go to the documentation of this file.
1// LIC// ====================================================================
2// LIC// This file forms part of oomph-lib, the object-oriented,
3// LIC// multi-physics finite-element library, available
4// LIC// at http://www.oomph-lib.org.
5// LIC//
6// LIC// Copyright (C) 2006-2022 Matthias Heil and Andrew Hazel
7// LIC//
8// LIC// This library is free software; you can redistribute it and/or
9// LIC// modify it under the terms of the GNU Lesser General Public
10// LIC// License as published by the Free Software Foundation; either
11// LIC// version 2.1 of the License, or (at your option) any later version.
12// LIC//
13// LIC// This library is distributed in the hope that it will be useful,
14// LIC// but WITHOUT ANY WARRANTY; without even the implied warranty of
15// LIC// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
16// LIC// Lesser General Public License for more details.
17// LIC//
18// LIC// You should have received a copy of the GNU Lesser General Public
19// LIC// License along with this library; if not, write to the Free Software
20// LIC// Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA
21// LIC// 02110-1301 USA.
22// LIC//
23// LIC// The authors may be contacted at oomph-lib@maths.man.ac.uk.
24// LIC//
25// LIC//====================================================================
26// Header file for a class that defines interfaces to Eigensolvers
27
28// Include guard to prevent multiple inclusions of the header
29#ifndef OOMPH_EIGEN_SOLVER_HEADER
30#define OOMPH_EIGEN_SOLVER_HEADER
31
32// Include the header generated by autoconfig
33#ifdef HAVE_CONFIG_H
34#include <oomph-lib-config.h>
35#endif
36
37#ifdef OOMPH_HAS_MPI
38#include "mpi.h"
39#endif
40
41#include <complex>
42#include "Vector.h"
43#include "complex_matrices.h"
44
45namespace oomph
46{
47 // Forward definition of problem class
48 class Problem;
49
50 // Forward definition of matrix class
51 class DoubleMatrixBase;
52
53 // Forward definition of linear solver class
54 class LinearSolver;
55
56 //=======================================================================
57 /// Base class for all EigenProblem solves. This simply defines standard
58 /// interfaces so that different solvers can be used easily.
59 //=======================================================================
61 {
62 protected:
63 /// Double value that represents the real part of the shift in
64 /// shifted eigensolvers
65 double Sigma_real;
66
67 public:
68 /// Empty constructor
70
71 /// Empty copy constructor
73
74 /// Empty destructor
75 virtual ~EigenSolver() {}
76
77 /// Actual eigensolver. This takes a pointer to a problem and returns
78 /// a vector of complex numbers representing the eigenvalues
79 /// and a corresponding vector of eigenvectors
80 virtual void solve_eigenproblem(Problem* const& problem_pt,
81 const int& n_eval,
82 Vector<std::complex<double>>& eigenvalue,
83 Vector<DoubleVector>& eigenvector) = 0;
84
85
86 /// Set the value of the shift
87 void set_shift(const double& shift_value)
88 {
89 Sigma_real = shift_value;
90 }
91
92 /// Return the value of the shift (const version)
93 const double& get_shift() const
94 {
95 return Sigma_real;
96 }
97 };
98
99
100 //=====================================================================
101 /// Class for the ARPACK eigensolver
102 //=====================================================================
103 class ARPACK : public EigenSolver
104 {
105 private:
106 /// Pointer to a linear solver
108
109 /// Pointer to a default linear solver
111
112 /// Integer to set whether the real, imaginary or magnitude is
113 /// required
114 /// to be small or large.
116
117
118 /// Number of Arnoldi vectors to compute
120
121
122 /// Boolean to set which part of the spectrum left (default) or right
123 /// of the shifted value.
124 bool Small;
125
126 /// Boolean to indicate whether or not to compute the eigenvectors
128
129
130 public:
131 /// Constructor
132 ARPACK();
133
134 /// Empty copy constructor
135 ARPACK(const ARPACK&) {}
136
137 /// Destructor, delete the linear solver
138 virtual ~ARPACK();
139
140 /// Access function for the number of Arnoldi vectors
141 int& narnoldi()
142 {
143 return NArnoldi;
144 }
145
146 /// Access function for the number of Arnoldi vectors (const version)
147 const int& narnoldi() const
148 {
149 return NArnoldi;
150 }
151
152 /// Set to enable the computation of the eigenvectors (default)
154 {
156 }
157
158 /// Set to disable the computation of the eigenvectors
160 {
161 Compute_eigenvectors = false;
162 }
163
164 /// Solve the eigen problem
165 void solve_eigenproblem(Problem* const& problem_pt,
166 const int& n_eval,
167 Vector<std::complex<double>>& eigenvalue,
168 Vector<DoubleVector>& eigenvector);
169
170 /// Use the eigensolver to find the eigenvalues of a given matrix
171 // void find_eigenvalues(const DoubleMatrixBase &A, const int &n_eval,
172 // Vector<std::complex<double> > &eigenvalue,
173 // Vector<Vector<double> > &eigenvector);
174
175
176 /// Set the desired eigenvalues to be left of the shift
178 {
179 Small = true;
180 }
181
182 /// Set the desired eigenvalues to be right of the shift
184 {
185 Small = false;
186 }
187
188 /// Set the real part to be the quantity of interest (default)
190 {
191 Spectrum = 1;
192 }
193
194 /// Set the imaginary part fo the quantity of interest
196 {
197 Spectrum = -1;
198 }
199
200 /// Set the magnitude to be the quantity of interest
202 {
203 Spectrum = 0;
204 }
205
206 /// Return a pointer to the linear solver object
208 {
209 return Linear_solver_pt;
210 }
211
212 /// Return a pointer to the linear solver object (const version)
214 {
215 return Linear_solver_pt;
216 }
217 };
218
219
220 //=====================================================================
221 /// Class for the LAPACK eigensolver
222 //=====================================================================
223 class LAPACK_QZ : public EigenSolver
224 {
225 public:
226 /// Empty constructor
228
229 /// Empty copy constructor
231
232 /// Empty desctructor
233 virtual ~LAPACK_QZ() {}
234
235 /// Solve the eigen problem
236 void solve_eigenproblem(Problem* const& problem_pt,
237 const int& n_eval,
238 Vector<std::complex<double>>& eigenvalue,
239 Vector<DoubleVector>& eigenvector);
240
241 /// Find the eigenvalues of a generalised eigenvalue problem
242 /// specified by \f$ Ax = \lambda Mx \f$
244 const ComplexMatrixBase& M,
245 Vector<std::complex<double>>& eigenvalue,
246 Vector<Vector<std::complex<double>>>& eigenvector);
247
248 /// Set the desired eigenvalues to be right of the shift
249 /// Dummy at the moment
251 };
252
253} // namespace oomph
254
255#endif
Class for the ARPACK eigensolver.
Definition: eigen_solver.h:104
ARPACK(const ARPACK &)
Empty copy constructor.
Definition: eigen_solver.h:135
void get_eigenvalues_right_of_shift()
Set the desired eigenvalues to be right of the shift.
Definition: eigen_solver.h:183
int Spectrum
Integer to set whether the real, imaginary or magnitude is required to be small or large.
Definition: eigen_solver.h:115
LinearSolver *const & linear_solver_pt() const
Return a pointer to the linear solver object (const version)
Definition: eigen_solver.h:213
int NArnoldi
Number of Arnoldi vectors to compute.
Definition: eigen_solver.h:119
void track_eigenvalue_imaginary_part()
Set the imaginary part fo the quantity of interest.
Definition: eigen_solver.h:195
void get_eigenvalues_left_of_shift()
Use the eigensolver to find the eigenvalues of a given matrix.
Definition: eigen_solver.h:177
LinearSolver * Linear_solver_pt
Pointer to a linear solver.
Definition: eigen_solver.h:107
virtual ~ARPACK()
Destructor, delete the linear solver.
Definition: eigen_solver.cc:62
LinearSolver *& linear_solver_pt()
Return a pointer to the linear solver object.
Definition: eigen_solver.h:207
const int & narnoldi() const
Access function for the number of Arnoldi vectors (const version)
Definition: eigen_solver.h:147
void track_eigenvalue_magnitude()
Set the magnitude to be the quantity of interest.
Definition: eigen_solver.h:201
void enable_compute_eigenvectors()
Set to enable the computation of the eigenvectors (default)
Definition: eigen_solver.h:153
int & narnoldi()
Access function for the number of Arnoldi vectors.
Definition: eigen_solver.h:141
bool Small
Boolean to set which part of the spectrum left (default) or right of the shifted value.
Definition: eigen_solver.h:124
void track_eigenvalue_real_part()
Set the real part to be the quantity of interest (default)
Definition: eigen_solver.h:189
void solve_eigenproblem(Problem *const &problem_pt, const int &n_eval, Vector< std::complex< double > > &eigenvalue, Vector< DoubleVector > &eigenvector)
Solve the eigen problem.
Definition: eigen_solver.cc:72
LinearSolver * Default_linear_solver_pt
Pointer to a default linear solver.
Definition: eigen_solver.h:110
void disable_compute_eigenvectors()
Set to disable the computation of the eigenvectors.
Definition: eigen_solver.h:159
ARPACK()
Constructor.
Definition: eigen_solver.cc:49
bool Compute_eigenvectors
Boolean to indicate whether or not to compute the eigenvectors.
Definition: eigen_solver.h:127
Abstract base class for matrices of complex doubles – adds abstract interfaces for solving,...
Base class for any linear algebra object that is distributable. Just contains storage for the LinearA...
Base class for all EigenProblem solves. This simply defines standard interfaces so that different sol...
Definition: eigen_solver.h:61
virtual ~EigenSolver()
Empty destructor.
Definition: eigen_solver.h:75
EigenSolver(const EigenSolver &)
Empty copy constructor.
Definition: eigen_solver.h:72
EigenSolver()
Empty constructor.
Definition: eigen_solver.h:69
virtual void solve_eigenproblem(Problem *const &problem_pt, const int &n_eval, Vector< std::complex< double > > &eigenvalue, Vector< DoubleVector > &eigenvector)=0
Actual eigensolver. This takes a pointer to a problem and returns a vector of complex numbers represe...
const double & get_shift() const
Return the value of the shift (const version)
Definition: eigen_solver.h:93
void set_shift(const double &shift_value)
Set the value of the shift.
Definition: eigen_solver.h:87
double Sigma_real
Double value that represents the real part of the shift in shifted eigensolvers.
Definition: eigen_solver.h:65
Class for the LAPACK eigensolver.
Definition: eigen_solver.h:224
void find_eigenvalues(const ComplexMatrixBase &A, const ComplexMatrixBase &M, Vector< std::complex< double > > &eigenvalue, Vector< Vector< std::complex< double > > > &eigenvector)
Find the eigenvalues of a generalised eigenvalue problem specified by .
virtual ~LAPACK_QZ()
Empty desctructor.
Definition: eigen_solver.h:233
void get_eigenvalues_right_of_shift()
Set the desired eigenvalues to be right of the shift Dummy at the moment.
Definition: eigen_solver.h:250
LAPACK_QZ()
Empty constructor.
Definition: eigen_solver.h:227
LAPACK_QZ(const LAPACK_QZ &)
Empty copy constructor.
Definition: eigen_solver.h:230
void solve_eigenproblem(Problem *const &problem_pt, const int &n_eval, Vector< std::complex< double > > &eigenvalue, Vector< DoubleVector > &eigenvector)
Solve the eigen problem.
Base class for all linear solvers. This merely defines standard interfaces for linear solvers,...
Definition: linear_solver.h:68
////////////////////////////////////////////////////////////////// //////////////////////////////////...
Definition: problem.h:151
A slight extension to the standard template vector class so that we can include "graceful" array rang...
Definition: Vector.h:58
//////////////////////////////////////////////////////////////////// ////////////////////////////////...