| /* Boost check_gmp.cpp test file |
| |
| Copyright 2010-2012 Mario Mulansky |
| Copyright 2011-2012 Karsten Ahnert |
| |
| This file tests the odeint library with the gmp arbitrary precision types |
| |
| Distributed under the Boost Software License, Version 1.0. |
| (See accompanying file LICENSE_1_0.txt or |
| copy at http://www.boost.org/LICENSE_1_0.txt) |
| */ |
| |
| #define BOOST_TEST_MODULE odeint_gmp |
| |
| #include <gmpxx.h> |
| |
| #include <boost/test/unit_test.hpp> |
| #include <boost/array.hpp> |
| |
| #include <boost/mpl/vector.hpp> |
| |
| #include <boost/numeric/odeint.hpp> |
| //#include <boost/numeric/odeint/algebra/vector_space_algebra.hpp> |
| |
| using namespace boost::unit_test; |
| using namespace boost::numeric::odeint; |
| |
| namespace mpl = boost::mpl; |
| |
| const int precision = 1024; |
| |
| typedef mpf_class value_type; |
| typedef mpf_class state_type; |
| |
| //provide min, max and pow functions for mpf types - required for controlled steppers |
| value_type min( const value_type a , const value_type b ) |
| { |
| if( a<b ) return a; |
| else return b; |
| } |
| value_type max( const value_type a , const value_type b ) |
| { |
| if( a>b ) return a; |
| else return b; |
| } |
| value_type pow( const value_type a , const value_type b ) |
| { |
| // do the calculation in double precision... |
| return value_type( std::pow( a.get_d() , b.get_d() ) ); |
| } |
| |
| |
| //provide vector_space reduce: |
| |
| namespace boost { namespace numeric { namespace odeint { |
| |
| template<> |
| struct vector_space_reduce< state_type > |
| { |
| template< class Op > |
| state_type operator()( state_type x , Op op , state_type init ) const |
| { |
| init = op( init , x ); |
| return init; |
| } |
| }; |
| |
| } } } |
| |
| |
| void constant_system( const state_type &x , state_type &dxdt , value_type t ) |
| { |
| dxdt = value_type( 1.0 , precision ); |
| } |
| |
| /* check runge kutta stepers */ |
| typedef mpl::vector< |
| euler< state_type , value_type , state_type , value_type , vector_space_algebra > , |
| modified_midpoint< state_type , value_type , state_type , value_type , vector_space_algebra > , |
| runge_kutta4< state_type , value_type , state_type , value_type , vector_space_algebra > , |
| runge_kutta4_classic< state_type , value_type , state_type , value_type , vector_space_algebra > , |
| runge_kutta_cash_karp54_classic< state_type , value_type , state_type , value_type , vector_space_algebra > , |
| runge_kutta_cash_karp54< state_type , value_type , state_type , value_type , vector_space_algebra > , |
| runge_kutta_dopri5< state_type , value_type , state_type , value_type , vector_space_algebra > , |
| runge_kutta_fehlberg78< state_type , value_type , state_type , value_type , vector_space_algebra > |
| > stepper_types; |
| |
| |
| template< class Stepper > |
| struct perform_integrate_const_test { |
| |
| void operator()( void ) |
| { |
| /* We have to specify the desired precision in advance! */ |
| mpf_set_default_prec( precision ); |
| |
| mpf_t eps_ , unity; |
| mpf_init( eps_ ); mpf_init( unity ); |
| mpf_set_d( unity , 1.0 ); |
| mpf_div_2exp( eps_ , unity , precision-1 ); // 2^(-precision+1) : smallest number that can be represented with used precision |
| value_type eps( eps_ ); |
| |
| Stepper stepper; |
| state_type x; |
| x = 0.0; |
| value_type t0( 0.0 ); |
| value_type tend( 1.0 ); |
| value_type dt(0.1); |
| |
| integrate_const( stepper , constant_system , x , t0 , tend , dt ); |
| |
| x = 0.0; |
| t0 = 0.0; |
| dt = 0.1; |
| size_t steps = 10; |
| |
| integrate_n_steps( stepper , constant_system , x , t0 , dt , steps ); |
| |
| BOOST_CHECK_MESSAGE( abs( x - 10*dt ) < eps , x ); |
| } |
| }; |
| |
| |
| BOOST_AUTO_TEST_CASE_TEMPLATE( integrate_const_test , Stepper , stepper_types ) |
| { |
| perform_integrate_const_test< Stepper > tester; |
| tester(); |
| } |
| |
| |
| typedef mpl::vector< |
| controlled_runge_kutta< runge_kutta_cash_karp54_classic< state_type , value_type , state_type , value_type , vector_space_algebra > > , |
| controlled_runge_kutta< runge_kutta_dopri5< state_type , value_type , state_type , value_type , vector_space_algebra > > , |
| controlled_runge_kutta< runge_kutta_fehlberg78< state_type , value_type , state_type , value_type , vector_space_algebra > > , |
| bulirsch_stoer< state_type , value_type , state_type , value_type , vector_space_algebra > |
| > controlled_stepper_types; |
| |
| |
| template< class Stepper > |
| struct perform_integrate_adaptive_test { |
| |
| void operator()( void ) |
| { |
| mpf_set_default_prec( precision ); |
| |
| mpf_t eps_ , unity; |
| mpf_init( eps_ ); mpf_init( unity ); |
| mpf_set_d( unity , 1.0 ); |
| mpf_div_2exp( eps_ , unity , precision-1 ); // 2^(-precision+1) : smallest number that can be represented with used precision |
| value_type eps( eps_ ); |
| |
| Stepper stepper; |
| state_type x; |
| x = 0.0; |
| value_type t0( 0.0 ); |
| value_type tend( 1.0 ); |
| value_type dt(0.1); |
| |
| integrate_adaptive( stepper , constant_system , x , t0 , tend , dt ); |
| |
| BOOST_CHECK_MESSAGE( abs( x - tend ) < eps , x - 0.1 ); |
| } |
| }; |
| |
| BOOST_AUTO_TEST_CASE_TEMPLATE( integrate_adaptive__test , Stepper , controlled_stepper_types ) |
| { |
| perform_integrate_adaptive_test< Stepper > tester; |
| tester(); |
| } |