1 /*
2 [auto_generated]
3 libs/numeric/odeint/test/integrate_implicit.cpp
4
5 [begin_description]
6 This file tests the integrate function and its variants with the implicit steppers.
7 [end_description]
8
9 Copyright 2011 Mario Mulansky
10 Copyright 2011-2013 Karsten Ahnert
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_integrate_functions_implicit
18
19 #include <vector>
20 #include <cmath>
21 #include <iostream>
22
23 #include <boost/numeric/odeint/config.hpp>
24
25 #include <boost/array.hpp>
26 #include <boost/ref.hpp>
27 #include <boost/iterator/counting_iterator.hpp>
28 #include <boost/numeric/ublas/vector.hpp>
29 #include <boost/numeric/ublas/matrix.hpp>
30
31 #include <boost/test/unit_test.hpp>
32
33 #include <boost/mpl/vector.hpp>
34
35 #ifndef ODEINT_INTEGRATE_ITERATOR
36 #include <boost/numeric/odeint/integrate/integrate_const.hpp>
37 #include <boost/numeric/odeint/integrate/integrate_adaptive.hpp>
38 #include <boost/numeric/odeint/integrate/integrate_times.hpp>
39 #include <boost/numeric/odeint/integrate/integrate_n_steps.hpp>
40 #else
41 #include <boost/numeric/odeint/iterator/integrate/integrate_const.hpp>
42 #include <boost/numeric/odeint/iterator/integrate/integrate_adaptive.hpp>
43 #include <boost/numeric/odeint/iterator/integrate/integrate_times.hpp>
44 #include <boost/numeric/odeint/iterator/integrate/integrate_n_steps.hpp>
45 #endif
46 #include <boost/numeric/odeint/stepper/rosenbrock4.hpp>
47 #include <boost/numeric/odeint/stepper/rosenbrock4_controller.hpp>
48 #include <boost/numeric/odeint/stepper/rosenbrock4_dense_output.hpp>
49
50
51 using namespace boost::unit_test;
52 using namespace boost::numeric::odeint;
53 namespace mpl = boost::mpl;
54 namespace ublas = boost::numeric::ublas;
55
56 typedef double value_type;
57 typedef ublas::vector< value_type > state_type;
58 typedef boost::numeric::ublas::matrix< value_type > matrix_type;
59
60 struct sys
61 {
operator ()sys62 void operator()( const state_type &x , state_type &dxdt , const value_type &t ) const
63 {
64 dxdt( 0 ) = x( 0 ) + 2 * x( 1 );
65 dxdt( 1 ) = x( 1 );
66 }
67 };
68
69 struct jacobi
70 {
operator ()jacobi71 void operator()( const state_type &x , matrix_type &jacobi , const value_type &t , state_type &dfdt ) const
72 {
73 jacobi( 0 , 0 ) = 1;
74 jacobi( 0 , 1 ) = 2;
75 jacobi( 1 , 0 ) = 0;
76 jacobi( 1 , 1 ) = 1;
77 dfdt( 0 ) = 0.0;
78 dfdt( 1 ) = 0.0;
79 }
80 };
81
82 struct push_back_time
83 {
84 std::vector< double >& m_times;
85
push_back_timepush_back_time86 push_back_time( std::vector< double > × )
87 : m_times( times ) { }
88
operator ()push_back_time89 void operator()( const state_type &x , double t )
90 {
91 m_times.push_back( t );
92 }
93 };
94
95 template< class Stepper >
96 struct perform_integrate_const_test
97 {
operator ()perform_integrate_const_test98 void operator()( void )
99 {
100 state_type x( 2 , 1.0 );
101 const value_type dt = 0.03;
102 const value_type t_end = 10.0;
103
104 std::vector< value_type > times;
105
106 integrate_const( Stepper() , std::make_pair( sys() , jacobi() ) , x , 0.0 , t_end ,
107 dt , push_back_time( times ) );
108
109 BOOST_CHECK_EQUAL( static_cast<int>(times.size()) , static_cast<int>(std::floor(t_end/dt))+1 );
110
111 for( size_t i=0 ; i<times.size() ; ++i )
112 {
113 //std::cout << i << std::endl;
114 // check if observer was called at times 0,1,2,...
115 BOOST_CHECK_SMALL( times[i] - static_cast< value_type >(i)*dt , (i+1) * 2.0e-16 );
116 }
117 }
118 };
119
120 template< class Stepper >
121 struct perform_integrate_adaptive_test
122 {
operator ()perform_integrate_adaptive_test123 void operator()( void )
124 {
125 state_type x( 2 , 1.0 );
126 const value_type dt = 0.03;
127 const value_type t_end = 10.0;
128
129 std::vector< value_type > times;
130
131 size_t steps = integrate_adaptive( Stepper() , std::make_pair( sys() , jacobi() ) , x , 0.0 , t_end ,
132 dt , push_back_time( times ) );
133
134 BOOST_CHECK_EQUAL( times.size() , steps+1 );
135
136 BOOST_CHECK_SMALL( times[0] - 0.0 , 2.0e-16 );
137 BOOST_CHECK_SMALL( times[times.size()-1] - t_end , times.size() * 3.0e-16 );
138 }
139 };
140
141
142 template< class Stepper >
143 struct perform_integrate_times_test
144 {
operator ()perform_integrate_times_test145 void operator()( void )
146 {
147 state_type x( 2 , 1.0 );
148
149 const value_type dt = 0.03;
150
151 std::vector< double > times;
152
153 // simple stepper
154 integrate_times( Stepper() , std::make_pair( sys() , jacobi() ) , x ,
155 boost::counting_iterator<int>(0) , boost::counting_iterator<int>(10) ,
156 dt , push_back_time( times ) );
157
158 BOOST_CHECK_EQUAL( static_cast<int>(times.size()) , 10 );
159
160 for( size_t i=0 ; i<times.size() ; ++i )
161 // check if observer was called at times 0,1,2,...
162 BOOST_CHECK_EQUAL( times[i] , static_cast<double>(i) );
163 }
164 };
165
166 template< class Stepper >
167 struct perform_integrate_n_steps_test
168 {
operator ()perform_integrate_n_steps_test169 void operator()( void )
170 {
171 state_type x( 2 , 1.0 );
172
173 const value_type dt = 0.03;
174 const int n = 200;
175
176 std::vector< double > times;
177
178 // simple stepper
179 value_type end_time = integrate_n_steps( Stepper() , std::make_pair( sys() , jacobi() ) , x , 0.0 , dt , n , push_back_time( times ) );
180
181 BOOST_CHECK_SMALL( end_time - n*dt , 3.0e-16 );
182 BOOST_CHECK_EQUAL( static_cast<int>(times.size()) , n+1 );
183
184 for( size_t i=0 ; i<times.size() ; ++i )
185 // check if observer was called at times 0,1,2,...
186 BOOST_CHECK_SMALL( times[i] - static_cast< value_type >(i)*dt , (i+1) * 2.0e-16 );
187 }
188 };
189
190
191
192 class stepper_methods : public mpl::vector<
193 rosenbrock4< value_type > ,
194 rosenbrock4_controller< rosenbrock4< value_type > > ,
195 rosenbrock4_dense_output< rosenbrock4_controller< rosenbrock4< value_type > > >
196 > { };
197
198
199
200 BOOST_AUTO_TEST_SUITE( integrate_test )
201
BOOST_AUTO_TEST_CASE_TEMPLATE(integrate_const_test_case,Stepper,stepper_methods)202 BOOST_AUTO_TEST_CASE_TEMPLATE( integrate_const_test_case , Stepper, stepper_methods )
203 {
204 perform_integrate_const_test< Stepper > tester;
205 tester();
206 }
207
208
BOOST_AUTO_TEST_CASE_TEMPLATE(integrate_adaptive_test_case,Stepper,stepper_methods)209 BOOST_AUTO_TEST_CASE_TEMPLATE( integrate_adaptive_test_case , Stepper, stepper_methods )
210 {
211 perform_integrate_adaptive_test< Stepper > tester;
212 tester();
213 }
214
215
BOOST_AUTO_TEST_CASE_TEMPLATE(integrate_times_test_case,Stepper,stepper_methods)216 BOOST_AUTO_TEST_CASE_TEMPLATE( integrate_times_test_case , Stepper, stepper_methods )
217 {
218 perform_integrate_times_test< Stepper > tester;
219 tester();
220 }
221
222
223 class simple_stepper_methods : public mpl::vector<
224 rosenbrock4< value_type >
225 > { };
226
BOOST_AUTO_TEST_CASE_TEMPLATE(integrate_n_steps_test_case,Stepper,simple_stepper_methods)227 BOOST_AUTO_TEST_CASE_TEMPLATE( integrate_n_steps_test_case , Stepper, simple_stepper_methods )
228 {
229 perform_integrate_n_steps_test< Stepper > tester;
230 tester();
231 }
232
233 BOOST_AUTO_TEST_SUITE_END()
234