You can not select more than 25 topics
Topics must start with a letter or number, can include dashes ('-') and can be up to 35 characters long.
94 lines
2.8 KiB
94 lines
2.8 KiB
#ifndef MEDUSA_BITS_INTEGRATORS_ADAMSBASHFORTH_HPP_
|
|
#define MEDUSA_BITS_INTEGRATORS_ADAMSBASHFORTH_HPP_
|
|
|
|
/**
|
|
* @file
|
|
* Adams-Bashforth integrator implementation.
|
|
*/
|
|
|
|
#include "AdamsBashforth_fwd.hpp"
|
|
#include "RKExplicit.hpp"
|
|
|
|
namespace mm {
|
|
|
|
template <typename Scalar, int num_steps>
|
|
template <typename func_t, typename initial_method_t>
|
|
typename AdamsBashforth<Scalar, num_steps>::template
|
|
Integrator<func_t, initial_method_t>::IterationStep&
|
|
AdamsBashforth<Scalar, num_steps>::Integrator<func_t,
|
|
initial_method_t>::IterationStep::operator++() {
|
|
cur_step++;
|
|
Eigen::VectorXd next;
|
|
if (cur_step < steps) { // don't have enough steps, do a RK method.
|
|
next = integrator.initial_method_.step(
|
|
integrator.func_, t, last_ys.col(steps - 1), integrator.dt_);
|
|
} else {
|
|
next = integrator.method_.step(integrator.func_, t, last_ys, integrator.dt_);
|
|
}
|
|
t += integrator.dt_;
|
|
// shift values
|
|
for (int i = 0; i < steps-1; ++i) {
|
|
last_ys.col(i) = last_ys.col(i+1);
|
|
}
|
|
last_ys.col(steps-1) = next;
|
|
return *this;
|
|
}
|
|
|
|
/// Output the method's tableau for debugging.
|
|
template <typename scalar_t, int num_steps>
|
|
std::ostream& operator<<(std::ostream& os, const AdamsBashforth<scalar_t, num_steps>& method) {
|
|
os << "AdamsBashforth with " << num_steps << " stages and "
|
|
<< "b = " << method.b_.transpose();
|
|
return os;
|
|
}
|
|
|
|
namespace integrators {
|
|
namespace ExplicitMultistep {
|
|
/// @cond
|
|
namespace of_order_internal {
|
|
|
|
template <int order, class scalar_t>
|
|
AdamsBashforth<scalar_t, order> _of_order_impl<order, scalar_t>::of_order() {
|
|
static_assert(order > 0, "Order must be positive.");
|
|
static_assert(order < 6, "Methods of so high orders are not implemented yet, but it's"
|
|
" very simple to add them!");
|
|
throw;
|
|
}
|
|
|
|
/// Specialization of order 1
|
|
template <class scalar_t>
|
|
struct _of_order_impl<1, scalar_t> {
|
|
static AdamsBashforth<scalar_t, 1> of_order() { return AB1(); }
|
|
};
|
|
|
|
/// Specialization of order 2
|
|
template <class scalar_t>
|
|
struct _of_order_impl<2, scalar_t> {
|
|
static AdamsBashforth<scalar_t, 2> of_order() { return AB2(); }
|
|
};
|
|
|
|
/// Specialization of order 3
|
|
template <class scalar_t>
|
|
struct _of_order_impl<3, scalar_t> {
|
|
static AdamsBashforth<scalar_t, 3> of_order() { return AB3(); }
|
|
};
|
|
|
|
/// Specialization of order 4
|
|
template <class scalar_t>
|
|
struct _of_order_impl<4, scalar_t> {
|
|
static AdamsBashforth<scalar_t, 4> of_order() { return AB4(); }
|
|
};
|
|
|
|
/// Specialization of order 5
|
|
template <class scalar_t>
|
|
struct _of_order_impl<5, scalar_t> {
|
|
static AdamsBashforth<scalar_t, 5> of_order() { return AB5(); }
|
|
};
|
|
|
|
} // namespace of_order_internal
|
|
/// @endcond
|
|
} // namespace ExplicitMultistep
|
|
} // namespace integrators
|
|
} // namespace mm
|
|
|
|
#endif // MEDUSA_BITS_INTEGRATORS_ADAMSBASHFORTH_HPP_
|
|
|