ceph/src/boost/boost/numeric/odeint/integrate/detail/integrate_n_steps.hpp

   1 /*
   2  [auto_generated]
   3  boost/numeric/odeint/integrate/detail/integrate_n_steps.hpp
   4
   5  [begin_description]
   6  integrate steps implementation
   7  [end_description]
   8
   9  Copyright 2012-2015 Mario Mulansky
  10  Copyright 2012 Christoph Koke
  11  Copyright 2012 Karsten Ahnert
  12
  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)
  16  */
  17
  18 #ifndef BOOST_NUMERIC_ODEINT_INTEGRATE_DETAIL_INTEGRATE_N_STEPS_HPP_INCLUDED
  19 #define BOOST_NUMERIC_ODEINT_INTEGRATE_DETAIL_INTEGRATE_N_STEPS_HPP_INCLUDED
  20
  21 #include <boost/numeric/odeint/util/unwrap_reference.hpp>
  22 #include <boost/numeric/odeint/stepper/stepper_categories.hpp>
  23 #include <boost/numeric/odeint/integrate/detail/integrate_adaptive.hpp>
  24 #include <boost/numeric/odeint/util/unit_helper.hpp>
  25
  26 #include <boost/numeric/odeint/util/detail/less_with_sign.hpp>
  27
  28 namespace boost {
  29 namespace numeric {
  30 namespace odeint {
  31 namespace detail {
  32
  33 // forward declaration
  34 template< class Stepper , class System , class State , class Time , class Observer >
  35 size_t integrate_adaptive_checked(
  36         Stepper stepper , System system , State &start_state ,
  37         Time &start_time , Time end_time , Time &dt ,
  38         Observer observer, controlled_stepper_tag
  39 );
  40
  41
  42 /* basic version */
  43 template< class Stepper , class System , class State , class Time , class Observer>
  44 Time integrate_n_steps(
  45         Stepper stepper , System system , State &start_state ,
  46         Time start_time , Time dt , size_t num_of_steps ,
  47         Observer observer , stepper_tag )
  48 {
  49     typename odeint::unwrap_reference< Observer >::type &obs = observer;
  50     typename odeint::unwrap_reference< Stepper >::type &st = stepper;
  51
  52     Time time = start_time;
  53
  54     for( size_t step = 0; step < num_of_steps ; ++step )
  55     {
  56         obs( start_state , time );
  57         st.do_step( system , start_state , time , dt );
  58         // direct computation of the time avoids error propagation happening when using time += dt
  59         // we need clumsy type analysis to get boost units working here
  60         time = start_time + static_cast< typename unit_value_type<Time>::type >( step+1 ) * dt;
  61     }
  62     obs( start_state , time );
  63
  64     return time;
  65 }
  66
  67
  68 /* controlled version */
  69 template< class Stepper , class System , class State , class Time , class Observer >
  70 Time integrate_n_steps(
  71         Stepper stepper , System system , State &start_state ,
  72         Time start_time , Time dt , size_t num_of_steps ,
  73         Observer observer , controlled_stepper_tag )
  74 {
  75     typename odeint::unwrap_reference< Observer >::type &obs = observer;
  76
  77     Time time = start_time;
  78     Time time_step = dt;
  79
  80     for( size_t step = 0; step < num_of_steps ; ++step )
  81     {
  82         obs( start_state , time );
  83         // integrate_adaptive_checked uses the given checker to throw if an overflow occurs
  84         detail::integrate_adaptive(stepper, system, start_state, time, static_cast<Time>(time + time_step), dt,
  85                                    null_observer(), controlled_stepper_tag());
  86         // direct computation of the time avoids error propagation happening when using time += dt
  87         // we need clumsy type analysis to get boost units working here
  88         time = start_time + static_cast< typename unit_value_type<Time>::type >(step+1) * time_step;
  89     }
  90     obs( start_state , time );
  91
  92     return time;
  93 }
  94
  95
  96 /* dense output version */
  97 template< class Stepper , class System , class State , class Time , class Observer >
  98 Time integrate_n_steps(
  99         Stepper stepper , System system , State &start_state ,
 100         Time start_time , Time dt , size_t num_of_steps ,
 101         Observer observer , dense_output_stepper_tag )
 102 {
 103     typename odeint::unwrap_reference< Observer >::type &obs = observer;
 104     typename odeint::unwrap_reference< Stepper >::type &st = stepper;
 105
 106     Time time = start_time;
 107     const Time end_time = start_time + static_cast< typename unit_value_type<Time>::type >(num_of_steps) * dt;
 108
 109     st.initialize( start_state , time , dt );
 110
 111     size_t step = 0;
 112
 113     while( step < num_of_steps )
 114     {
 115         while( less_with_sign( time , st.current_time() , st.current_time_step() ) )
 116         {
 117             st.calc_state( time , start_state );
 118             obs( start_state , time );
 119             ++step;
 120             // direct computation of the time avoids error propagation happening when using time += dt
 121             // we need clumsy type analysis to get boost units working here
 122             time = start_time + static_cast< typename unit_value_type<Time>::type >(step) * dt;
 123         }
 124
 125         // we have not reached the end, do another real step
 126         if( less_with_sign( static_cast<Time>(st.current_time()+st.current_time_step()) ,
 127                             end_time ,
 128                             st.current_time_step() ) )
 129         {
 130             st.do_step( system );
 131         }
 132         else if( less_with_sign( st.current_time() , end_time , st.current_time_step() ) )
 133         { // do the last step ending exactly on the end point
 134             st.initialize( st.current_state() , st.current_time() , static_cast<Time>(end_time - st.current_time()) );
 135             st.do_step( system );
 136         }
 137     }
 138
 139     // make sure we really end exactly where we should end
 140     while( st.current_time() < end_time )
 141     {
 142         if( less_with_sign( end_time ,
 143                             static_cast<Time>(st.current_time()+st.current_time_step()) ,
 144                             st.current_time_step() ) )
 145             st.initialize( st.current_state() , st.current_time() , static_cast<Time>(end_time - st.current_time()) );
 146         st.do_step( system );
 147     }
 148
 149     // observation at final point
 150     obs( st.current_state() , end_time );
 151
 152     return time;
 153 }
 154
 155
 156 }
 157 }
 158 }
 159 }
 160
 161 #endif /* BOOST_NUMERIC_ODEINT_INTEGRATE_DETAIL_INTEGRATE_N_STEPS_HPP_INCLUDED */