libMesh
eigen_time_solver.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2017 Benjamin S. Kirk, John W. Peterson, Roy H. Stogner
3 
4 // This library is free software; you can redistribute it and/or
5 // modify it under the terms of the GNU Lesser General Public
6 // License as published by the Free Software Foundation; either
7 // version 2.1 of the License, or (at your option) any later version.
8 
9 // This library is distributed in the hope that it will be useful,
10 // but WITHOUT ANY WARRANTY; without even the implied warranty of
11 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
12 // Lesser General Public License for more details.
13 
14 // You should have received a copy of the GNU Lesser General Public
15 // License along with this library; if not, write to the Free Software
16 // Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307 USA
17 
18 
19 
20 #ifndef LIBMESH_EIGEN_TIME_SOLVER_H
21 #define LIBMESH_EIGEN_TIME_SOLVER_H
22 
23 #include "libmesh/libmesh_config.h"
24 #ifdef LIBMESH_HAVE_SLEPC
25 
26 // Local includes
27 #include "libmesh/time_solver.h"
28 
29 // C++ includes
30 
31 namespace libMesh
32 {
33 
34 // Forward declarations
35 template <typename T> class EigenSolver;
36 
37 
67 {
68 public:
73 
77  typedef TimeSolver Parent;
78 
83  explicit
84  EigenTimeSolver (sys_type & s);
85 
89  virtual ~EigenTimeSolver ();
90 
95  virtual void init () libmesh_override;
96 
101  virtual void reinit () libmesh_override;
102 
107  virtual void solve () libmesh_override;
108 
112  virtual void advance_timestep () libmesh_override {}
113 
118  Real error_order() const { return 0.; }
119 
124  virtual bool element_residual (bool get_jacobian,
125  DiffContext &) libmesh_override;
126 
130  virtual bool side_residual (bool get_jacobian,
131  DiffContext &) libmesh_override;
132 
136  virtual bool nonlocal_residual (bool get_jacobian,
137  DiffContext &) libmesh_override;
138 
144  virtual Real du (const SystemNorm &) const libmesh_override { return 0.; }
145 
149  virtual bool is_steady() const libmesh_override { return true; }
150 
156 
162 
166  unsigned int maxits;
167 
172 
183 
189 
194  unsigned int n_iterations_reqd;
195 
196 private:
197 
203 
208 
213  };
214 
219 };
220 
221 } // namespace libMesh
222 
223 
224 #endif // LIBMESH_HAVE_SLEPC
225 #endif // LIBMESH_EIGEN_TIME_SOLVER_H
This class provides all data required for a physics package (e.g.
Definition: diff_context.h:54
The enum is in an invalid state.
This is a generic class that defines a solver to handle time integration of DifferentiableSystems.
Definition: time_solver.h:58
virtual void solve() libmesh_override
Implements the assembly of both matrices A and B, and calls the EigenSolver to compute the eigenvalue...
virtual void reinit() libmesh_override
The reinitialization function.
This class defines a norm/seminorm to be applied to a NumericVector which contains coefficients in a ...
Definition: system_norm.h:43
The libMesh namespace provides an interface to certain functionality in the library.
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
This class provides a specific system class.
Definition: diff_system.h:53
virtual bool is_steady() const libmesh_override
This is effectively a steady-state solver.
virtual bool nonlocal_residual(bool get_jacobian, DiffContext &) libmesh_override
Forms the jacobian of the nonlocal terms.
Real error_order() const
error convergence order against deltat is not applicable to an eigenvalue problem.
TimeSolver Parent
The parent class.
virtual Real du(const SystemNorm &) const libmesh_override
UniquePtr< EigenSolver< Number > > eigen_solver
The EigenSolver object.
virtual void init() libmesh_override
The initialization function.
virtual bool element_residual(bool get_jacobian, DiffContext &) libmesh_override
Forms either the spatial (Jacobian) or mass matrix part of the operator, depending on which is reques...
unsigned int maxits
The maximum number of iterations allowed to solve the problem.
unsigned int n_iterations_reqd
After a solve, holds the number of iterations required to converge the requested number of eigenpairs...
virtual bool side_residual(bool get_jacobian, DiffContext &) libmesh_override
Forms the jacobian of the boundary terms.
virtual ~EigenTimeSolver()
Destructor.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual void advance_timestep() libmesh_override
It doesn&#39;t make sense to advance the timestep, so we shouldn&#39;t call this.
NowAssembling now_assembling
Flag which controls the internals of element_residual() and side_residual().
unsigned int n_basis_vectors_to_use
The number of basis vectors to use in the computation.
Real tol
The linear solver tolerance to be used when solving the eigenvalue problem.
DifferentiableSystem sys_type
The type of system.
EigenTimeSolver(sys_type &s)
Constructor.
The name of this class is confusing...it&#39;s meant to refer to the base class (TimeSolver) while still ...
unsigned int n_converged_eigenpairs
After a solve, holds the number of eigenpairs successfully converged.
The matrix associated with the spatial part of the operator.
unsigned int n_eigenpairs_to_compute
The number of eigenvectors/values to be computed.
The matrix associated with the time derivative (mass matrix).