diff options
Diffstat (limited to 'src/boost/libs/numeric/odeint/examples/stepper_details.cpp')
| -rw-r--r-- | src/boost/libs/numeric/odeint/examples/stepper_details.cpp | 200 |
1 files changed, 200 insertions, 0 deletions
diff --git a/src/boost/libs/numeric/odeint/examples/stepper_details.cpp b/src/boost/libs/numeric/odeint/examples/stepper_details.cpp new file mode 100644 index 00000000..d4ae8bd1 --- /dev/null +++ b/src/boost/libs/numeric/odeint/examples/stepper_details.cpp @@ -0,0 +1,200 @@ +/* + * stepper_details.cpp + * + * This example demonstrates some details about the steppers in odeint. + * + * + * Copyright 2011-2012 Karsten Ahnert + * Copyright 2012 Mario Mulansky + * Copyright 2013 Pascal Germroth + * + * 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) + */ + +#include <iostream> +#include <boost/array.hpp> +#include <boost/bind.hpp> +#include <boost/numeric/odeint.hpp> + +using namespace std; +using namespace boost::numeric::odeint; + +const size_t N = 3; + +typedef boost::array< double , N > state_type; + +//[ system_function_structure +void sys( const state_type & /*x*/ , state_type & /*dxdt*/ , const double /*t*/ ) +{ + // ... +} +//] + +void sys1( const state_type &/*x*/ , state_type &/*dxdt*/ , const double /*t*/ ) +{ +} + +void sys2( const state_type &/*x*/ , state_type &/*dxdt*/ , const double /*t*/ ) +{ +} + + +//[ symplectic_stepper_detail_system_function +typedef boost::array< double , 1 > vector_type; + + +struct harm_osc_f1 +{ + void operator()( const vector_type &p , vector_type &dqdt ) + { + dqdt[0] = p[0]; + } +}; + +struct harm_osc_f2 +{ + void operator()( const vector_type &q , vector_type &dpdt ) + { + dpdt[0] = -q[0]; + } +}; +//] + +//[ symplectic_stepper_detail_system_class +struct harm_osc +{ + void f1( const vector_type &p , vector_type &dqdt ) const + { + dqdt[0] = p[0]; + } + + void f2( const vector_type &q , vector_type &dpdt ) const + { + dpdt[0] = -q[0]; + } +}; +//] + +int main( int argc , char **argv ) +{ + using namespace std; + + // Explicit stepper example + { + double t( 0.0 ) , dt( 0.1 ); + state_type in , out , dxdtin , inout; + //[ explicit_stepper_detail_example + runge_kutta4< state_type > rk; + rk.do_step( sys1 , inout , t , dt ); // In-place transformation of inout + rk.do_step( sys2 , inout , t , dt ); // call with different system: Ok + rk.do_step( sys1 , in , t , out , dt ); // Out-of-place transformation + rk.do_step( sys1 , inout , dxdtin , t , dt ); // In-place tranformation of inout + rk.do_step( sys1 , in , dxdtin , t , out , dt ); // Out-of-place transformation + //] + } + + + + // FSAL stepper example + { + double t( 0.0 ) , dt( 0.1 ); + state_type in , in2 , in3 , out , dxdtin , dxdtout , inout , dxdtinout; + //[ fsal_stepper_detail_example + runge_kutta_dopri5< state_type > rk; + rk.do_step( sys1 , in , t , out , dt ); + rk.do_step( sys2 , in , t , out , dt ); // DONT do this, sys1 is assumed + + rk.do_step( sys2 , in2 , t , out , dt ); + rk.do_step( sys2 , in3 , t , out , dt ); // DONT do this, in2 is assumed + + rk.do_step( sys1 , inout , dxdtinout , t , dt ); + rk.do_step( sys2 , inout , dxdtinout , t , dt ); // Ok, internal derivative is not used, dxdtinout is updated + + rk.do_step( sys1 , in , dxdtin , t , out , dxdtout , dt ); + rk.do_step( sys2 , in , dxdtin , t , out , dxdtout , dt ); // Ok, internal derivative is not used + //] + } + + + // Symplectic harmonic oscillator example + { + double t( 0.0 ) , dt( 0.1 ); + //[ symplectic_stepper_detail_example + pair< vector_type , vector_type > x; + x.first[0] = 1.0; x.second[0] = 0.0; + symplectic_rkn_sb3a_mclachlan< vector_type > rkn; + rkn.do_step( make_pair( harm_osc_f1() , harm_osc_f2() ) , x , t , dt ); + //] + + //[ symplectic_stepper_detail_system_class_example + harm_osc h; + rkn.do_step( make_pair( boost::bind( &harm_osc::f1 , h , _1 , _2 ) , boost::bind( &harm_osc::f2 , h , _1 , _2 ) ) , + x , t , dt ); + //] + } + + // Simplified harmonic oscillator example + { + double t = 0.0, dt = 0.1; + //[ simplified_symplectic_stepper_example + pair< vector_type , vector_type > x; + x.first[0] = 1.0; x.second[0] = 0.0; + symplectic_rkn_sb3a_mclachlan< vector_type > rkn; + rkn.do_step( harm_osc_f1() , x , t , dt ); + //] + + vector_type q = {{ 1.0 }} , p = {{ 0.0 }}; + //[ symplectic_stepper_detail_ref_usage + rkn.do_step( harm_osc_f1() , make_pair( boost::ref( q ) , boost::ref( p ) ) , t , dt ); + rkn.do_step( harm_osc_f1() , q , p , t , dt ); + rkn.do_step( make_pair( harm_osc_f1() , harm_osc_f2() ) , q , p , t , dt ); + //] + } + + // adams_bashforth_moulton stepper example + { + double t = 0.0 , dt = 0.1; + state_type inout; + //[ multistep_detail_example + adams_bashforth_moulton< 5 , state_type > abm; + abm.initialize( sys , inout , t , dt ); + abm.do_step( sys , inout , t , dt ); + //] + + //[ multistep_detail_own_stepper_initialization + abm.initialize( runge_kutta_fehlberg78< state_type >() , sys , inout , t , dt ); + //] + } + + + + // dense output stepper examples + { + double t = 0.0 , dt = 0.1; + state_type in; + //[ dense_output_detail_example + dense_output_runge_kutta< controlled_runge_kutta< runge_kutta_dopri5< state_type > > > dense; + dense.initialize( in , t , dt ); + pair< double , double > times = dense.do_step( sys ); + (void)times; + //] + + state_type inout; + double t_start = 0.0 , t_end = 1.0; + //[ dense_output_detail_generation1 + typedef boost::numeric::odeint::result_of::make_dense_output< + runge_kutta_dopri5< state_type > >::type dense_stepper_type; + dense_stepper_type dense2 = make_dense_output( 1.0e-6 , 1.0e-6 , runge_kutta_dopri5< state_type >() ); + (void)dense2; + //] + + //[ dense_output_detail_generation2 + integrate_const( make_dense_output( 1.0e-6 , 1.0e-6 , runge_kutta_dopri5< state_type >() ) , sys , inout , t_start , t_end , dt ); + //] + } + + + return 0; +} |