|
|
Boost-Commit : |
Subject: [Boost-commit] svn:boost r69465 - in sandbox/odeint/branches/karsten: . boost/numeric/odeint/integrate/detail boost/numeric/odeint/stepper boost/numeric/odeint/stepper/base libs/numeric/odeint/regression_test libs/numeric/odeint/test
From: karsten.ahnert_at_[hidden]
Date: 2011-03-02 02:33:08
Author: karsten
Date: 2011-03-02 02:33:04 EST (Wed, 02 Mar 2011)
New Revision: 69465
URL: http://svn.boost.org/trac/boost/changeset/69465
Log:
range for dense output and controlled stepper
Text files modified:
sandbox/odeint/branches/karsten/TODO | 15 ++++---
sandbox/odeint/branches/karsten/boost/numeric/odeint/integrate/detail/integrate_adaptive.hpp | 12 +++++-
sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/base/symplectic_rkn_stepper_base.hpp | 27 ++++++++------
sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_controlled_explicit_fsal.hpp | 11 ++++++
sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_explicit.hpp | 9 +++++
sandbox/odeint/branches/karsten/libs/numeric/odeint/regression_test/integrate_functions.cpp | 2
sandbox/odeint/branches/karsten/libs/numeric/odeint/test/stepper_with_ranges.cpp | 72 ++++++++++++++++++++++++++++++++++++++-
7 files changed, 123 insertions(+), 25 deletions(-)
Modified: sandbox/odeint/branches/karsten/TODO
==============================================================================
--- sandbox/odeint/branches/karsten/TODO (original)
+++ sandbox/odeint/branches/karsten/TODO 2011-03-02 02:33:04 EST (Wed, 02 Mar 2011)
@@ -16,13 +16,7 @@
OK * move error_checker into controlled_stepper
* rename controlled_stepper to a more specific name
* ranges:
- * test ranges in symplectic_rkn_stepper
- * dense_output
* integrate functions
- OK * explicit_stepper_and_error_stepper_fsal_base
- OK * controlled_error_stepper
- OK * check comments (spelling and if the comment is true, in some versions dxdt is already const)
- OK * check names of the impl functions
* general:
* check if everywhere static_cast< value_type > is used
* check header guards
@@ -71,4 +65,11 @@
OK * test fusion_algebra
OK * test, if copy construct of stepper_base is called when explicit_euler is used
OK * test units with dense output
-OK * include rosenbrock4 in trunk
\ No newline at end of file
+OK * include rosenbrock4 in trunk
+* ranges:
+ OK * test ranges in symplectic_rkn_stepper
+ OK * dense_output
+ OK * explicit_stepper_and_error_stepper_fsal_base
+ OK * controlled_error_stepper
+ OK * check comments (spelling and if the comment is true, in some versions dxdt is already const)
+ OK * check names of the impl functions
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint/integrate/detail/integrate_adaptive.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/integrate/detail/integrate_adaptive.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/integrate/detail/integrate_adaptive.hpp 2011-03-02 02:33:04 EST (Wed, 02 Mar 2011)
@@ -14,9 +14,17 @@
namespace detail {
template< class Stepper , class System , class State , class Time , class Observer >
-size_t integrate_adaptive( Stepper stepper , System system , State &start_state , const Time &start_time , const Time &end_time , const Time &dt , Observer observer , stepper_tag )
+size_t integrate_adaptive( Stepper stepper , System system , State &start_state , Time start_time , const Time &end_time , Time dt , Observer observer , stepper_tag )
{
- return 0;
+ size_t count = 0;
+ while( start_time < end_time )
+ {
+ stepper.do_step( system , start_state , start_time , dt );
+ observer( start_state , start_time );
+ start_time += dt;
+ ++count;
+ }
+ return count;
}
template< class Stepper , class System , class State , class Time , class Observer >
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/base/symplectic_rkn_stepper_base.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/base/symplectic_rkn_stepper_base.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/base/symplectic_rkn_stepper_base.hpp 2011-03-02 02:33:04 EST (Wed, 02 Mar 2011)
@@ -111,7 +111,7 @@
/*
* Version 1 : do_step( system , x , t , dt )
*
- * The overloads are needed in order to solve the forwarding problem
+ * This version does not solve the forwarding problem, boost.range can not be used.
*/
template< class System , class StateInOut >
void do_step( System system , const StateInOut &state , const time_type &t , const time_type &dt )
@@ -128,7 +128,11 @@
}
- // for convenience
+ /*
+ * Version 2 : do_step( system , q , p , t , dt );
+ *
+ * The two overloads are needed in order to solve the forwarding problem.
+ */
template< class System , class CoorInOut , class MomentumInOut >
void do_step( System system , CoorInOut &q , MomentumInOut &p , const time_type &t , const time_type &dt )
{
@@ -187,9 +191,9 @@
typedef typename boost::unwrap_reference< StateIn >::type state_in_type;
typedef typename boost::unwrap_reference< typename state_in_type::first_type >::type coor_in_type;
typedef typename boost::unwrap_reference< typename state_in_type::second_type >::type momentum_in_type;
- state_in_type &state_in = in;
- coor_in_type &coor_in = state_in.first;
- momentum_in_type &momentum_in = state_in.first;
+ const state_in_type &state_in = in;
+ const coor_in_type &coor_in = state_in.first;
+ const momentum_in_type &momentum_in = state_in.second;
typedef typename boost::unwrap_reference< StateOut >::type state_out_type;
typedef typename boost::unwrap_reference< typename state_out_type::first_type >::type coor_out_type;
@@ -209,19 +213,19 @@
if( l == 0 )
{
coor_func( momentum_in , m_dqdt );
- algebra_type::for_each2( coor_out , coor_in , m_dqdt ,
+ algebra_type::for_each3( coor_out , coor_in , m_dqdt ,
typename operations_type::template scale_sum2< value_type , time_type >( 1.0 , m_coef_a[l] * dt ) );
- momentum_func( coor_out , m_dqdt );
- algebra_type::for_each2( momentum_out , momentum_in , m_dqdt ,
+ momentum_func( coor_out , m_dpdt );
+ algebra_type::for_each3( momentum_out , momentum_in , m_dpdt ,
typename operations_type::template scale_sum2< value_type , time_type >( 1.0 , m_coef_b[l] * dt ) );
}
else
{
coor_func( momentum_out , m_dqdt );
- algebra_type::for_each2( coor_out , coor_out , m_dqdt ,
+ algebra_type::for_each3( coor_out , coor_out , m_dqdt ,
typename operations_type::template scale_sum2< value_type , time_type >( 1.0 , m_coef_a[l] * dt ) );
- momentum_func( coor_out , m_dqdt );
- algebra_type::for_each2( momentum_out , momentum_out , m_dqdt ,
+ momentum_func( coor_out , m_dpdt );
+ algebra_type::for_each3( momentum_out , momentum_out , m_dpdt ,
typename operations_type::template scale_sum2< value_type , time_type >( 1.0 , m_coef_b[l] * dt ) );
}
}
@@ -277,7 +281,6 @@
}
}
- // ToDo : check if a and b are appropriate names for the coefficients
const coef_type m_coef_a;
const coef_type m_coef_b;
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_controlled_explicit_fsal.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_controlled_explicit_fsal.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_controlled_explicit_fsal.hpp 2011-03-02 02:33:04 EST (Wed, 02 Mar 2011)
@@ -167,12 +167,23 @@
return std::make_pair( m_t_old , m_t );
}
+
+ /*
+ * The two overloads are needed in order to solve the forwarding problem.
+ */
template< class StateOut >
void calc_state( const time_type &t , StateOut &x )
{
m_stepper.stepper().calc_state( t , x , *m_old_state , *m_old_deriv , m_t_old , *m_current_state , *m_current_deriv , m_t );
}
+ template< class StateOut >
+ void calc_state( const time_type &t , const StateOut &x )
+ {
+ m_stepper.stepper().calc_state( t , x , *m_old_state , *m_old_deriv , m_t_old , *m_current_state , *m_current_deriv , m_t );
+ }
+
+
template< class StateType >
void adjust_size( const StateType &x )
{
Modified: sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_explicit.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_explicit.hpp (original)
+++ sandbox/odeint/branches/karsten/boost/numeric/odeint/stepper/dense_output_explicit.hpp 2011-03-02 02:33:04 EST (Wed, 02 Mar 2011)
@@ -123,12 +123,21 @@
return std::make_pair( m_t_old , m_dt );
}
+ /*
+ * The next two overloads are needed to solve the forwarding problem
+ */
template< class StateOut >
void calc_state( const time_type &t , StateOut &x )
{
m_stepper.calc_state( x , t , *m_old_state , m_t_old , *m_current_state , m_t );
}
+ template< class StateOut >
+ void calc_state( const time_type &t , const StateOut &x )
+ {
+ m_stepper.calc_state( x , t , *m_old_state , m_t_old , *m_current_state , m_t );
+ }
+
template< class StateType >
void adjust_size( const StateType &x )
{
Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/regression_test/integrate_functions.cpp
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/regression_test/integrate_functions.cpp (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/regression_test/integrate_functions.cpp 2011-03-02 02:33:04 EST (Wed, 02 Mar 2011)
@@ -69,7 +69,7 @@
state_type x1;
vector_type x2( 3 );
- integrate( implicit_euler< double >() , make_pair( lorenz() , lorenz_jacobi() ) , x2 , 0.0 , 10.0 , 0.1 , do_nothing_observer() );
+// integrate( implicit_euler< double >() , make_pair( lorenz() , lorenz_jacobi() ) , x2 , 0.0 , 10.0 , 0.1 , do_nothing_observer() );
// integrate_n_steps( implicit_euler< double >() , make_pair( lorenz() , lorenz_jacobi() ) , x2 , 0.0 , 1000 , 0.1 );
// integrate_adaptive( implicit_euler< double >() , make_pair( lorenz() , lorenz_jacobi() ) , x2 , 0.0 , 10.0 , 0.1 );
Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/test/stepper_with_ranges.cpp
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/test/stepper_with_ranges.cpp (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/test/stepper_with_ranges.cpp 2011-03-02 02:33:04 EST (Wed, 02 Mar 2011)
@@ -14,12 +14,15 @@
#include <utility>
#include <boost/range.hpp>
+#include <boost/ref.hpp>
#include <boost/numeric/odeint/stepper/explicit_euler.hpp>
#include <boost/numeric/odeint/stepper/explicit_error_rk54_ck.hpp>
#include <boost/numeric/odeint/stepper/explicit_error_dopri5.hpp>
#include <boost/numeric/odeint/stepper/controlled_error_stepper.hpp>
#include <boost/numeric/odeint/stepper/symplectic_euler.hpp>
+#include <boost/numeric/odeint/stepper/dense_output_explicit.hpp>
+#include <boost/numeric/odeint/stepper/dense_output_controlled_explicit_fsal.hpp>
typedef std::vector< double > state_type;
typedef std::tr1::array< double , 3 > state_type2;
@@ -83,10 +86,37 @@
};
+/*
+ * Useful for Hamiltonian systems
+ */
+struct ham_sys
+{
+ template< class State , class Deriv >
+ void operator()( const State &x_ , Deriv &dxdt_ )
+ {
+ typename boost::range_iterator< Deriv >::type dxdt = boost::begin( dxdt_ );
+ dxdt[0] = 1.0;
+ dxdt[1] = 2.0;
+ dxdt[2] = 3.0;
+ }
+
+ template< class State , class Deriv >
+ void operator()( const State &x_ , const Deriv &dxdt_ )
+ {
+ typename boost::range_iterator< Deriv >::type dxdt = boost::begin( dxdt_ );
+ dxdt[0] = 1.0;
+ dxdt[1] = 2.0;
+ dxdt[2] = 3.0;
+ }
+};
+
+
struct vector_fixture
{
const static size_t dim = 6;
std::tr1::array< double , dim > in;
+ std::tr1::array< double , dim > q;
+ std::tr1::array< double , dim > p;
state_type err;
vector_fixture( void )
@@ -95,7 +125,7 @@
{
for( size_t i=0 ; i<dim ; ++i )
{
- in[ i ] = double( i );
+ in[ i ] = q[i] = p[i] = double( i );
}
for( size_t i=0 ; i<3 ; ++i )
{
@@ -183,16 +213,52 @@
{
vector_fixture f;
boost::numeric::odeint::symplectic_euler< state_type > euler;
-// euler.do_step( )
+ euler.do_step( ham_sys() ,
+ std::make_pair( f.q.begin() + 1 , f.q.begin() + 4 ) ,
+ std::make_pair( f.p.begin() + 3 , f.p.begin() + 6 ) ,
+ 0.0 , 0.1 );
+ CHECK_VALUES( f.q , 0.0 , 1.3 , 2.4 , 3.5 , 4.0 , 5.0 );
+ CHECK_VALUES( f.p , 0.0 , 1.0 , 2.0 , 3.1 , 4.2 , 5.3 );
}
BOOST_AUTO_TEST_CASE( symplectic_euler_coor_and_mom_func )
{
vector_fixture f;
boost::numeric::odeint::symplectic_euler< state_type > euler;
-// euler.do_step( )
+ euler.do_step( std::make_pair( ham_sys() , ham_sys() ) ,
+ std::make_pair( f.q.begin() + 1 , f.q.begin() + 4 ) ,
+ std::make_pair( f.p.begin() + 3 , f.p.begin() + 6 ) ,
+ 0.0 , 0.1 );
+ CHECK_VALUES( f.q , 0.0 , 1.1 , 2.2 , 3.3 , 4.0 , 5.0 );
+ CHECK_VALUES( f.p , 0.0 , 1.0 , 2.0 , 3.1 , 4.2 , 5.3 );
}
+BOOST_AUTO_TEST_CASE( dense_output_euler_with_ranges )
+{
+ using namespace boost::numeric::odeint;
+ vector_fixture f;
+ dense_output_explicit< explicit_euler< state_type > > stepper;
+ stepper.initialize( std::make_pair( f.in.begin() + 1, f.in.begin() + 4 ) , 0.0 , 0.1 );
+ stepper.do_step( system1() );
+ stepper.calc_state( 0.05 , std::make_pair( f.in.begin() + 1 ,f.in.begin() +4 ) );
+ CHECK_VALUES( f.in , 0.0 , 1.05 , 2.1 , 3.15 , 4.0 , 5.0 );
+}
+
+BOOST_AUTO_TEST_CASE( dense_output_dopri5_with_ranges )
+{
+ using namespace boost::numeric::odeint;
+ vector_fixture f;
+ dense_output_controlled_explicit_fsal<
+ controlled_error_stepper<
+ explicit_error_dopri5< state_type >
+ > > stepper;
+ stepper.initialize( std::make_pair( f.in.begin() + 1, f.in.begin() + 4 ) , 0.0 , 0.1 );
+ stepper.do_step( system2() );
+ stepper.calc_state( 0.05 , std::make_pair( f.in.begin() + 1 ,f.in.begin() +4 ) );
+ CHECK_VALUES( f.in , 0.0 , 1.05 , 2.1 , 3.15 , 4.0 , 5.0 );
+}
+
+
BOOST_AUTO_TEST_SUITE_END()
Boost-Commit list run by bdawes at acm.org, david.abrahams at rcn.com, gregod at cs.rpi.edu, cpdaniel at pacbell.net, john at johnmaddock.co.uk