156 lines
4.2 KiB
C++
156 lines
4.2 KiB
C++
/*
|
|
[auto_generated]
|
|
libs/numeric/odeint/test/bulirsch_stoer.cpp
|
|
|
|
[begin_description]
|
|
This file tests the Bulirsch-Stoer stepper.
|
|
[end_description]
|
|
|
|
Copyright 2011 Mario Mulansky
|
|
Copyright 2012 Karsten Ahnert
|
|
|
|
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)
|
|
*/
|
|
|
|
|
|
// disable checked iterator warning for msvc
|
|
#include <boost/config.hpp>
|
|
#ifdef BOOST_MSVC
|
|
#pragma warning(disable:4996)
|
|
#endif
|
|
|
|
#define BOOST_TEST_MODULE odeint_bulirsch_stoer
|
|
|
|
#include <utility>
|
|
#include <iostream>
|
|
|
|
#include <boost/array.hpp>
|
|
|
|
#include <boost/test/unit_test.hpp>
|
|
|
|
#include <boost/numeric/odeint/stepper/bulirsch_stoer.hpp>
|
|
#include <boost/numeric/odeint/stepper/bulirsch_stoer_dense_out.hpp>
|
|
|
|
#include <boost/numeric/odeint/integrate/integrate_adaptive.hpp>
|
|
|
|
using namespace boost::unit_test;
|
|
using namespace boost::numeric::odeint;
|
|
|
|
typedef double value_type;
|
|
typedef boost::array< value_type , 3 > state_type;
|
|
|
|
const double sigma = 10.0;
|
|
const double R = 28.0;
|
|
const double b = 8.0 / 3.0;
|
|
|
|
struct lorenz
|
|
{
|
|
template< class State , class Deriv >
|
|
void operator()( const State &x , Deriv &dxdt , double t ) const
|
|
{
|
|
dxdt[0] = sigma * ( x[1] - x[0] );
|
|
dxdt[1] = R * x[0] - x[1] - x[0] * x[2];
|
|
dxdt[2] = -b * x[2] + x[0] * x[1];
|
|
}
|
|
};
|
|
|
|
struct const_system
|
|
{
|
|
template< class State , class Deriv >
|
|
void operator()( const State &x , Deriv &dxdt , double t ) const
|
|
{
|
|
dxdt[0] = 1.0;
|
|
dxdt[1] = 1.0;
|
|
dxdt[2] = 1.0;
|
|
}
|
|
};
|
|
|
|
struct sin_system
|
|
{
|
|
template< class State , class Deriv >
|
|
void operator()( const State &x , Deriv &dxdt , double t ) const
|
|
{
|
|
dxdt[0] = sin( x[0] );
|
|
dxdt[1] = cos( x[1] );
|
|
dxdt[2] = sin( x[2] ) + cos( x[2] );
|
|
}
|
|
};
|
|
|
|
BOOST_AUTO_TEST_SUITE( bulirsch_stoer_test )
|
|
|
|
BOOST_AUTO_TEST_CASE( test_bulirsch_stoer )
|
|
{
|
|
typedef bulirsch_stoer< state_type > stepper_type;
|
|
stepper_type stepper( 1E-9 , 1E-9 , 1.0 , 0.0 );
|
|
|
|
state_type x;
|
|
x[0] = 10.0 ; x[1] = 10.0 ; x[2] = 5.0;
|
|
|
|
double dt = 0.1;
|
|
|
|
//stepper.try_step( lorenz() , x , t , dt );
|
|
|
|
std::cout << "starting integration..." << std::endl;
|
|
|
|
size_t steps = integrate_adaptive( stepper , lorenz() , x , 0.0 , 10.0 , dt );
|
|
|
|
std::cout << "required steps: " << steps << std::endl;
|
|
|
|
bulirsch_stoer_dense_out< state_type > bs_do( 1E-9 , 1E-9 , 1.0 , 0.0 );
|
|
x[0] = 10.0 ; x[1] = 10.0 ; x[2] = 5.0;
|
|
double t = 0.0;
|
|
dt = 1E-1;
|
|
bs_do.initialize( x , t , dt );
|
|
bs_do.do_step( sin_system() );
|
|
std::cout << "one step successful, new time: " << bs_do.current_time() << " (" << t << ")" << std::endl;
|
|
|
|
x = bs_do.current_state();
|
|
std::cout << "x( " << bs_do.current_time() << " ) = [ " << x[0] << " , " << x[1] << " , " << x[2] << " ]" << std::endl;
|
|
|
|
bs_do.calc_state( bs_do.current_time()/3 , x );
|
|
std::cout << "x( " << bs_do.current_time()/3 << " ) = [ " << x[0] << " , " << x[1] << " , " << x[2] << " ]" << std::endl;
|
|
|
|
std::cout << std::endl << "=======================================================================" << std::endl << std::endl;
|
|
|
|
x[0] = 10.0 ; x[1] = 10.0 ; x[2] = 5.0;
|
|
t = 0.0; dt /= 3;
|
|
bs_do.initialize( x , t , dt );
|
|
bs_do.do_step( sin_system() );
|
|
x = bs_do.current_state();
|
|
std::cout << "x( " << bs_do.current_time() << " ) = [ " << x[0] << " , " << x[1] << " , " << x[2] << " ]" << std::endl;
|
|
|
|
t = dt;
|
|
bs_do.initialize( x , t , dt );
|
|
bs_do.do_step( sin_system() );
|
|
x = bs_do.current_state();
|
|
|
|
t = 2*dt;
|
|
bs_do.initialize( x , t , dt );
|
|
bs_do.do_step( sin_system() );
|
|
x = bs_do.current_state();
|
|
|
|
std::cout << "x( " << bs_do.current_time() << " ) = [ " << x[0] << " , " << x[1] << " , " << x[2] << " ]" << std::endl << std::endl << std::endl;
|
|
}
|
|
|
|
BOOST_AUTO_TEST_CASE( test_bulirsch_stoer_adjust_size )
|
|
{
|
|
typedef bulirsch_stoer< state_type > stepper_type;
|
|
stepper_type stepper( 1E-9 , 1E-9 , 1.0 , 0.0 );
|
|
|
|
state_type x; x[0] = 10.0 ; x[1] = 10.0 ; x[2] = 5.0;
|
|
|
|
stepper.adjust_size( x );
|
|
|
|
|
|
double dt = 0.1;
|
|
|
|
size_t steps = integrate_adaptive( stepper , lorenz() , x , 0.0 , 10.0 , dt );
|
|
|
|
std::cout << "required steps: " << steps << std::endl;
|
|
}
|
|
|
|
|
|
BOOST_AUTO_TEST_SUITE_END()
|