1 /*
2 [auto_generated]
3 libs/numeric/odeint/test/stepper_with_units.cpp
4
5 [begin_description]
6 This file tests if the steppers play well with Boost.Units.
7 [end_description]
8
9 Copyright 2011-2012 Karsten Ahnert
10 Copyright 2011-2013 Mario Mulansky
11
12 Distributed under the Boost Software License, Version 1.0.
13 (See accompanying file LICENSE_1_0.txt or
14 copy at http://www.boost.org/LICENSE_1_0.txt)
15 */
16
17 #define BOOST_TEST_MODULE odeint_stepper_with_units
18
19 // the runge-kutta 78 stepper invoked with boost units requires increased fusion macro variables!
20 // note that by default the rk78 + units test case is disabled as it requires enormous memory when compiling (5 GB)
21 #define BOOST_FUSION_INVOKE_MAX_ARITY 15
22 #define BOOST_RESULT_OF_NUM_ARGS 15
23
24
25 #include <boost/numeric/odeint/config.hpp>
26
27 #include <boost/test/unit_test.hpp>
28
29 #include <boost/units/systems/si/length.hpp>
30 #include <boost/units/systems/si/time.hpp>
31 #include <boost/units/systems/si/velocity.hpp>
32 #include <boost/units/systems/si/acceleration.hpp>
33 #include <boost/units/systems/si/io.hpp>
34
35 #include <boost/fusion/container.hpp>
36 #include <boost/mpl/vector.hpp>
37
38 #include <boost/numeric/odeint/stepper/euler.hpp>
39 #include <boost/numeric/odeint/stepper/runge_kutta4.hpp>
40 #include <boost/numeric/odeint/stepper/runge_kutta4_classic.hpp>
41 #include <boost/numeric/odeint/stepper/runge_kutta_cash_karp54.hpp>
42 #include <boost/numeric/odeint/stepper/runge_kutta_cash_karp54_classic.hpp>
43 #include <boost/numeric/odeint/stepper/runge_kutta_dopri5.hpp>
44 #include <boost/numeric/odeint/stepper/runge_kutta_fehlberg78.hpp>
45 #include <boost/numeric/odeint/stepper/controlled_runge_kutta.hpp>
46 #include <boost/numeric/odeint/stepper/dense_output_runge_kutta.hpp>
47 #include <boost/numeric/odeint/stepper/bulirsch_stoer.hpp>
48 #include <boost/numeric/odeint/stepper/bulirsch_stoer_dense_out.hpp>
49 #include <boost/numeric/odeint/algebra/fusion_algebra.hpp>
50 #include <boost/numeric/odeint/algebra/fusion_algebra_dispatcher.hpp>
51
52
53 using namespace boost::numeric::odeint;
54 using namespace boost::unit_test;
55 namespace mpl = boost::mpl;
56 namespace fusion = boost::fusion;
57 namespace units = boost::units;
58 namespace si = boost::units::si;
59
60 typedef double value_type;
61 typedef units::quantity< si::time , value_type > time_type;
62 typedef units::quantity< si::length , value_type > length_type;
63 typedef units::quantity< si::velocity , value_type > velocity_type;
64 typedef units::quantity< si::acceleration , value_type > acceleration_type;
65 typedef fusion::vector< length_type , velocity_type > state_type;
66 typedef fusion::vector< velocity_type , acceleration_type > deriv_type;
67
oscillator(const state_type & x,deriv_type & dxdt,time_type t)68 void oscillator( const state_type &x , deriv_type &dxdt , time_type t )
69 {
70 const units::quantity< si::frequency , value_type > omega = 1.0 * si::hertz;
71 fusion::at_c< 0 >( dxdt ) = fusion::at_c< 1 >( x );
72 fusion::at_c< 1 >( dxdt ) = - omega * omega * fusion::at_c< 0 >( x );
73 }
74
75 template< class Stepper >
check_stepper(Stepper & stepper)76 void check_stepper( Stepper &stepper )
77 {
78 typedef Stepper stepper_type;
79 typedef typename stepper_type::state_type state_type;
80 typedef typename stepper_type::value_type value_type;
81 typedef typename stepper_type::deriv_type deriv_type;
82 typedef typename stepper_type::time_type time_type;
83 typedef typename stepper_type::order_type order_type;
84 typedef typename stepper_type::algebra_type algebra_type;
85 typedef typename stepper_type::operations_type operations_type;
86
87 const time_type t( 0.0 * si::second );
88 time_type dt( 0.1 * si::second );
89 state_type x( 1.0 * si::meter , 0.0 * si::meter_per_second );
90
91 // test call method one
92 stepper.do_step( oscillator , x , t , dt );
93
94 // test call method two
95 stepper.do_step( oscillator , x , t , x , dt );
96
97 // test call method three
98 deriv_type dxdt;
99 oscillator( x , dxdt , t );
100 stepper.do_step( oscillator , x , dxdt , t , dt );
101
102 // test call method four
103 oscillator( x , dxdt , t );
104 stepper.do_step( oscillator , x , dxdt , t , x , dt );
105 }
106
107 template< class Stepper >
check_fsal_stepper(Stepper & stepper)108 void check_fsal_stepper( Stepper &stepper )
109 {
110 typedef Stepper stepper_type;
111 typedef typename stepper_type::state_type state_type;
112 typedef typename stepper_type::value_type value_type;
113 typedef typename stepper_type::deriv_type deriv_type;
114 typedef typename stepper_type::time_type time_type;
115 typedef typename stepper_type::order_type order_type;
116 typedef typename stepper_type::algebra_type algebra_type;
117 typedef typename stepper_type::operations_type operations_type;
118
119 const time_type t( 0.0 * si::second );
120 time_type dt( 0.1 * si::second );
121 state_type x( 1.0 * si::meter , 0.0 * si::meter_per_second );
122
123 // test call method one
124 stepper.do_step( oscillator , x , t , dt );
125
126 // test call method two
127 stepper.do_step( oscillator , x , t , x , dt );
128
129 // test call method three
130 deriv_type dxdt;
131 oscillator( x , dxdt , t );
132 stepper.do_step( oscillator , x , dxdt , t , dt );
133
134 // test call method four
135 stepper.do_step( oscillator , x , dxdt , t , x , dxdt , dt );
136 }
137
138 template< class Stepper >
check_error_stepper(Stepper & stepper)139 void check_error_stepper( Stepper &stepper )
140 {
141 typedef Stepper stepper_type;
142 typedef typename stepper_type::state_type state_type;
143 typedef typename stepper_type::value_type value_type;
144 typedef typename stepper_type::deriv_type deriv_type;
145 typedef typename stepper_type::time_type time_type;
146 typedef typename stepper_type::order_type order_type;
147 typedef typename stepper_type::algebra_type algebra_type;
148 typedef typename stepper_type::operations_type operations_type;
149
150 const time_type t( 0.0 * si::second );
151 time_type dt( 0.1 * si::second );
152 state_type x( 1.0 * si::meter , 0.0 * si::meter_per_second ) , xerr;
153
154 // test call method one
155 stepper.do_step( oscillator , x , t , dt , xerr );
156
157 // test call method two
158 stepper.do_step( oscillator , x , t , x , dt , xerr );
159
160 // test call method three
161 deriv_type dxdt;
162 oscillator( x , dxdt , t );
163 stepper.do_step( oscillator , x , dxdt , t , dt , xerr );
164
165 // test call method four
166 stepper.do_step( oscillator , x , dxdt , t , x , dt , xerr );
167 }
168
169 template< class Stepper >
check_fsal_error_stepper(Stepper & stepper)170 void check_fsal_error_stepper( Stepper &stepper )
171 {
172 typedef Stepper stepper_type;
173 typedef typename stepper_type::state_type state_type;
174 typedef typename stepper_type::value_type value_type;
175 typedef typename stepper_type::deriv_type deriv_type;
176 typedef typename stepper_type::time_type time_type;
177 typedef typename stepper_type::order_type order_type;
178 typedef typename stepper_type::algebra_type algebra_type;
179 typedef typename stepper_type::operations_type operations_type;
180
181 const time_type t( 0.0 * si::second );
182 time_type dt( 0.1 * si::second );
183 state_type x( 1.0 * si::meter , 0.0 * si::meter_per_second ) , xerr;
184
185 // test call method one
186 stepper.do_step( oscillator , x , t , dt , xerr );
187
188 // test call method two
189 stepper.do_step( oscillator , x , t , x , dt , xerr );
190
191 // test call method three
192 deriv_type dxdt;
193 oscillator( x , dxdt , t );
194 stepper.do_step( oscillator , x , dxdt , t , dt , xerr );
195
196 // test call method four
197 stepper.do_step( oscillator , x , dxdt , t , x , dxdt , dt , xerr );
198 }
199
200 template< class Stepper >
check_controlled_stepper(Stepper & stepper)201 void check_controlled_stepper( Stepper &stepper )
202 {
203 typedef Stepper stepper_type;
204 typedef typename stepper_type::state_type state_type;
205 typedef typename stepper_type::value_type value_type;
206 typedef typename stepper_type::deriv_type deriv_type;
207 typedef typename stepper_type::time_type time_type;
208
209 time_type t( 0.0 * si::second );
210 time_type dt( 0.1 * si::second );
211 state_type x( 1.0 * si::meter , 0.0 * si::meter_per_second );
212
213 // test call method one
214 stepper.try_step( oscillator , x , t , dt );
215 }
216
217
218 template< class Stepper >
check_dense_output_stepper(Stepper & stepper)219 void check_dense_output_stepper( Stepper &stepper )
220 {
221 typedef Stepper stepper_type;
222 typedef typename stepper_type::state_type state_type;
223 typedef typename stepper_type::value_type value_type;
224 typedef typename stepper_type::deriv_type deriv_type;
225 typedef typename stepper_type::time_type time_type;
226 // typedef typename stepper_type::order_type order_type;
227
228 time_type t( 0.0 * si::second );
229 time_type dt( 0.1 * si::second );
230 state_type x( 1.0 * si::meter , 0.0 * si::meter_per_second ) , x2;
231
232 stepper.initialize( x , t , dt );
233 stepper.do_step( oscillator );
234 stepper.calc_state( dt / 2.0 , x2 );
235 }
236
237
238
239
240
241
242 class stepper_types : public mpl::vector
243 <
244 euler< state_type , value_type , deriv_type , time_type >,
245 runge_kutta4< state_type , value_type , deriv_type , time_type > ,
246 runge_kutta4_classic< state_type , value_type , deriv_type , time_type > ,
247 runge_kutta_cash_karp54< state_type , value_type , deriv_type , time_type >,
248 runge_kutta_cash_karp54_classic< state_type , value_type , deriv_type , time_type >
249 // don't run rk78 test - gcc requires > 5GB RAM to compile this
250 //, runge_kutta_fehlberg78< state_type , value_type , deriv_type , time_type >
251 > { };
252
253 class fsal_stepper_types : public mpl::vector
254 <
255 runge_kutta_dopri5< state_type , value_type , deriv_type , time_type >
256 > { };
257
258 class error_stepper_types : public mpl::vector
259 <
260 runge_kutta_cash_karp54_classic< state_type , value_type , deriv_type , time_type >
261 //, runge_kutta_fehlberg78< state_type , value_type , deriv_type , time_type >
262 > { };
263
264 class fsal_error_stepper_types : public mpl::vector
265 <
266 runge_kutta_dopri5< state_type , value_type , deriv_type , time_type >
267 > { };
268
269 class controlled_stepper_types : public mpl::vector
270 <
271 controlled_runge_kutta< runge_kutta_cash_karp54_classic< state_type , value_type , deriv_type , time_type > > ,
272 controlled_runge_kutta< runge_kutta_dopri5< state_type , value_type , deriv_type , time_type > >
273 , bulirsch_stoer< state_type , value_type , deriv_type , time_type >
274 // rk78 with units needs up to 3GB memory to compile - disable testing...
275 //, controlled_runge_kutta< runge_kutta_fehlberg78< state_type , value_type , deriv_type , time_type > >
276 > { };
277
278 class dense_output_stepper_types : public mpl::vector
279 <
280 dense_output_runge_kutta< euler< state_type , value_type , deriv_type , time_type > > ,
281 dense_output_runge_kutta<
282 controlled_runge_kutta< runge_kutta_dopri5< state_type , value_type , deriv_type , time_type > > >
283 //, bulirsch_stoer_dense_out< state_type , value_type , deriv_type , time_type >
284 > { };
285
286
287
288
289 BOOST_AUTO_TEST_SUITE( stepper_with_units )
290
BOOST_AUTO_TEST_CASE_TEMPLATE(stepper_test,Stepper,stepper_types)291 BOOST_AUTO_TEST_CASE_TEMPLATE( stepper_test , Stepper , stepper_types )
292 {
293 Stepper stepper;
294 check_stepper( stepper );
295 }
296
BOOST_AUTO_TEST_CASE_TEMPLATE(fsl_stepper_test,Stepper,fsal_stepper_types)297 BOOST_AUTO_TEST_CASE_TEMPLATE( fsl_stepper_test , Stepper , fsal_stepper_types )
298 {
299 Stepper stepper;
300 check_fsal_stepper( stepper );
301 }
302
BOOST_AUTO_TEST_CASE_TEMPLATE(error_stepper_test,Stepper,error_stepper_types)303 BOOST_AUTO_TEST_CASE_TEMPLATE( error_stepper_test , Stepper , error_stepper_types )
304 {
305 Stepper stepper;
306 check_error_stepper( stepper );
307 }
308
BOOST_AUTO_TEST_CASE_TEMPLATE(fsal_error_stepper_test,Stepper,fsal_error_stepper_types)309 BOOST_AUTO_TEST_CASE_TEMPLATE( fsal_error_stepper_test , Stepper , fsal_error_stepper_types )
310 {
311 Stepper stepper;
312 check_fsal_error_stepper( stepper );
313 }
314
BOOST_AUTO_TEST_CASE_TEMPLATE(controlled_stepper_test,Stepper,controlled_stepper_types)315 BOOST_AUTO_TEST_CASE_TEMPLATE( controlled_stepper_test , Stepper , controlled_stepper_types )
316 {
317 Stepper stepper;
318 check_controlled_stepper( stepper );
319 }
320
BOOST_AUTO_TEST_CASE_TEMPLATE(dense_ouput_test,Stepper,dense_output_stepper_types)321 BOOST_AUTO_TEST_CASE_TEMPLATE( dense_ouput_test , Stepper , dense_output_stepper_types )
322 {
323 Stepper stepper;
324 check_dense_output_stepper( stepper );
325 }
326
327
328 BOOST_AUTO_TEST_SUITE_END()
329