libMesh
unsteady_solver.h
Go to the documentation of this file.
1 // The libMesh Finite Element Library.
2 // Copyright (C) 2002-2024 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_UNSTEADY_SOLVER_H
21 #define LIBMESH_UNSTEADY_SOLVER_H
22 
23 // Local includes
24 #include "libmesh/libmesh_common.h"
25 #include "libmesh/time_solver.h"
26 
27 // C++ includes
28 #include <memory>
29 
30 namespace libMesh
31 {
32 
33 // Forward declarations
34 template <typename T> class NumericVector;
35 
50 class UnsteadySolver : public TimeSolver
51 {
52 public:
57  explicit
59 
63  virtual ~UnsteadySolver ();
64 
69  virtual void init () override;
70 
75  virtual void init_adjoints () override;
76 
82  virtual void init_data () override;
83 
88  virtual void reinit () override;
89 
95  virtual void solve () override;
96 
103  virtual void advance_timestep () override;
104 
105  void update();
106 
111  virtual std::pair<unsigned int, Real> adjoint_solve (const QoISet & qoi_indices) override;
112 
118  virtual void adjoint_advance_timestep () override;
119 
124  virtual void retrieve_timestep () override;
125 
129  virtual void integrate_qoi_timestep() override;
130 
136  virtual void integrate_adjoint_sensitivity(const QoISet & qois, const ParameterVector & parameter_vector, SensitivityData & sensitivities) override;
137 
138 #ifdef LIBMESH_ENABLE_AMR
139 
146  virtual void integrate_adjoint_refinement_error_estimate(AdjointRefinementEstimator & /*adjoint_refinement_error_estimator*/, ErrorVector & /*QoI_elementwise_error*/) override;
147 #endif // LIBMESH_ENABLE_AMR
148 
156  virtual Real error_order () const = 0;
157 
165  virtual unsigned int time_order () const = 0;
166 
171  Number old_nonlinear_solution (const dof_id_type global_dof_number) const;
172 
178  std::shared_ptr<NumericVector<Number>> old_local_nonlinear_solution;
179 
189  virtual Real du(const SystemNorm & norm) const override;
190 
194  virtual bool is_steady() const override { return false; }
195 
199  void set_first_adjoint_step(bool first_adjoint_step_setting)
200  {
201  first_adjoint_step = first_adjoint_step_setting;
202  }
203 
204  /*
205  * A setter for the first_solver boolean. Useful for example if we are using
206  * a nested time solver, and the outer solver wants to tell the inner one that
207  * the initial conditions have already been handled.
208  */
209  void set_first_solve(bool first_solve_setting)
210  {
211  first_solve = first_solve_setting;
212  }
213 
214 protected:
215 
221 
227 
231  std::vector< std::unique_ptr<NumericVector<Number>> > old_adjoints;
232 
239 };
240 
241 
242 
243 } // namespace libMesh
244 
245 
246 #endif // LIBMESH_UNSTEADY_SOLVER_H
virtual void reinit() override
The reinitialization function.
virtual void integrate_adjoint_refinement_error_estimate(AdjointRefinementEstimator &, ErrorVector &) override
A method to compute the adjoint refinement error estimate at the current timestep.
virtual std::pair< unsigned int, Real > adjoint_solve(const QoISet &qoi_indices) override
This method solves for the adjoint solution at the next adjoint timestep (or a steady state adjoint s...
virtual bool is_steady() const override
This is not a steady-state solver.
UnsteadySolver(sys_type &s)
Constructor.
Data structure for specifying which Parameters should be independent variables in a parameter sensiti...
This is a generic class that defines a solver to handle time integration of DifferentiableSystems.
Definition: time_solver.h:63
This class implements a "brute force" goal-oriented error estimator which computes an estimate of err...
void set_first_adjoint_step(bool first_adjoint_step_setting)
A setter for the first_adjoint_step boolean.
Data structure for specifying which Quantities of Interest should be calculated in an adjoint or a pa...
Definition: qoi_set.h:45
The ErrorVector is a specialization of the StatisticsVector for error data computed on a finite eleme...
Definition: error_vector.h:50
virtual void init_data() override
The data initialization function.
virtual void integrate_adjoint_sensitivity(const QoISet &qois, const ParameterVector &parameter_vector, SensitivityData &sensitivities) override
A method to integrate the adjoint sensitivity w.r.t a given parameter vector.
Number old_nonlinear_solution(const dof_id_type global_dof_number) const
virtual Real error_order() const =0
This method should return the expected convergence order of the (non-local) error of the time discret...
This class defines a norm/seminorm to be applied to a NumericVector which contains coefficients in a ...
Definition: system_norm.h:45
The libMesh namespace provides an interface to certain functionality in the library.
void set_first_solve(bool first_solve_setting)
This class provides a specific system class.
Definition: diff_system.h:54
std::vector< std::unique_ptr< NumericVector< Number > > > old_adjoints
A vector of pointers to vectors holding the adjoint solution at the last time step.
virtual void adjoint_advance_timestep() override
This method advances the adjoint solution to the previous timestep, after an adjoint_solve() has been...
virtual unsigned int time_order() const =0
virtual ~UnsteadySolver()
Destructor.
std::shared_ptr< NumericVector< Number > > old_local_nonlinear_solution
Serial vector of _system.get_vector("_old_nonlinear_solution") This is a shared_ptr so that it can be...
Data structure for holding completed parameter sensitivity calculations.
virtual void retrieve_timestep() override
This method retrieves all the stored solutions at the current system.time.
auto norm(const T &a) -> decltype(std::abs(a))
Definition: tensor_tools.h:74
virtual void integrate_qoi_timestep() override
A method to integrate the system::QoI functionals.
virtual void solve() override
This method solves for the solution at the next timestep.
virtual void advance_timestep() override
This method advances the solution to the next timestep, after a solve() has been performed.
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
virtual void init_adjoints() override
Add adjoint vectors and old_adjoint_vectors as per the indices of QoISet.
virtual Real du(const SystemNorm &norm) const override
Computes the size of ||u^{n+1} - u^{n}|| in some norm.
Real last_step_deltat
We will need to move the system.time around to ensure that residuals are built with the right deltat ...
virtual void init() override
The initialization function.
bool first_adjoint_step
A bool that will be true the first time adjoint_advance_timestep() is called, (when the primal soluti...
uint8_t dof_id_type
Definition: id_types.h:67
bool first_solve
A bool that will be true the first time solve() is called, and false thereafter.