20 #ifndef LIBMESH_UNSTEADY_SOLVER_H 21 #define LIBMESH_UNSTEADY_SOLVER_H 24 #include "libmesh/libmesh_common.h" 25 #include "libmesh/time_solver.h" 34 template <
typename T>
class NumericVector;
69 virtual void init ()
override;
88 virtual void reinit ()
override;
95 virtual void solve ()
override;
138 #ifdef LIBMESH_ENABLE_AMR 147 #endif // LIBMESH_ENABLE_AMR 194 virtual bool is_steady()
const override {
return false; }
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.
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...
The ErrorVector is a specialization of the StatisticsVector for error data computed on a finite eleme...
virtual void init_data() override
The data initialization function.
virtual void integrate_adjoint_sensitivity(const QoISet &qois, const ParameterVector ¶meter_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 ...
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.
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.
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...
bool first_solve
A bool that will be true the first time solve() is called, and false thereafter.