Medusa  1.1
Coordinate Free Mehless Method implementation
AdamsBashforth.hpp
Go to the documentation of this file.
1 #ifndef MEDUSA_BITS_INTEGRATORS_ADAMSBASHFORTH_HPP_
2 #define MEDUSA_BITS_INTEGRATORS_ADAMSBASHFORTH_HPP_
3 
9 #include "AdamsBashforth_fwd.hpp"
10 #include "RKExplicit.hpp"
11 
12 namespace mm {
13 
14 template <typename Scalar, int num_steps>
15 template <typename func_t, typename initial_method_t>
16 typename AdamsBashforth<Scalar, num_steps>::template
17  Integrator<func_t, initial_method_t>::IterationStep&
18 AdamsBashforth<Scalar, num_steps>::Integrator<func_t,
19  initial_method_t>::IterationStep::operator++() {
20  cur_step++;
21  Eigen::VectorXd next;
22  if (cur_step < steps) { // don't have enough steps, do a RK method.
23  next = integrator.initial_method_.step(
24  integrator.func_, t, last_ys.col(steps - 1), integrator.dt_);
25  } else {
26  next = integrator.method_.step(integrator.func_, t, last_ys, integrator.dt_);
27  }
28  t += integrator.dt_;
29  // shift values
30  for (int i = 0; i < steps-1; ++i) {
31  last_ys.col(i) = last_ys.col(i+1);
32  }
33  last_ys.col(steps-1) = next;
34  return *this;
35 }
36 
38 template <typename scalar_t, int num_steps>
39 std::ostream& operator<<(std::ostream& os, const AdamsBashforth<scalar_t, num_steps>& method) {
40  os << "AdamsBashforth with " << num_steps << " stages and "
41  << "b = " << method.b_.transpose();
42  return os;
43 }
44 
45 namespace integrators {
46 namespace ExplicitMultistep {
48 namespace of_order_internal {
49 
50 template <int order, class scalar_t>
52  static_assert(order > 0, "Order must be positive.");
53  static_assert(order < 6, "Methods of so high orders are not implemented yet, but it's"
54  " very simple to add them!");
55  throw;
56 }
57 
59 template <class scalar_t>
60 struct _of_order_impl<1, scalar_t> {
61  static AdamsBashforth<scalar_t, 1> of_order() { return AB1(); }
62 };
63 
65 template <class scalar_t>
66 struct _of_order_impl<2, scalar_t> {
67  static AdamsBashforth<scalar_t, 2> of_order() { return AB2(); }
68 };
69 
71 template <class scalar_t>
72 struct _of_order_impl<3, scalar_t> {
73  static AdamsBashforth<scalar_t, 3> of_order() { return AB3(); }
74 };
75 
77 template <class scalar_t>
78 struct _of_order_impl<4, scalar_t> {
79  static AdamsBashforth<scalar_t, 4> of_order() { return AB4(); }
80 };
81 
83 template <class scalar_t>
84 struct _of_order_impl<5, scalar_t> {
85  static AdamsBashforth<scalar_t, 5> of_order() { return AB5(); }
86 };
87 
88 } // namespace of_order_internal
90 } // namespace ExplicitMultistep
91 } // namespace integrators
92 } // namespace mm
93 
94 #endif // MEDUSA_BITS_INTEGRATORS_ADAMSBASHFORTH_HPP_
mm::integrators::ExplicitMultistep::AB2
static AdamsBashforth< scalar_t, 2 > AB2()
Two step AB method.
Definition: AdamsBashforth_fwd.hpp:255
mm
Root namespace for the whole library.
Definition: Gaussian.hpp:14
scalar_t
Scalar scalar_t
Type of the elements, alias of Scalar.
Definition: MatrixBaseAddons.hpp:16
mm::integrators::ExplicitMultistep::AB1
static AdamsBashforth< scalar_t, 1 > AB1()
Standard Euler's method.
Definition: AdamsBashforth_fwd.hpp:248
mm::operator<<
std::ostream & operator<<(std::ostream &os, const Gaussian< S > &b)
Output basic information about given Gaussian RBF.
Definition: Gaussian.hpp:37
mm::integrators::ExplicitMultistep::of_order
static AdamsBashforth< scalar_t, order > of_order()
Returns Adams-Bashforth explicit method of requested order with given floating point type.
Definition: AdamsBashforth_fwd.hpp:296
mm::AdamsBashforth::b_
Eigen::Matrix< scalar_t, steps, 1 > b_
coefficients of the multistep method
Definition: AdamsBashforth_fwd.hpp:36
mm::integrators::ExplicitMultistep::AB4
static AdamsBashforth< scalar_t, 4 > AB4()
Four step AB method.
Definition: AdamsBashforth_fwd.hpp:269
mm::integrators::ExplicitMultistep::AB5
static AdamsBashforth< scalar_t, 5 > AB5()
Five step AB method.
Definition: AdamsBashforth_fwd.hpp:276
AdamsBashforth_fwd.hpp
mm::integrators::ExplicitMultistep::AB3
static AdamsBashforth< scalar_t, 3 > AB3()
Three step AB method.
Definition: AdamsBashforth_fwd.hpp:262
mm::AdamsBashforth< scalar_t, num_steps >
RKExplicit.hpp