//! \brief Base class for computing a step of the Range Kutta algorithm
template <int dim, class butcher_tableau_t>
class EmbeddedRungeKuttaStep : public OptionsHandler<EmbeddedRungeKuttaStepOptions>
{
public:
using vec_t = vector_t<dim>;
//! \brief Constructor with default options
EmbeddedRungeKuttaStep(
const rhs_f<dim>& rhs,
const butcher_tableau_t& butcher):
m_rhs(rhs),
m_butcher(butcher)
{}
//! \brief Constructor with custom options
EmbeddedRungeKuttaStep(
const rhs_f<dim>& rhs,
const butcher_tableau_t& butcher,
EmbeddedRungeKuttaStepOptions& options
):
OptionsHandler(options),
m_rhs(rhs),
m_butcher(butcher)
{}
//! \brief compute a step
//!
//! \param y initial value
//! \param dydx initial value of the right hand site
//! \param x value of the integration variable
//! \param timestep the integration step
//! \param y_out value at the end of the step (y+timestep)
//! \param y_err estimation of the truncation error
void run_step(
const vec_t& y,
const vec_t& dydx,
scalar_t x,
scalar_t timestep,
vec_t& y_out,
vec_t& y_err);
protected:
//! Return the value of the right hand side at x for value y
void get_rhs(
scalar_t x,
const vector_t<dim>& y,
vec_t& dydx
) {return m_rhs(x, y, dydx);}
//! Return the Butcher tableau
butcher_tableau_t& butcher() {return m_butcher;}
private:
rhs_f<dim> m_rhs;
butcher_tableau_t m_butcher;
};
//! Embedded Runge Kutta of fifth order
template <int dim, class butcher_tableau_t>
class EmbeddedRungeKuttaStep45: EmbeddedRungeKuttaStep<dim, butcher_tableau_t>
{
public:
using base = EmbeddedRungeKuttaStep<dim, butcher_tableau_t>;
using vec_t = typename base::vec_t;
using base::butcher;
using base::get_rhs;
using base::get_options;
//check that template argument are consistent
static_assert(butcher_tableau_t::RK_order == 5, "The butcher tableau must correspond to a Runge-Kutta method of order 5");
static_assert(butcher_tableau_t::RK_evaluations == 6 || butcher_tableau_t::RK_evaluations == 7, "The number of functions evaluations must be 6 or 7 for an embedded Runge-Kutta method");
//! \brief Constructor with default options
EmbeddedRungeKuttaStep45(const rhs_f<dim>& rhs,
const butcher_tableau_t& butcher):
base(rhs, butcher)
{
}
//! \brief Constructor with custom options
EmbeddedRungeKuttaStep45(
const rhs_f<dim>& rhs,
const butcher_tableau_t& butcher,
EmbeddedRungeKuttaStepOptions options
):
base(rhs, butcher, options)
{
}
//! \brief compute a step
//!
//! \param y initial value
//! \param dydx initial value of the right hand site
//! \param x value of the integration variable
//! \param timestep the integration step
//! \param y_out value at the end of the step (y+timestep)
//! \param y_err estimation of the truncation error