annotate DEPENDENCIES/generic/include/boost/numeric/odeint/integrate/detail/integrate_times.hpp @ 125:34e428693f5d vext

Vext -> Repoint
author Chris Cannam
date Thu, 14 Jun 2018 11:15:39 +0100
parents c530137014c0
children
rev   line source
Chris@16 1 /*
Chris@16 2 [auto_generated]
Chris@16 3 boost/numeric/odeint/integrate/detail/integrate_times.hpp
Chris@16 4
Chris@16 5 [begin_description]
Chris@16 6 Default integrate times implementation.
Chris@16 7 [end_description]
Chris@16 8
Chris@101 9 Copyright 2011-2012 Mario Mulansky
Chris@101 10 Copyright 2012 Karsten Ahnert
Chris@101 11 Copyright 2012 Christoph Koke
Chris@16 12
Chris@16 13 Distributed under the Boost Software License, Version 1.0.
Chris@16 14 (See accompanying file LICENSE_1_0.txt or
Chris@16 15 copy at http://www.boost.org/LICENSE_1_0.txt)
Chris@16 16 */
Chris@16 17
Chris@16 18
Chris@16 19 #ifndef BOOST_NUMERIC_ODEINT_INTEGRATE_DETAIL_INTEGRATE_TIMES_HPP_INCLUDED
Chris@16 20 #define BOOST_NUMERIC_ODEINT_INTEGRATE_DETAIL_INTEGRATE_TIMES_HPP_INCLUDED
Chris@16 21
Chris@16 22 #include <stdexcept>
Chris@16 23
Chris@16 24 #include <boost/config.hpp>
Chris@101 25 #include <boost/throw_exception.hpp>
Chris@16 26 #include <boost/numeric/odeint/util/unwrap_reference.hpp>
Chris@16 27 #include <boost/numeric/odeint/stepper/controlled_step_result.hpp>
Chris@16 28 #include <boost/numeric/odeint/util/detail/less_with_sign.hpp>
Chris@16 29
Chris@16 30
Chris@16 31 namespace boost {
Chris@16 32 namespace numeric {
Chris@16 33 namespace odeint {
Chris@16 34 namespace detail {
Chris@16 35
Chris@16 36
Chris@16 37
Chris@16 38 /*
Chris@16 39 * integrate_times for simple stepper
Chris@16 40 */
Chris@16 41 template< class Stepper , class System , class State , class TimeIterator , class Time , class Observer >
Chris@16 42 size_t integrate_times(
Chris@16 43 Stepper stepper , System system , State &start_state ,
Chris@16 44 TimeIterator start_time , TimeIterator end_time , Time dt ,
Chris@16 45 Observer observer , stepper_tag
Chris@16 46 )
Chris@16 47 {
Chris@16 48 typename odeint::unwrap_reference< Observer >::type &obs = observer;
Chris@101 49 typename odeint::unwrap_reference< Stepper >::type &st = stepper;
Chris@101 50 typedef typename unit_value_type<Time>::type time_type;
Chris@16 51
Chris@16 52 size_t steps = 0;
Chris@16 53 Time current_dt = dt;
Chris@16 54 while( true )
Chris@16 55 {
Chris@16 56 Time current_time = *start_time++;
Chris@16 57 obs( start_state , current_time );
Chris@16 58 if( start_time == end_time )
Chris@16 59 break;
Chris@101 60 while( less_with_sign( current_time , static_cast<time_type>(*start_time) , current_dt ) )
Chris@16 61 {
Chris@101 62 current_dt = min_abs( dt , *start_time - current_time );
Chris@101 63 st.do_step( system , start_state , current_time , current_dt );
Chris@16 64 current_time += current_dt;
Chris@16 65 steps++;
Chris@16 66 }
Chris@16 67 }
Chris@16 68 return steps;
Chris@16 69 }
Chris@16 70
Chris@16 71 /*
Chris@16 72 * integrate_times for controlled stepper
Chris@16 73 */
Chris@16 74 template< class Stepper , class System , class State , class TimeIterator , class Time , class Observer >
Chris@16 75 size_t integrate_times(
Chris@16 76 Stepper stepper , System system , State &start_state ,
Chris@16 77 TimeIterator start_time , TimeIterator end_time , Time dt ,
Chris@16 78 Observer observer , controlled_stepper_tag
Chris@16 79 )
Chris@16 80 {
Chris@16 81 typename odeint::unwrap_reference< Observer >::type &obs = observer;
Chris@101 82 typename odeint::unwrap_reference< Stepper >::type &st = stepper;
Chris@101 83 typedef typename unit_value_type<Time>::type time_type;
Chris@16 84
Chris@16 85 const size_t max_attempts = 1000;
Chris@16 86 const char *error_string = "Integrate adaptive : Maximal number of iterations reached. A step size could not be found.";
Chris@16 87 size_t steps = 0;
Chris@16 88 while( true )
Chris@16 89 {
Chris@16 90 size_t fail_steps = 0;
Chris@16 91 Time current_time = *start_time++;
Chris@16 92 obs( start_state , current_time );
Chris@16 93 if( start_time == end_time )
Chris@16 94 break;
Chris@101 95 while( less_with_sign( current_time , static_cast<time_type>(*start_time) , dt ) )
Chris@16 96 {
Chris@101 97 // adjust stepsize to end up exactly at the observation point
Chris@101 98 Time current_dt = min_abs( dt , *start_time - current_time );
Chris@101 99 if( st.try_step( system , start_state , current_time , current_dt ) == success )
Chris@16 100 {
Chris@16 101 ++steps;
Chris@101 102 // continue with the original step size if dt was reduced due to observation
Chris@101 103 dt = max_abs( dt , current_dt );
Chris@16 104 }
Chris@16 105 else
Chris@16 106 {
Chris@16 107 ++fail_steps;
Chris@101 108 dt = current_dt;
Chris@16 109 }
Chris@101 110 if( fail_steps == max_attempts ) BOOST_THROW_EXCEPTION( std::overflow_error( error_string ));
Chris@16 111 }
Chris@16 112 }
Chris@16 113 return steps;
Chris@16 114 }
Chris@16 115
Chris@16 116 /*
Chris@16 117 * integrate_times for dense output stepper
Chris@16 118 */
Chris@16 119 template< class Stepper , class System , class State , class TimeIterator , class Time , class Observer >
Chris@16 120 size_t integrate_times(
Chris@16 121 Stepper stepper , System system , State &start_state ,
Chris@16 122 TimeIterator start_time , TimeIterator end_time , Time dt ,
Chris@16 123 Observer observer , dense_output_stepper_tag
Chris@16 124 )
Chris@16 125 {
Chris@16 126 typename odeint::unwrap_reference< Observer >::type &obs = observer;
Chris@101 127 typename odeint::unwrap_reference< Stepper >::type &st = stepper;
Chris@101 128 typedef typename unit_value_type<Time>::type time_type;
Chris@16 129
Chris@16 130 if( start_time == end_time )
Chris@16 131 return 0;
Chris@16 132
Chris@101 133 Time last_time_point = static_cast<time_type>(*(end_time-1));
Chris@16 134
Chris@101 135 st.initialize( start_state , *start_time , dt );
Chris@16 136 obs( start_state , *start_time++ );
Chris@16 137
Chris@16 138 size_t count = 0;
Chris@16 139 while( start_time != end_time )
Chris@16 140 {
Chris@101 141 while( ( start_time != end_time ) && less_eq_with_sign( static_cast<time_type>(*start_time) , st.current_time() , st.current_time_step() ) )
Chris@16 142 {
Chris@101 143 st.calc_state( *start_time , start_state );
Chris@16 144 obs( start_state , *start_time );
Chris@16 145 start_time++;
Chris@16 146 }
Chris@16 147
Chris@16 148 // we have not reached the end, do another real step
Chris@101 149 if( less_eq_with_sign( st.current_time() + st.current_time_step() ,
Chris@16 150 last_time_point ,
Chris@101 151 st.current_time_step() ) )
Chris@16 152 {
Chris@101 153 st.do_step( system );
Chris@16 154 ++count;
Chris@16 155 }
Chris@16 156 else if( start_time != end_time )
Chris@16 157 { // do the last step ending exactly on the end point
Chris@101 158 st.initialize( st.current_state() , st.current_time() , last_time_point - st.current_time() );
Chris@101 159 st.do_step( system );
Chris@16 160 ++count;
Chris@16 161 }
Chris@16 162 }
Chris@16 163 return count;
Chris@16 164 }
Chris@16 165
Chris@16 166
Chris@16 167 } // namespace detail
Chris@16 168 } // namespace odeint
Chris@16 169 } // namespace numeric
Chris@16 170 } // namespace boost
Chris@16 171
Chris@16 172
Chris@16 173 #endif // BOOST_NUMERIC_ODEINT_INTEGRATE_DETAIL_INTEGRATE_ADAPTIVE_HPP_INCLUDED