3 boost/numeric/odeint/integrate/detail/integrate_adaptive.hpp
6 Default Integrate adaptive implementation.
9 Copyright 2011-2013 Karsten Ahnert
10 Copyright 2011-2015 Mario Mulansky
11 Copyright 2012 Christoph Koke
13 Distributed under the Boost Software License, Version 1.0.
14 (See accompanying file LICENSE_1_0.txt or
15 copy at http://www.boost.org/LICENSE_1_0.txt)
19 #ifndef BOOST_NUMERIC_ODEINT_INTEGRATE_DETAIL_INTEGRATE_ADAPTIVE_HPP_INCLUDED
20 #define BOOST_NUMERIC_ODEINT_INTEGRATE_DETAIL_INTEGRATE_ADAPTIVE_HPP_INCLUDED
24 #include <boost/throw_exception.hpp>
26 #include <boost/numeric/odeint/stepper/stepper_categories.hpp>
27 #include <boost/numeric/odeint/stepper/controlled_step_result.hpp>
28 #include <boost/numeric/odeint/integrate/max_step_checker.hpp>
29 #include <boost/numeric/odeint/integrate/detail/integrate_const.hpp>
30 #include <boost/numeric/odeint/util/bind.hpp>
31 #include <boost/numeric/odeint/util/unwrap_reference.hpp>
32 #include <boost/numeric/odeint/util/copy.hpp>
34 #include <boost/numeric/odeint/util/detail/less_with_sign.hpp>
44 // forward declaration
45 template< class Stepper , class System , class State , class Time , class Observer >
46 size_t integrate_const(
47 Stepper stepper , System system , State &start_state ,
48 Time start_time , Time end_time , Time dt ,
49 Observer observer , stepper_tag );
52 * integrate_adaptive for simple stepper is basically an integrate_const + some last step
54 template< class Stepper , class System , class State , class Time , class Observer >
55 size_t integrate_adaptive(
56 Stepper stepper , System system , State &start_state ,
57 Time start_time , Time end_time , Time dt ,
58 Observer observer , stepper_tag
61 size_t steps = detail::integrate_const( stepper , system , start_state , start_time ,
62 end_time , dt , observer , stepper_tag() );
63 typename odeint::unwrap_reference< Observer >::type &obs = observer;
64 typename odeint::unwrap_reference< Stepper >::type &st = stepper;
66 Time end = start_time + dt*steps;
67 if( less_with_sign( end , end_time , dt ) )
68 { //make a last step to end exactly at end_time
69 st.do_step( system , start_state , end , end_time - end );
71 obs( start_state , end_time );
78 * integrate adaptive for controlled stepper
80 template< class Stepper , class System , class State , class Time , class Observer >
81 size_t integrate_adaptive(
82 Stepper stepper , System system , State &start_state ,
83 Time &start_time , Time end_time , Time &dt ,
84 Observer observer , controlled_stepper_tag
87 typename odeint::unwrap_reference< Observer >::type &obs = observer;
88 typename odeint::unwrap_reference< Stepper >::type &st = stepper;
90 failed_step_checker fail_checker; // to throw a runtime_error if step size adjustment fails
92 while( less_with_sign( start_time , end_time , dt ) )
94 obs( start_state , start_time );
95 if( less_with_sign( end_time , static_cast<Time>(start_time + dt) , dt ) )
97 dt = end_time - start_time;
100 controlled_step_result res;
103 res = st.try_step( system , start_state , start_time , dt );
104 fail_checker(); // check number of failed steps
106 while( res == fail );
107 fail_checker.reset(); // if we reach here, the step was successful -> reset fail checker
111 obs( start_state , start_time );
117 * integrate adaptive for dense output steppers
119 * step size control is used if the stepper supports it
121 template< class Stepper , class System , class State , class Time , class Observer >
122 size_t integrate_adaptive(
123 Stepper stepper , System system , State &start_state ,
124 Time start_time , Time end_time , Time dt ,
125 Observer observer , dense_output_stepper_tag )
127 typename odeint::unwrap_reference< Observer >::type &obs = observer;
128 typename odeint::unwrap_reference< Stepper >::type &st = stepper;
131 st.initialize( start_state , start_time , dt );
133 while( less_with_sign( st.current_time() , end_time , st.current_time_step() ) )
135 while( less_eq_with_sign( static_cast<Time>(st.current_time() + st.current_time_step()) ,
137 st.current_time_step() ) )
138 { //make sure we don't go beyond the end_time
139 obs( st.current_state() , st.current_time() );
140 st.do_step( system );
143 // calculate time step to arrive exactly at end time
144 st.initialize( st.current_state() , st.current_time() , static_cast<Time>(end_time - st.current_time()) );
146 obs( st.current_state() , st.current_time() );
147 // overwrite start_state with the final point
148 boost::numeric::odeint::copy( st.current_state() , start_state );
155 } // namespace detail
156 } // namespace odeint
157 } // namespace numeric
161 #endif // BOOST_NUMERIC_ODEINT_INTEGRATE_DETAIL_INTEGRATE_ADAPTIVE_HPP_INCLUDED