3 boost/numeric/odeint/stepper/adams_moulton.hpp
6 Implementation of the Adams-Moulton method. This is method is not a real stepper, it is more a helper class
7 which computes the corrector step in the Adams-Bashforth-Moulton method.
10 Copyright 2011-2012 Karsten Ahnert
11 Copyright 2011-2013 Mario Mulansky
12 Copyright 2012 Christoph Koke
14 Distributed under the Boost Software License, Version 1.0.
15 (See accompanying file LICENSE_1_0.txt or
16 copy at http://www.boost.org/LICENSE_1_0.txt)
20 #ifndef BOOST_NUMERIC_ODEINT_STEPPER_ADAMS_MOULTON_HPP_INCLUDED
21 #define BOOST_NUMERIC_ODEINT_STEPPER_ADAMS_MOULTON_HPP_INCLUDED
24 #include <boost/numeric/odeint/util/bind.hpp>
26 #include <boost/numeric/odeint/algebra/range_algebra.hpp>
27 #include <boost/numeric/odeint/algebra/default_operations.hpp>
28 #include <boost/numeric/odeint/algebra/algebra_dispatcher.hpp>
29 #include <boost/numeric/odeint/algebra/operations_dispatcher.hpp>
31 #include <boost/numeric/odeint/util/state_wrapper.hpp>
32 #include <boost/numeric/odeint/util/is_resizeable.hpp>
33 #include <boost/numeric/odeint/util/resizer.hpp>
35 #include <boost/numeric/odeint/stepper/stepper_categories.hpp>
36 #include <boost/numeric/odeint/stepper/runge_kutta4_classic.hpp>
38 #include <boost/numeric/odeint/stepper/detail/adams_moulton_call_algebra.hpp>
39 #include <boost/numeric/odeint/stepper/detail/adams_moulton_coefficients.hpp>
40 #include <boost/numeric/odeint/stepper/detail/rotating_buffer.hpp>
51 * Static implicit Adams-Moulton multistep-solver without step size control and without dense output.
56 class Value = double ,
59 class Algebra = typename algebra_dispatcher< State >::algebra_type ,
60 class Operations = typename operations_dispatcher< State >::operations_type ,
61 class Resizer = initially_resizer
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 Algebra algebra_type;
77 typedef Operations operations_type;
78 typedef Resizer resizer_type;
79 typedef stepper_tag stepper_category;
81 typedef adams_moulton< Steps , State , Value , Deriv , Time , Algebra , Operations , Resizer > stepper_type;
83 static const size_t steps = Steps;
85 typedef unsigned short order_type;
86 static const order_type order_value = steps + 1;
88 typedef detail::rotating_buffer< wrapped_deriv_type , steps > step_storage_type;
91 : m_coefficients() , m_dxdt() , m_resizer() ,
92 m_algebra_instance() , m_algebra( m_algebra_instance )
95 adams_moulton( algebra_type &algebra )
96 : m_coefficients() , m_dxdt() , m_resizer() ,
97 m_algebra_instance() , m_algebra( algebra )
100 adams_moulton& operator=( const adams_moulton &stepper )
102 m_dxdt = stepper.m_dxdt;
103 m_resizer = stepper.m_resizer;
104 m_algebra = stepper.m_algebra;
108 order_type order( void ) const { return order_value; }
112 * Version 1 : do_step( system , x , t , dt , buf );
114 * solves the forwarding problem
116 template< class System , class StateInOut , class StateIn , class ABBuf >
117 void do_step( System system , StateInOut &x , StateIn const & pred , time_type t , time_type dt , const ABBuf &buf )
119 do_step( system , x , pred , t , x , dt , buf );
122 template< class System , class StateInOut , class StateIn , class ABBuf >
123 void do_step( System system , const StateInOut &x , StateIn const & pred , time_type t , time_type dt , const ABBuf &buf )
125 do_step( system , x , pred , t , x , dt , buf );
131 * Version 2 : do_step( system , in , t , out , dt , buf );
133 * solves the forwarding problem
135 template< class System , class StateIn , class PredIn , class StateOut , class ABBuf >
136 void do_step( System system , const StateIn &in , const PredIn &pred , time_type t , StateOut &out , time_type dt , const ABBuf &buf )
138 do_step_impl( system , in , pred , t , out , dt , buf );
141 template< class System , class StateIn , class PredIn , class StateOut , class ABBuf >
142 void do_step( System system , const StateIn &in , const PredIn &pred , time_type t , const StateOut &out , time_type dt , const ABBuf &buf )
144 do_step_impl( system , in , pred , t , out , dt , buf );
149 template< class StateType >
150 void adjust_size( const StateType &x )
155 algebra_type& algebra()
156 { return m_algebra; }
158 const algebra_type& algebra() const
159 { return m_algebra; }
165 template< class System , class StateIn , class PredIn , class StateOut , class ABBuf >
166 void do_step_impl( System system , const StateIn &in , const PredIn &pred , time_type t , StateOut &out , time_type dt , const ABBuf &buf )
168 typename odeint::unwrap_reference< System >::type &sys = system;
169 m_resizer.adjust_size( in , detail::bind( &stepper_type::template resize_impl<StateIn> , detail::ref( *this ) , detail::_1 ) );
170 sys( pred , m_dxdt.m_v , t );
171 detail::adams_moulton_call_algebra< steps , algebra_type , operations_type >()( m_algebra , in , out , m_dxdt.m_v , buf , m_coefficients , dt );
175 template< class StateIn >
176 bool resize_impl( const StateIn &x )
178 return adjust_size_by_resizeability( m_dxdt , x , typename is_resizeable<deriv_type>::type() );
182 const detail::adams_moulton_coefficients< value_type , steps > m_coefficients;
183 wrapped_deriv_type m_dxdt;
184 resizer_type m_resizer;
188 algebra_type m_algebra_instance;
189 algebra_type &m_algebra;
201 #endif // BOOST_NUMERIC_ODEINT_STEPPER_ADAMS_MOULTON_HPP_INCLUDED