adams_bashforth.hpp 15 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418
  1. /*
  2. [auto_generated]
  3. boost/numeric/odeint/stepper/adams_bashforth.hpp
  4. [begin_description]
  5. Implementaton of the Adam-Bashforth method a multistep method used for the predictor step in the
  6. Adams-Bashforth-Moulton method.
  7. [end_description]
  8. Copyright 2011-2013 Karsten Ahnert
  9. Copyright 2011-2013 Mario Mulansky
  10. Copyright 2012 Christoph Koke
  11. Copyright 2013 Pascal Germroth
  12. Distributed under the Boost Software License, Version 1.0.
  13. (See accompanying file LICENSE_1_0.txt or
  14. copy at http://www.boost.org/LICENSE_1_0.txt)
  15. */
  16. #ifndef BOOST_NUMERIC_ODEINT_STEPPER_ADAMS_BASHFORTH_HPP_INCLUDED
  17. #define BOOST_NUMERIC_ODEINT_STEPPER_ADAMS_BASHFORTH_HPP_INCLUDED
  18. #include <boost/static_assert.hpp>
  19. #include <boost/numeric/odeint/util/unwrap_reference.hpp>
  20. #include <boost/numeric/odeint/algebra/range_algebra.hpp>
  21. #include <boost/numeric/odeint/algebra/default_operations.hpp>
  22. #include <boost/numeric/odeint/algebra/algebra_dispatcher.hpp>
  23. #include <boost/numeric/odeint/algebra/operations_dispatcher.hpp>
  24. #include <boost/numeric/odeint/util/state_wrapper.hpp>
  25. #include <boost/numeric/odeint/util/is_resizeable.hpp>
  26. #include <boost/numeric/odeint/util/resizer.hpp>
  27. #include <boost/numeric/odeint/stepper/stepper_categories.hpp>
  28. #include <boost/numeric/odeint/stepper/runge_kutta4.hpp>
  29. #include <boost/numeric/odeint/stepper/extrapolation_stepper.hpp>
  30. #include <boost/numeric/odeint/stepper/base/algebra_stepper_base.hpp>
  31. #include <boost/numeric/odeint/stepper/detail/adams_bashforth_coefficients.hpp>
  32. #include <boost/numeric/odeint/stepper/detail/adams_bashforth_call_algebra.hpp>
  33. #include <boost/numeric/odeint/stepper/detail/rotating_buffer.hpp>
  34. #include <boost/mpl/arithmetic.hpp>
  35. #include <boost/mpl/min_max.hpp>
  36. #include <boost/mpl/equal_to.hpp>
  37. namespace mpl = boost::mpl;
  38. namespace boost {
  39. namespace numeric {
  40. namespace odeint {
  41. using mpl::int_;
  42. /* if N >= 4, returns the smallest even number > N, otherwise returns 4 */
  43. template < int N >
  44. struct order_helper
  45. : mpl::max< typename mpl::eval_if<
  46. mpl::equal_to< mpl::modulus< int_< N >, int_< 2 > >,
  47. int_< 0 > >,
  48. int_< N >, int_< N + 1 > >::type,
  49. int_< 4 > >::type
  50. { };
  51. template<
  52. size_t Steps ,
  53. class State ,
  54. class Value = double ,
  55. class Deriv = State ,
  56. class Time = Value ,
  57. class Algebra = typename algebra_dispatcher< State >::algebra_type ,
  58. class Operations = typename operations_dispatcher< State >::operations_type ,
  59. class Resizer = initially_resizer ,
  60. class InitializingStepper = extrapolation_stepper< order_helper<Steps>::value,
  61. State, Value, Deriv, Time,
  62. Algebra, Operations, Resizer >
  63. >
  64. class adams_bashforth : public algebra_stepper_base< Algebra , Operations >
  65. {
  66. #ifndef DOXYGEN_SKIP
  67. static_assert(( Steps > 0 && Steps < 9 ), "Must have between 1 and 8 steps inclusive");
  68. #endif
  69. public :
  70. typedef State state_type;
  71. typedef state_wrapper< state_type > wrapped_state_type;
  72. typedef Value value_type;
  73. typedef Deriv deriv_type;
  74. typedef state_wrapper< deriv_type > wrapped_deriv_type;
  75. typedef Time time_type;
  76. typedef Resizer resizer_type;
  77. typedef stepper_tag stepper_category;
  78. typedef InitializingStepper initializing_stepper_type;
  79. typedef algebra_stepper_base< Algebra , Operations > algebra_stepper_base_type;
  80. typedef typename algebra_stepper_base_type::algebra_type algebra_type;
  81. typedef typename algebra_stepper_base_type::operations_type operations_type;
  82. #ifndef DOXYGEN_SKIP
  83. typedef adams_bashforth< Steps , State , Value , Deriv , Time , Algebra , Operations , Resizer , InitializingStepper > stepper_type;
  84. #endif
  85. static const size_t steps = Steps;
  86. typedef unsigned short order_type;
  87. static const order_type order_value = steps;
  88. typedef detail::rotating_buffer< wrapped_deriv_type , steps > step_storage_type;
  89. order_type order( void ) const { return order_value; }
  90. adams_bashforth( const algebra_type &algebra = algebra_type() )
  91. : algebra_stepper_base_type( algebra ) ,
  92. m_step_storage() , m_resizer() , m_coefficients() ,
  93. m_steps_initialized( 0 ) , m_initializing_stepper()
  94. { }
  95. /*
  96. * Version 1 : do_step( system , x , t , dt );
  97. *
  98. * solves the forwarding problem
  99. */
  100. template< class System , class StateInOut >
  101. void do_step( System system , StateInOut &x , time_type t , time_type dt )
  102. {
  103. do_step( system , x , t , x , dt );
  104. }
  105. /**
  106. * \brief Second version to solve the forwarding problem, can be called with Boost.Range as StateInOut.
  107. */
  108. template< class System , class StateInOut >
  109. void do_step( System system , const StateInOut &x , time_type t , time_type dt )
  110. {
  111. do_step( system , x , t , x , dt );
  112. }
  113. /*
  114. * Version 2 : do_step( system , in , t , out , dt );
  115. *
  116. * solves the forwarding problem
  117. */
  118. template< class System , class StateIn , class StateOut >
  119. void do_step( System system , const StateIn &in , time_type t , StateOut &out , time_type dt )
  120. {
  121. do_step_impl( system , in , t , out , dt );
  122. }
  123. /**
  124. * \brief Second version to solve the forwarding problem, can be called with Boost.Range as StateOut.
  125. */
  126. template< class System , class StateIn , class StateOut >
  127. void do_step( System system , const StateIn &in , time_type t , const StateOut &out , time_type dt )
  128. {
  129. do_step_impl( system , in , t , out , dt );
  130. }
  131. template< class StateType >
  132. void adjust_size( const StateType &x )
  133. {
  134. resize_impl( x );
  135. }
  136. const step_storage_type& step_storage( void ) const
  137. {
  138. return m_step_storage;
  139. }
  140. step_storage_type& step_storage( void )
  141. {
  142. return m_step_storage;
  143. }
  144. template< class ExplicitStepper , class System , class StateIn >
  145. void initialize( ExplicitStepper explicit_stepper , System system , StateIn &x , time_type &t , time_type dt )
  146. {
  147. typename odeint::unwrap_reference< ExplicitStepper >::type &stepper = explicit_stepper;
  148. typename odeint::unwrap_reference< System >::type &sys = system;
  149. m_resizer.adjust_size(x, [this](auto&& arg) { return this->resize_impl<StateIn>(std::forward<decltype(arg)>(arg)); });
  150. for( size_t i=0 ; i+1<steps ; ++i )
  151. {
  152. if( i != 0 ) m_step_storage.rotate();
  153. sys( x , m_step_storage[0].m_v , t );
  154. stepper.do_step_dxdt_impl( system, x, m_step_storage[0].m_v, t,
  155. dt );
  156. t += dt;
  157. }
  158. m_steps_initialized = steps;
  159. }
  160. template< class System , class StateIn >
  161. void initialize( System system , StateIn &x , time_type &t , time_type dt )
  162. {
  163. initialize( std::ref( m_initializing_stepper ) , system , x , t , dt );
  164. }
  165. void reset( void )
  166. {
  167. m_steps_initialized = 0;
  168. }
  169. bool is_initialized( void ) const
  170. {
  171. return m_steps_initialized >= ( steps - 1 );
  172. }
  173. const initializing_stepper_type& initializing_stepper( void ) const { return m_initializing_stepper; }
  174. initializing_stepper_type& initializing_stepper( void ) { return m_initializing_stepper; }
  175. private:
  176. template< class System , class StateIn , class StateOut >
  177. void do_step_impl( System system , const StateIn &in , time_type t , StateOut &out , time_type dt )
  178. {
  179. typename odeint::unwrap_reference< System >::type &sys = system;
  180. if( m_resizer.adjust_size(in, [this](auto&& arg) { return this->resize_impl<StateIn>(std::forward<decltype(arg)>(arg)); } ) )
  181. {
  182. m_steps_initialized = 0;
  183. }
  184. if( m_steps_initialized + 1 < steps )
  185. {
  186. if( m_steps_initialized != 0 ) m_step_storage.rotate();
  187. sys( in , m_step_storage[0].m_v , t );
  188. m_initializing_stepper.do_step_dxdt_impl(
  189. system, in, m_step_storage[0].m_v, t, out, dt );
  190. ++m_steps_initialized;
  191. }
  192. else
  193. {
  194. m_step_storage.rotate();
  195. sys( in , m_step_storage[0].m_v , t );
  196. detail::adams_bashforth_call_algebra< steps , algebra_type , operations_type >()( this->m_algebra , in , out , m_step_storage , m_coefficients , dt );
  197. }
  198. }
  199. template< class StateIn >
  200. bool resize_impl( const StateIn &x )
  201. {
  202. bool resized( false );
  203. for( size_t i=0 ; i<steps ; ++i )
  204. {
  205. resized |= adjust_size_by_resizeability( m_step_storage[i] , x , typename is_resizeable<deriv_type>::type() );
  206. }
  207. return resized;
  208. }
  209. step_storage_type m_step_storage;
  210. resizer_type m_resizer;
  211. detail::adams_bashforth_coefficients< value_type , steps > m_coefficients;
  212. size_t m_steps_initialized;
  213. initializing_stepper_type m_initializing_stepper;
  214. };
  215. /***** DOXYGEN *****/
  216. /**
  217. * \class adams_bashforth
  218. * \brief The Adams-Bashforth multistep algorithm.
  219. *
  220. * The Adams-Bashforth method is a multi-step algorithm with configurable step
  221. * number. The step number is specified as template parameter Steps and it
  222. * then uses the result from the previous Steps steps. See also
  223. * <a href="http://en.wikipedia.org/wiki/Linear_multistep_method">en.wikipedia.org/wiki/Linear_multistep_method</a>.
  224. * Currently, a maximum of Steps=8 is supported.
  225. * The method is explicit and fulfills the Stepper concept. Step size control
  226. * or continuous output are not provided.
  227. *
  228. * This class derives from algebra_base and inherits its interface via
  229. * CRTP (current recurring template pattern). For more details see
  230. * algebra_stepper_base.
  231. *
  232. * \tparam Steps The number of steps (maximal 8).
  233. * \tparam State The state type.
  234. * \tparam Value The value type.
  235. * \tparam Deriv The type representing the time derivative of the state.
  236. * \tparam Time The time representing the independent variable - the time.
  237. * \tparam Algebra The algebra type.
  238. * \tparam Operations The operations type.
  239. * \tparam Resizer The resizer policy type.
  240. * \tparam InitializingStepper The stepper for the first two steps.
  241. */
  242. /**
  243. * \fn adams_bashforth::adams_bashforth( const algebra_type &algebra )
  244. * \brief Constructs the adams_bashforth class. This constructor can be used as a default
  245. * constructor if the algebra has a default constructor.
  246. * \param algebra A copy of algebra is made and stored.
  247. */
  248. /**
  249. * \fn order_type adams_bashforth::order( void ) const
  250. * \brief Returns the order of the algorithm, which is equal to the number of steps.
  251. * \return order of the method.
  252. */
  253. /**
  254. * \fn void adams_bashforth::do_step( System system , StateInOut &x , time_type t , time_type dt )
  255. * \brief This method performs one step. It transforms the result in-place.
  256. *
  257. * \param system The system function to solve, hence the r.h.s. of the ordinary differential equation. It must fulfill the
  258. * Simple System concept.
  259. * \param x The state of the ODE which should be solved. After calling do_step the result is updated in x.
  260. * \param t The value of the time, at which the step should be performed.
  261. * \param dt The step size.
  262. */
  263. /**
  264. * \fn void adams_bashforth::do_step( System system , const StateIn &in , time_type t , StateOut &out , time_type dt )
  265. * \brief The method performs one step with the stepper passed by Stepper. The state of the ODE is updated out-of-place.
  266. *
  267. * \param system The system function to solve, hence the r.h.s. of the ODE. It must fulfill the
  268. * Simple System concept.
  269. * \param in The state of the ODE which should be solved. in is not modified in this method
  270. * \param t The value of the time, at which the step should be performed.
  271. * \param out The result of the step is written in out.
  272. * \param dt The step size.
  273. */
  274. /**
  275. * \fn void adams_bashforth::adjust_size( const StateType &x )
  276. * \brief Adjust the size of all temporaries in the stepper manually.
  277. * \param x A state from which the size of the temporaries to be resized is deduced.
  278. */
  279. /**
  280. * \fn const step_storage_type& adams_bashforth::step_storage( void ) const
  281. * \brief Returns the storage of intermediate results.
  282. * \return The storage of intermediate results.
  283. */
  284. /**
  285. * \fn step_storage_type& adams_bashforth::step_storage( void )
  286. * \brief Returns the storage of intermediate results.
  287. * \return The storage of intermediate results.
  288. */
  289. /**
  290. * \fn void adams_bashforth::initialize( ExplicitStepper explicit_stepper , System system , StateIn &x , time_type &t , time_type dt )
  291. * \brief Initialized the stepper. Does Steps-1 steps with the explicit_stepper to fill the buffer.
  292. * \param explicit_stepper the stepper used to fill the buffer of previous step results
  293. * \param system The system function to solve, hence the r.h.s. of the ordinary differential equation. It must fulfill the
  294. * Simple System concept.
  295. * \param x The state of the ODE which should be solved. After calling do_step the result is updated in x.
  296. * \param t The value of the time, at which the step should be performed.
  297. * \param dt The step size.
  298. */
  299. /**
  300. * \fn void adams_bashforth::initialize( System system , StateIn &x , time_type &t , time_type dt )
  301. * \brief Initialized the stepper. Does Steps-1 steps with an internal instance of InitializingStepper to fill the buffer.
  302. * \note The state x and time t are updated to the values after Steps-1 initial steps.
  303. * \param system The system function to solve, hence the r.h.s. of the ordinary differential equation. It must fulfill the
  304. * Simple System concept.
  305. * \param x The initial state of the ODE which should be solved, updated in this method.
  306. * \param t The initial value of the time, updated in this method.
  307. * \param dt The step size.
  308. */
  309. /**
  310. * \fn void adams_bashforth::reset( void )
  311. * \brief Resets the internal buffer of the stepper.
  312. */
  313. /**
  314. * \fn bool adams_bashforth::is_initialized( void ) const
  315. * \brief Returns true if the stepper has been initialized.
  316. * \return bool true if stepper is initialized, false otherwise
  317. */
  318. /**
  319. * \fn const initializing_stepper_type& adams_bashforth::initializing_stepper( void ) const
  320. * \brief Returns the internal initializing stepper instance.
  321. * \return initializing_stepper
  322. */
  323. /**
  324. * \fn const initializing_stepper_type& adams_bashforth::initializing_stepper( void ) const
  325. * \brief Returns the internal initializing stepper instance.
  326. * \return initializing_stepper
  327. */
  328. /**
  329. * \fn initializing_stepper_type& adams_bashforth::initializing_stepper( void )
  330. * \brief Returns the internal initializing stepper instance.
  331. * \return initializing_stepper
  332. */
  333. } // odeint
  334. } // numeric
  335. } // boost
  336. #endif // BOOST_NUMERIC_ODEINT_STEPPER_ADAMS_BASHFORTH_HPP_INCLUDED