19 #include "libmesh/twostep_time_solver.h" 21 #include "libmesh/adjoint_refinement_estimator.h" 22 #include "libmesh/diff_system.h" 23 #include "libmesh/enum_norm_type.h" 24 #include "libmesh/error_vector.h" 25 #include "libmesh/euler_solver.h" 26 #include "libmesh/int_range.h" 27 #include "libmesh/numeric_vector.h" 28 #include "libmesh/parameter_vector.h" 29 #include "libmesh/sensitivity_data.h" 30 #include "libmesh/solution_history.h" 69 bool max_tolerance_met =
false;
72 Real single_norm(0.), double_norm(0.), error_norm(0.),
81 while (!max_tolerance_met)
90 libMesh::out <<
"\n === Computing adaptive timestep === " 100 std::unique_ptr<NumericVector<Number>> double_solution =
102 std::unique_ptr<NumericVector<Number>> old_solution =
108 libMesh::out <<
"Double norm = " << double_norm << std::endl;
140 libMesh::out <<
"Single norm = " << single_norm << std::endl;
150 relative_error = error_norm / old_deltat /
151 std::max(double_norm, single_norm);
154 if (!double_norm && !single_norm)
159 libMesh::out <<
"Error norm = " << error_norm << std::endl;
162 std::max(double_norm, single_norm))
165 << (error_norm / old_deltat /
166 std::max(double_norm, single_norm))
168 libMesh::out <<
"old delta t = " << old_deltat << std::endl;
212 max_tolerance_met =
true;
232 const Real global_shrink_or_growth_factor =
236 const Real local_shrink_or_growth_factor =
238 (error_norm/std::max(double_norm, single_norm)),
244 << global_shrink_or_growth_factor << std::endl;
246 << local_shrink_or_growth_factor << std::endl;
255 Real shrink_or_growth_factor =
257 local_shrink_or_growth_factor;
263 libMesh::out <<
"delta t is constrained by max_growth" << std::endl;
275 libMesh::out <<
"delta t is constrained by maximum-allowable delta t." 286 libMesh::out <<
"delta t is constrained by minimum-allowable delta t." 327 std::pair<unsigned int, Real> full_adjoint_output =
core_time_solver->adjoint_solve(qoi_indices);
338 return full_adjoint_output;
385 core_time_solver->integrate_adjoint_sensitivity(qois, parameter_vector, sensitivities_first_half);
390 core_time_solver->integrate_adjoint_sensitivity(qois, parameter_vector, sensitivities_second_half);
393 const auto pv_size = parameter_vector.
size();
396 sensitivities[i][j] = sensitivities_first_half[i][j] + sensitivities_second_half[i][j];
399 #ifdef LIBMESH_ENABLE_AMR 405 std::vector<Number> qoi_error_estimates_first_half(
_system.
n_qois());
406 std::vector<Number> qoi_error_estimates_second_half(
_system.
n_qois());
411 core_time_solver->integrate_adjoint_refinement_error_estimate(adjoint_refinement_error_estimator, QoI_elementwise_error_first_half);
426 core_time_solver->integrate_adjoint_refinement_error_estimate(adjoint_refinement_error_estimator, QoI_elementwise_error_second_half);
440 QoI_elementwise_error[i] = QoI_elementwise_error_first_half[i] + QoI_elementwise_error_second_half[i];
451 #endif // LIBMESH_ENABLE_AMR bool quiet
Print extra debugging information if quiet == false.
virtual Real calculate_norm(System &, NumericVector< Number > &)
A helper function to calculate error norms.
Real time
For time-dependent problems, this is the time t at the beginning of the current timestep.
QoISet & qoi_set()
Access to the QoISet (default: weight all QoIs equally) to use when computing errors.
Real completed_timestep_size
The adaptive time solver's have two notions of deltat.
Data structure for specifying which Parameters should be independent variables in a parameter sensiti...
This class implements a "brute force" goal-oriented error estimator which computes an estimate of err...
Data structure for specifying which Quantities of Interest should be calculated in an adjoint or a pa...
Number get_qoi_value(unsigned int qoi_index) const
unsigned int n_qois() const
Number of currently active quantities of interest.
The ErrorVector is a specialization of the StatisticsVector for error data computed on a finite eleme...
Real target_tolerance
This tolerance is the target relative error between an exact time integration and a single time step ...
virtual std::unique_ptr< NumericVector< T > > clone() const =0
The libMesh namespace provides an interface to certain functionality in the library.
virtual void integrate_qoi_timestep() override
A method to integrate the system::QoI functionals.
virtual void integrate_adjoint_refinement_error_estimate(AdjointRefinementEstimator &adjoint_refinement_error_estimator, ErrorVector &QoI_elementwise_error) override
A method to compute the adjoint refinement error estimate at the current timestep.
This class provides a specific system class.
bool has_index(std::size_t) const
Return whether or not this index is in the set to be calculated.
bool global_tolerance
This flag, which is true by default, grows (shrinks) the timestep based on the expected global accura...
sys_type & _system
A reference to the system we are solving.
Real max_deltat
Do not allow the adaptive time solver to select deltat > max_deltat.
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...
Data structure for holding completed parameter sensitivity calculations.
Real deltat
For time-dependent problems, this is the amount delta t to advance the solution in time...
Real upper_tolerance
This tolerance is the maximum relative error between an exact time integration and a single time step...
std::unique_ptr< NumericVector< Number > > solution
Data structure to hold solution values.
Real calculate_norm(const NumericVector< Number > &v, unsigned int var, FEMNormType norm_type, std::set< unsigned int > *skip_dimensions=nullptr) const
This class wraps another UnsteadySolver derived class, and compares the results of timestepping with ...
unsigned int reduce_deltat_on_diffsolver_failure
This value (which defaults to zero) is the number of times the TimeSolver is allowed to halve deltat ...
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.
std::size_t size(const System &sys) const
TwostepTimeSolver(sys_type &s)
Constructor.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void set_qoi_error_estimate(unsigned int qoi_index, Number qoi_error_estimate)
~TwostepTimeSolver()
Destructor.
Real min_deltat
Do not allow the adaptive time solver to select deltat < min_deltat.
IntRange< T > make_range(T beg, T end)
The 2-parameter make_range() helper function returns an IntRange<T> when both input parameters are of...
NumericVector< Number > & get_adjoint_solution(unsigned int i=0)
void set_qoi(unsigned int qoi_index, Number qoi_value)
Real last_deltat
The deltat for the last completed timestep before the current one.
virtual void advance_timestep() override
This method advances the solution to the next timestep, after a solve() has been performed.
Number get_qoi_error_estimate_value(unsigned int qoi_index) const
virtual void solve() override
This method solves for the solution at the next timestep.
unsigned int n_vars() const
std::unique_ptr< UnsteadySolver > core_time_solver
This object is used to take timesteps.
auto index_range(const T &sizable)
Helper function that returns an IntRange<std::size_t> representing all the indices of the passed-in v...
const NumericVector< Number > & get_vector(std::string_view vec_name) const
bool first_solve
A bool that will be true the first time solve() is called, and false thereafter.
Real max_growth
Do not allow the adaptive time solver to select a new deltat greater than max_growth times the old de...