ceph/src/boost/libs/numeric/odeint/test/integrate_implicit.cpp

   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 {
  62     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 {
  71     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
  86     push_back_time( std::vector< double > &times )
  87     :  m_times( times ) { }
  88
  89     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 {
  98     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 {
 123     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 {
 145     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 {
 169     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
 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
 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
 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
 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()