Boost logo

Boost-Commit :

Subject: [Boost-commit] svn:boost r66772 - in sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas: butcher generic_stepper
From: karsten.ahnert_at_[hidden]
Date: 2010-11-26 11:49:47


Author: karsten
Date: 2010-11-26 11:49:44 EST (Fri, 26 Nov 2010)
New Revision: 66772
URL: http://svn.boost.org/trac/boost/changeset/66772

Log:
* generic stepper typedefs and type construction
Text files modified:
   sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/butcher/algebra.hpp | 12
   sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/butcher/butcher_performance.cpp | 4
   sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/generic_stepper/Jamfile | 6
   sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/generic_stepper/rk_test.cpp | 67 ++--------
   sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/generic_stepper/runge_kutta_stepper.hpp | 241 +++++++++++++++++++++++----------------
   5 files changed, 165 insertions(+), 165 deletions(-)

Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/butcher/algebra.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/butcher/algebra.hpp (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/butcher/algebra.hpp 2010-11-26 11:49:44 EST (Fri, 26 Nov 2010)
@@ -38,14 +38,14 @@
         {
                 const static double a1 = dt * convert_value< typename mpl::at< coef_type , mpl::int_< 0 > >::type >::get_value();
 
- iterator first1 = x_tmp.begin() , last1 = x_tmp.end();
- const_iterator first2 = x.begin() , first3 = k_vector[0].begin();
- while( first1 != last1 )
- *first1++ = *first2++ + a1 * *first3++ ;
+// iterator first1 = x_tmp.begin() , last1 = x_tmp.end();
+// const_iterator first2 = x.begin() , first3 = k_vector[0].begin();
+// while( first1 != last1 )
+// *first1++ = *first2++ + a1 * *first3++ ;
 
 
-// std_algebra::for_each3( x_tmp , x , k_vector[0] ,
-// std_op::scale_sum2( 1.0 , a1 ) );
+ std_algebra::for_each3( x_tmp , x , k_vector[0] ,
+ std_op::scale_sum2( 1.0 , a1 ) );
 
 
 // for( size_t i=0 ; i<x.size() ; ++i )

Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/butcher/butcher_performance.cpp
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/butcher/butcher_performance.cpp (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/butcher/butcher_performance.cpp 2010-11-26 11:49:44 EST (Fri, 26 Nov 2010)
@@ -37,8 +37,8 @@
 
 
 typedef std::tr1::array< double , 3 > state_type;
-typedef mpl_rk4_stepper< state_type > stepper_type;
-typedef boost::numeric::odeint::explicit_rk4< state_type > stepper_std_type;
+typedef mpl_euler_stepper< state_type > stepper_type;
+typedef boost::numeric::odeint::explicit_euler< state_type > stepper_std_type;
 
 
 void lorenz( const state_type &x , state_type &dxdt , double t )

Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/generic_stepper/Jamfile
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/generic_stepper/Jamfile (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/generic_stepper/Jamfile 2010-11-26 11:49:44 EST (Fri, 26 Nov 2010)
@@ -15,9 +15,9 @@
       <include>$(BOOST_ROOT)
     ;
 
-exe test
- : test.cpp
- ;
+#exe test
+# : test.cpp
+# ;
 
 exe rk_test
     : rk_test.cpp

Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/generic_stepper/rk_test.cpp
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/generic_stepper/rk_test.cpp (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/generic_stepper/rk_test.cpp 2010-11-26 11:49:44 EST (Fri, 26 Nov 2010)
@@ -1,58 +1,21 @@
-/*
- * rk_test.cpp
- *
- * Created on: Nov 23, 2010
- * Author: mario
- */
-
-#include <iostream>
-#include <vector>
-#include <tr1/array>
-#include <boost/array.hpp>
-#include <boost/fusion/container.hpp>
-
 #include "runge_kutta_stepper.hpp"
 
-using namespace std;
-
-typedef tr1::array< double , 3 > state_type;
-typedef runge_kutta_stepper< state_type , 1 > euler_stepper;
-typedef runge_kutta_stepper< state_type , 2 > midpoint_stepper;
-
-
-const double sigma = 10.0;
-const double R = 28.0;
-const double b = 8.0 / 3.0;
-
-void lorenz( const state_type &x , state_type &dxdt , double t )
+int main( int argc , char **argv )
 {
- dxdt[0] = sigma * ( x[1] - x[0] );
- dxdt[1] = R * x[0] - x[1] - x[0] * x[2];
- dxdt[2] = x[0]*x[1] - b * x[2];
-}
-
-const double dt = 0.001;
-
-int main( void )
-{
-
- boost::array< double , 1 > b = {{ 1.0 }};
- boost::array< double , 1 > c = {{ 0.0 }};
-
- euler_stepper euler( fusion::vector0<>() , b , c );
- euler.print_vals();
- cout << typeid(euler_stepper::stage_vector_base).name() << endl;
+ typedef runge_kutta_stepper< 4 > rk_type;
+ typedef rk_type::coef_a_type coef_a_type;
+ typedef rk_type::coef_b_type coef_b_type;
+ typedef rk_type::coef_c_type coef_c_type;
+
+ const boost::array< double , 1 > a1 = {{ 5.1 }};
+ const boost::array< double , 2 > a2 = {{ 5.2 , 5.3 }};
+ const boost::array< double , 3 > a3 = {{ 5.4 , 5.5 , 5.6 }};
+
+ const coef_a_type a = fusion::make_vector( a1 , a2 , a3 );
+ const coef_b_type b = {{ 6.1 , 6.2 , 6.3 , 6.4 }};
+ const coef_c_type c = {{ 1.1 , 1.2 , 1.3 , 1.4 }};
 
-// boost::array< double , 1 > a2 = {{ 0.5 }};
-// boost::array< double , 2 > b2 = {{ 0.0 , 1.0 }};
-// boost::array< double , 2 > c2 = {{ 0.0 , 0.5 }};
-//
-// cout << typeid(midpoint_stepper::coef_a_type).name() << endl;
-// cout << typeid(midpoint_stepper::coef_b_type).name() << endl;
-// cout << typeid(midpoint_stepper::coef_c_type).name() << endl;
-// cout << typeid(midpoint_stepper::stage_vector_base).name() << endl;
-//
-// midpoint_stepper midpoint( fusion::make_vector( a2 ) , b2 , c2 );
-// midpoint.print_vals();
+ rk_type rk( a , b , c );
+ rk.print_vals();
 
 }

Modified: sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/generic_stepper/runge_kutta_stepper.hpp
==============================================================================
--- sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/generic_stepper/runge_kutta_stepper.hpp (original)
+++ sandbox/odeint/branches/karsten/libs/numeric/odeint/ideas/generic_stepper/runge_kutta_stepper.hpp 2010-11-26 11:49:44 EST (Fri, 26 Nov 2010)
@@ -1,4 +1,3 @@
-
 #include <boost/mpl/vector.hpp>
 #include <boost/mpl/size.hpp>
 #include <boost/mpl/push_back.hpp>
@@ -9,49 +8,64 @@
 #include <boost/mpl/at.hpp>
 #include <boost/mpl/copy.hpp>
 #include <boost/mpl/insert_range.hpp>
-
-#include <vector>
-#include <algorithm>
-#include <iostream>
+#include <boost/mpl/size_t.hpp>
 
 #include <boost/fusion/container.hpp>
 #include <boost/fusion/algorithm.hpp>
-#include <boost/fusion/include/mpl.hpp>
-#include <boost/fusion/view/zip_view.hpp>
-#include <boost/fusion/include/zip_view.hpp>
-#include <boost/fusion/algorithm/transformation/push_back.hpp>
-#include <boost/fusion/include/push_back.hpp>
+#include <boost/fusion/sequence.hpp>
+#include <boost/fusion/adapted.hpp>
 
 
-#include <boost/array.hpp>
+#include <vector>
+#include <algorithm>
+#include <iostream>
+#include <string>
 
+#include <boost/array.hpp>
 #include <typeinfo>
 
 
-#include "algebra.hpp"
-
 namespace mpl = boost::mpl;
 namespace fusion = boost::fusion;
 
 using namespace std;
 
 
+struct first_stage {};
+struct intermediate_stage {};
+struct last_stage {};
+
+ostream& operator<<( ostream& out , const first_stage& )
+{
+ out << "first stage";
+ return out;
+}
+
+ostream& operator<<( ostream& out , const intermediate_stage& )
+{
+ out << "intermediate stage";
+ return out;
+}
+
+ostream& operator<<( ostream& out , const last_stage& )
+{
+ out << "last stage";
+ return out;
+}
+
+
+
 template< class T , class Constant >
 struct array_wrapper
 {
     typedef typename boost::array< T , Constant::value > type;
 };
 
-struct first_stage {};
-struct intermediate_stage {};
-struct last_stage {};
 
 template< typename T , size_t N , typename StageCategory >
 struct stage_fusion
 {
- //typedef typename fusion::result_of::as_vector< mpl::vector< size_t , T , boost::array< T , N > , StageCategory > >::type type;
     typedef fusion::vector4< size_t , T , boost::array< T , N > , StageCategory > type;
- typedef fusion::vector4< size_t , T , boost::array< T , N > , StageCategory >& ref_type;
 };
 
 template< class T , class Constant , class StageCategory >
@@ -72,14 +86,56 @@
     }
 };
 
+struct print_sequence
+{
+ size_t m_indent;
+
+ print_sequence( size_t indent = 0 )
+ : m_indent( indent ) { }
 
-/* Runge Kutta Stepper - consisting of several stages */
+ template< typename T >
+ void operator()( T ) const
+ {
+ string in;
+ for( size_t i=0 ; i<m_indent ; ++i ) in += " ";
+ cout << in << typeid(T).name() << endl;
+ }
+};
 
-template< typename State , size_t stage_count >
+struct print_values
+{
+ size_t m_indent;
+
+ print_values( size_t indent = 0 )
+ : m_indent( indent ) { }
+
+ template< typename T >
+ void operator()( const T &t ) const
+ {
+ typedef typename fusion::traits::is_sequence< T >::type seq_type;
+ print_val( t , seq_type() );
+ }
+
+ template< class R >
+ void print_val( const R &r , mpl::false_ ) const
+ {
+ string in;
+ for( size_t i=0 ; i<m_indent ; ++i ) in += " ";
+ cout << in << r << endl;
+ }
+
+ template< class S >
+ void print_val( const S &s , mpl::true_ ) const
+ {
+ fusion::for_each( s , print_values( m_indent + 1 ) );
+ }
+};
+
+
+template< size_t stage_count >
 class runge_kutta_stepper
 {
 public:
- typedef State state_type;
 
     typedef mpl::range_c< size_t , 1 , stage_count > stage_indices;
 
@@ -97,89 +153,60 @@
>::type coef_a_type;
 
     typedef boost::array< double , stage_count > coef_b_type;
-
     typedef boost::array< double , stage_count > coef_c_type;
 
 
- typedef //typename fusion::result_of::as_vector<
- typename mpl::copy<
+ typedef typename mpl::copy
+ <
             stage_indices,
             mpl::inserter
- <
+ <
                 mpl::vector0<> ,
                 mpl::push_back< mpl::_1 , stage_fusion_wrapper< double , mpl::_2 , intermediate_stage > >
- >
- >::type stage_vector_base0;
- //>::type
-
- typedef typename fusion::result_of::as_vector<
- typename mpl::push_back<
- stage_vector_base0 ,
- typename stage_fusion_wrapper< double , mpl::integral_c< size_t , stage_count > , last_stage >::type
+ >
+ >::type base_of_stage_vector_base;
+
+ typedef typename fusion::result_of::as_vector
+ <
+ typename mpl::push_back
+ <
+ base_of_stage_vector_base ,
+ typename stage_fusion_wrapper< double , mpl::size_t< stage_count > , last_stage >::type
>::type
>::type stage_vector_base;
 
-
- struct stage_vector_type : public stage_vector_base
+ struct stage_vector : public stage_vector_base
     {
- struct init_stage{
-
- const coef_c_type &m_c;
-
- init_stage( const coef_b_type &b , const coef_c_type &c )
- : m_c( c ) {}
-
-
- template< typename T , size_t N , typename StageCategory>
- void operator() (
- fusion::vector2< fusion::vector4< size_t , T , boost::array< T , N > , StageCategory >& ,
- //stage_fusion< T , N , StageCategory >::ref_type ,
- boost::array< T , N > > stage ) const
- {
- std::cout << N << std::endl;
- fusion::at_c<0>( fusion::at_c<0>(stage) ) = N;
- fusion::at_c<1>( fusion::at_c<0>(stage) ) = m_c[N-1];
- fusion::at_c<2>( fusion::at_c<0>(stage) ) = fusion::at_c<1>( stage );
- }
-
- };
-
-
- stage_vector_type( const coef_a_type &a , const coef_b_type &b , const coef_c_type &c )
- {
- //typedef typename fusion::result_of::push_back< const coef_a_type , const coef_b_type >::type ab_type;
- typedef fusion::joint_view< const coef_a_type , const fusion::single_view< coef_b_type > > ab_type;
- const ab_type ab = fusion::push_back( a , b );
- // iterate stage and c in parallel
- typedef fusion::vector< stage_vector_base& , const ab_type& > zipped_stage;
-
- fusion::for_each( fusion::zip_view<zipped_stage>( zipped_stage( static_cast< stage_vector_base& >( *this ) , ab ) ) ,
- init_stage( b , c ) );
- }
+ struct do_insertion
+ {
+ stage_vector_base &m_base;
+ const coef_a_type &m_a;
+ const coef_c_type &m_c;
+
+ do_insertion( stage_vector_base &base , const coef_a_type &a , const coef_c_type &c )
+ : m_base( base ) , m_a( a ) , m_c( c ) { }
+
+ template< class Index >
+ void operator()( Index ) const
+ {
+ fusion::at_c< 0 >( fusion::at< Index >( m_base ) ) = Index::value;
+ fusion::at_c< 1 >( fusion::at< Index >( m_base ) ) = m_c[ Index::value ];
+ fusion::at_c< 2 >( fusion::at< Index >( m_base ) ) = fusion::at< Index >( m_a );
+ }
+ };
+
+ stage_vector( const coef_a_type &a , const coef_b_type &b , const coef_c_type &c )
+ {
+ typedef mpl::range_c< size_t , 0 , stage_count - 1 > indices;
+ mpl::for_each< indices >( do_insertion( *this , a , c ) );
+ fusion::at_c< 0 >( fusion::at_c< stage_count - 1 >( *this ) ) = stage_count - 1 ;
+ fusion::at_c< 1 >( fusion::at_c< stage_count - 1 >( *this ) ) = c[ stage_count - 1 ];
+ fusion::at_c< 2 >( fusion::at_c< stage_count - 1 >( *this ) ) = b;
+ }
     };
 
 
-// template< class System >
-// struct calculate_stage
-// {
-// System &system;
-// state_type &x , &x_tmp;
-// state_type *k_vector;
-// const double t;
-// const double dt;
-//
-// calculate_stage( System &_system , state_type &_x , state_type &_x_tmp , state_type *_k_vector ,
-// const double _t , const double _dt )
-// : system( _system ) , x( _x ) , x_tmp( _x_tmp ) , k_vector( _k_vector ) , t( _t ) , dt( _dt )
-// {}
-//
-// template< class Stage >
-// void operator()( Stage const &stage ) const
-// {
-// stage( system , x_tmp , x , k_vector , t , dt );
-// }
-//
-// };
+
 
 
 public:
@@ -187,26 +214,36 @@
     runge_kutta_stepper( const coef_a_type &a ,
                             const coef_b_type &b ,
                             const coef_c_type &c )
- : m_stages( a , b , c )
+ : m_a( a ) , m_b( b ) , m_c( c ) ,
+ m_stages( a , b , c )
+
     { }
 
- template< class System >
- void do_step( System &system , state_type &x , double t , const double dt )
- {
-// fusion::for_each( m_stages , calculate_stage< System >( system , x , m_x_tmp , m_k_vector , t , dt ) );
-// m_last_stage( system , m_x_tmp , x , m_k_vector , t , dt );
- }
 
     void print_vals()
     {
- fusion::for_each( m_stages , print_stage() );
+ cout << "Typeof state_vector_base : " << endl;
+ mpl::for_each< stage_vector_base >( print_sequence( 1 ) );
+
+ cout << "m_a : " << endl;
+ fusion::for_each( m_a , print_values( 1 ) );
+
+ cout << "m_b : " << endl;
+ fusion::for_each( m_b , print_values( 1 ) );
+
+ cout << "m_c : " << endl;
+ fusion::for_each( m_c , print_values( 1 ) );
+
+ cout << "m_stages : " << endl;
+ fusion::for_each( m_stages , print_values( 1 ) );
     }
 
 
+
 private:
 
- state_type m_k_vector[stage_count];
- state_type m_x_tmp;
- const stage_vector_type m_stages;
-// last_stage_type m_last_stage;
+ const coef_a_type m_a;
+ const coef_b_type m_b;
+ const coef_c_type m_c;
+ const stage_vector m_stages;
 };


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