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