libMesh
adaptive_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_ADAPTIVE_TIME_SOLVER_H
21 #define LIBMESH_ADAPTIVE_TIME_SOLVER_H
22 
23 // Local includes
24 #include "libmesh/system_norm.h"
25 #include "libmesh/first_order_unsteady_solver.h"
26 
27 // C++ includes
28 
29 namespace libMesh
30 {
31 
32 // Forward declarations
33 class System;
34 
50 {
51 public:
56 
61  explicit
63 
67  virtual ~AdaptiveTimeSolver ();
68 
69  virtual void init() libmesh_override;
70 
71  virtual void reinit() libmesh_override;
72 
73  virtual void solve() libmesh_override = 0;
74 
75  virtual void advance_timestep() libmesh_override;
76 
80  virtual Real error_order () const libmesh_override;
81 
85  virtual bool element_residual (bool get_jacobian,
86  DiffContext &) libmesh_override;
87 
91  virtual bool side_residual (bool get_jacobian,
92  DiffContext &) libmesh_override;
93 
97  virtual bool nonlocal_residual (bool get_jacobian,
98  DiffContext &) libmesh_override;
99 
103  virtual UniquePtr<DiffSolver> & diff_solver() libmesh_override;
104 
109  virtual UniquePtr<LinearSolver<Number>> & linear_solver() libmesh_override;
110 
115 
120 
127  std::vector<float> component_scale;
128 
145 
162 
168 
174 
182 
199 
200 protected:
201 
208 
213 };
214 
215 
216 } // namespace libMesh
217 
218 
219 #endif // LIBMESH_ADAPTIVE_TIME_SOLVER_H
virtual Real calculate_norm(System &, NumericVector< Number > &)
A helper function to calculate error norms.
std::vector< float > component_scale
If component_norms is non-empty, each variable&#39;s contribution to the error of a system will also be s...
virtual bool side_residual(bool get_jacobian, DiffContext &) libmesh_override
This method is passed on to the core_time_solver.
This class provides all data required for a physics package (e.g.
Definition: diff_context.h:54
SystemNorm component_norm
Error calculations are done in this norm, DISCRETE_L2 by default.
Real last_deltat
We need to store the value of the last deltat used, so that advance_timestep() will increment the sys...
Real target_tolerance
This tolerance is the target relative error between an exact time integration and a single time step ...
Numeric vector.
Definition: dof_map.h:66
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.
virtual UniquePtr< LinearSolver< Number > > & linear_solver() libmesh_override
An implicit linear solver to use for adjoint and sensitivity problems.
virtual void advance_timestep() libmesh_override
This method advances the solution to the next timestep, after a solve() has been performed.
FirstOrderUnsteadySolver Parent
The parent class.
AdaptiveTimeSolver(sys_type &s)
Constructor.
std::unique_ptr< T > UniquePtr
Definition: auto_ptr.h:46
virtual bool element_residual(bool get_jacobian, DiffContext &) libmesh_override
This method is passed on to the core_time_solver.
This class provides a specific system class.
Definition: diff_system.h:53
This is a generic class that defines a solver to handle ImplicitSystem classes, including NonlinearIm...
Definition: diff_solver.h:70
bool global_tolerance
This flag, which is true by default, grows (shrinks) the timestep based on the expected global accura...
This base class can be inherited from to provide interfaces to linear solvers from different packages...
Definition: linear_solver.h:58
Real max_deltat
Do not allow the adaptive time solver to select deltat > max_deltat.
Generic class from which first order UnsteadySolvers should subclass.
This is the base class for classes which contain information related to any physical process that mig...
Definition: system.h:76
Real upper_tolerance
This tolerance is the maximum relative error between an exact time integration and a single time step...
This class wraps another UnsteadySolver derived class, and compares the results of timestepping with ...
virtual void init() libmesh_override
The initialization function.
virtual void reinit() libmesh_override
The reinitialization function.
virtual void solve() libmesh_override=0
This method solves for the solution at the next timestep.
This is a generic class that defines a solver to handle time integration of DifferentiableSystems.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
Real min_deltat
Do not allow the adaptive time solver to select deltat < min_deltat.
virtual ~AdaptiveTimeSolver()
Destructor.
UniquePtr< UnsteadySolver > core_time_solver
This object is used to take timesteps.
virtual UniquePtr< DiffSolver > & diff_solver() libmesh_override
An implicit linear or nonlinear solver to use at each timestep.
virtual Real error_order() const libmesh_override
This method is passed on to the core_time_solver.
virtual bool nonlocal_residual(bool get_jacobian, DiffContext &) libmesh_override
This method is passed on to the core_time_solver.
Real max_growth
Do not allow the adaptive time solver to select a new deltat greater than max_growth times the old de...