comparison DEPENDENCIES/generic/include/boost/numeric/odeint/integrate/detail/integrate_const.hpp @ 101:c530137014c0

Update Boost headers (1.58.0)
author Chris Cannam
date Mon, 07 Sep 2015 11:12:49 +0100
parents 2665513ce2d3
children
comparison
equal deleted inserted replaced
100:793467b5e61c 101:c530137014c0
4 4
5 [begin_description] 5 [begin_description]
6 integrate const implementation 6 integrate const implementation
7 [end_description] 7 [end_description]
8 8
9 Copyright 2009-2012 Karsten Ahnert 9 Copyright 2012 Mario Mulansky
10 Copyright 2009-2012 Mario Mulansky 10 Copyright 2012 Christoph Koke
11 Copyright 2012 Karsten Ahnert
11 12
12 Distributed under the Boost Software License, Version 1.0. 13 Distributed under the Boost Software License, Version 1.0.
13 (See accompanying file LICENSE_1_0.txt or 14 (See accompanying file LICENSE_1_0.txt or
14 copy at http://www.boost.org/LICENSE_1_0.txt) 15 copy at http://www.boost.org/LICENSE_1_0.txt)
15 */ 16 */
44 Time start_time , Time end_time , Time dt , 45 Time start_time , Time end_time , Time dt ,
45 Observer observer , stepper_tag 46 Observer observer , stepper_tag
46 ) 47 )
47 { 48 {
48 typename odeint::unwrap_reference< Observer >::type &obs = observer; 49 typename odeint::unwrap_reference< Observer >::type &obs = observer;
50 typename odeint::unwrap_reference< Stepper >::type &st = stepper;
49 51
50 Time time = start_time; 52 Time time = start_time;
51 int step = 0; 53 int step = 0;
52 54 // cast time+dt explicitely in case of expression templates (e.g. multiprecision)
53 while( less_eq_with_sign( time+dt , end_time , dt ) ) 55 while( less_eq_with_sign( static_cast<Time>(time+dt) , end_time , dt ) )
54 { 56 {
55 obs( start_state , time ); 57 obs( start_state , time );
56 stepper.do_step( system , start_state , time , dt ); 58 st.do_step( system , start_state , time , dt );
57 // direct computation of the time avoids error propagation happening when using time += dt 59 // direct computation of the time avoids error propagation happening when using time += dt
58 // we need clumsy type analysis to get boost units working here 60 // we need clumsy type analysis to get boost units working here
59 ++step; 61 ++step;
60 time = start_time + static_cast< typename unit_value_type<Time>::type >(step) * dt; 62 time = start_time + static_cast< typename unit_value_type<Time>::type >(step) * dt;
61 } 63 }
75 { 77 {
76 typename odeint::unwrap_reference< Observer >::type &obs = observer; 78 typename odeint::unwrap_reference< Observer >::type &obs = observer;
77 79
78 Time time = start_time; 80 Time time = start_time;
79 const Time time_step = dt; 81 const Time time_step = dt;
82 int real_steps = 0;
80 int step = 0; 83 int step = 0;
81 84
82 while( less_eq_with_sign( time+time_step , end_time , dt ) ) 85 while( less_eq_with_sign( static_cast<Time>(time+time_step) , end_time , dt ) )
83 { 86 {
84 obs( start_state , time ); 87 obs( start_state , time );
85 detail::integrate_adaptive( stepper , system , start_state , time , time+time_step , dt , 88 real_steps += detail::integrate_adaptive( stepper , system , start_state , time , time+time_step , dt ,
86 null_observer() , controlled_stepper_tag() ); 89 null_observer() , controlled_stepper_tag() );
87 // direct computation of the time avoids error propagation happening when using time += dt 90 // direct computation of the time avoids error propagation happening when using time += dt
88 // we need clumsy type analysis to get boost units working here 91 // we need clumsy type analysis to get boost units working here
89 ++step; 92 step++;
90 time = start_time + static_cast< typename unit_value_type<Time>::type >(step) * time_step; 93 time = start_time + static_cast< typename unit_value_type<Time>::type >(step) * time_step;
91 } 94 }
92 obs( start_state , time ); 95 obs( start_state , time );
93 96
94 return step; 97 return real_steps;
95 } 98 }
96 99
97 100
98 template< class Stepper , class System , class State , class Time , class Observer > 101 template< class Stepper , class System , class State , class Time , class Observer >
99 size_t integrate_const( 102 size_t integrate_const(
101 Time start_time , Time end_time , Time dt , 104 Time start_time , Time end_time , Time dt ,
102 Observer observer , dense_output_stepper_tag 105 Observer observer , dense_output_stepper_tag
103 ) 106 )
104 { 107 {
105 typename odeint::unwrap_reference< Observer >::type &obs = observer; 108 typename odeint::unwrap_reference< Observer >::type &obs = observer;
109 typename odeint::unwrap_reference< Stepper >::type &st = stepper;
106 110
107 Time time = start_time; 111 Time time = start_time;
108 112
109 stepper.initialize( start_state , time , dt ); 113 st.initialize( start_state , time , dt );
110 obs( start_state , time ); 114 obs( start_state , time );
111 time += dt; 115 time += dt;
112 116
113 int obs_step( 1 ); 117 int obs_step( 1 );
114 int real_step( 0 ); 118 int real_step( 0 );
115 119
116 while( less_with_sign( time+dt , end_time , dt ) ) 120 while( less_with_sign( static_cast<Time>(time+dt) , end_time , dt ) )
117 { 121 {
118 while( less_eq_with_sign( time , stepper.current_time() , dt ) ) 122 while( less_eq_with_sign( time , st.current_time() , dt ) )
119 { 123 {
120 stepper.calc_state( time , start_state ); 124 st.calc_state( time , start_state );
121 obs( start_state , time ); 125 obs( start_state , time );
122 ++obs_step; 126 ++obs_step;
123 // direct computation of the time avoids error propagation happening when using time += dt 127 // direct computation of the time avoids error propagation happening when using time += dt
124 // we need clumsy type analysis to get boost units working here 128 // we need clumsy type analysis to get boost units working here
125 time = start_time + static_cast< typename unit_value_type<Time>::type >(obs_step) * dt; 129 time = start_time + static_cast< typename unit_value_type<Time>::type >(obs_step) * dt;
126 } 130 }
127 // we have not reached the end, do another real step 131 // we have not reached the end, do another real step
128 if( less_with_sign( stepper.current_time()+stepper.current_time_step() , 132 if( less_with_sign( static_cast<Time>(st.current_time()+st.current_time_step()) ,
129 end_time , 133 end_time ,
130 stepper.current_time_step() ) ) 134 st.current_time_step() ) )
131 { 135 {
132 while( less_eq_with_sign( stepper.current_time() , time , dt ) ) 136 while( less_eq_with_sign( st.current_time() , time , dt ) )
133 { 137 {
134 stepper.do_step( system ); 138 st.do_step( system );
135 ++real_step; 139 ++real_step;
136 } 140 }
137 } 141 }
138 else if( less_with_sign( stepper.current_time() , end_time , stepper.current_time_step() ) ) 142 else if( less_with_sign( st.current_time() , end_time , st.current_time_step() ) )
139 { // do the last step ending exactly on the end point 143 { // do the last step ending exactly on the end point
140 stepper.initialize( stepper.current_state() , stepper.current_time() , end_time - stepper.current_time() ); 144 st.initialize( st.current_state() , st.current_time() , end_time - st.current_time() );
141 stepper.do_step( system ); 145 st.do_step( system );
142 ++real_step; 146 ++real_step;
143 } 147 }
144 148
145 } 149 }
146 // last observation, if we are still in observation interval 150 // last observation, if we are still in observation interval
147 if( less_eq_with_sign( time , end_time , dt ) ) 151 if( less_eq_with_sign( time , end_time , dt ) )
148 { 152 {
149 stepper.calc_state( time , start_state ); 153 st.calc_state( time , start_state );
150 obs( start_state , time ); 154 obs( start_state , time );
151 } 155 }
152 156
153 return real_step; 157 return real_step;
154 } 158 }