diff --git a/boost/numeric/odeint.hpp b/boost/numeric/odeint.hpp index d3b3c569..49903028 100644 --- a/boost/numeric/odeint.hpp +++ b/boost/numeric/odeint.hpp @@ -33,6 +33,7 @@ #include #include #include +#include /* * Including this algebra slows down the compilation time diff --git a/boost/numeric/odeint/stepper/rosenbrock4.hpp b/boost/numeric/odeint/stepper/rosenbrock4.hpp index 8b4527d9..82bfa836 100644 --- a/boost/numeric/odeint/stepper/rosenbrock4.hpp +++ b/boost/numeric/odeint/stepper/rosenbrock4.hpp @@ -162,6 +162,7 @@ public: rosenbrock4& operator=( const rosenbrock4 &rb ) { copy( rb ); + return *this; } template< class System > diff --git a/boost/numeric/odeint/stepper/rosenbrock4_controller.hpp b/boost/numeric/odeint/stepper/rosenbrock4_controller.hpp index 3fbdcce1..fa860bf9 100644 --- a/boost/numeric/odeint/stepper/rosenbrock4_controller.hpp +++ b/boost/numeric/odeint/stepper/rosenbrock4_controller.hpp @@ -21,6 +21,25 @@ namespace odeint { template< class Stepper > class rosenbrock4_controller { +private: + + void initialize_variables( void ) + { + m_state_adjuster.register_state( 0 , m_xerr ); + } + + void copy_variables( const rosenbrock4_controller &rb ) + { + m_stepper = rb.m_stepper; + m_xerr = rb.m_xerr; + m_atol = rb.m_atol; + m_rtol = rb.m_rtol; + m_first_step = rb.m_first_step; + m_err_old = rb.m_err_old; + m_dt_old = rb.m_dt_old; + m_last_rejected = rb.m_last_rejected; + } + public: typedef Stepper stepper_type; @@ -28,14 +47,33 @@ public: typedef typename stepper_type::state_type state_type; typedef typename stepper_type::time_type time_type; typedef typename stepper_type::deriv_type deriv_type; + typedef typename stepper_type::adjust_size_policy adjust_size_policy; typedef controlled_stepper_tag stepper_category; rosenbrock4_controller( value_type atol = 1.0e-6 , value_type rtol = 1.0e-6 , const stepper_type &stepper = stepper_type() ) - : m_stepper() , m_atol( atol ) , m_rtol( rtol ) , + : m_stepper() , m_state_adjuster() , m_xerr() , + m_atol( atol ) , m_rtol( rtol ) , m_first_step( true ) , m_err_old( 0.0 ) , m_dt_old( 0.0 ) , m_last_rejected( false ) { + initialize_variables(); + } + + rosenbrock4_controller( const rosenbrock4_controller &rb ) + : m_stepper() , m_state_adjuster() , m_xerr() , + m_atol( 1.0e-6l ) , m_rtol( 1.0e-6 ) , + m_first_step( true ) , m_err_old( 0.0 ) , m_dt_old( 0.0 ) , + m_last_rejected( false ) + { + initialize_variables(); + copy_variables( rb ); + } + + rosenbrock4_controller& operator=( const rosenbrock4_controller &rb ) + { + copy_variables( rb ); + return *this; } value_type error( const state_type &x , const state_type &xold , const state_type &xerr ) @@ -82,10 +120,10 @@ public: { static const value_type safe = 0.9 , fac1 = 5.0 , fac2 = 1.0 / 6.0; - const size_t n = x.size(); - state_type xerr( n ); - m_stepper.do_step( sys , x , t , xout , dt , xerr ); - value_type err = error( xout , x , xerr ); + m_state_adjuster.adjust_size_by_policy( x , adjust_size_policy() ); + + m_stepper.do_step( sys , x , t , xout , dt , m_xerr ); + value_type err = error( xout , x , m_xerr ); value_type fac = std::max( fac2 ,std::min( fac1 , std::pow( err , 0.25 ) / safe ) ); value_type dt_new = dt / fac; @@ -121,6 +159,15 @@ public: } + template< class StateType > + void adjust_size( const StateType &x ) + { + m_stepper.adjust_size( x ); + m_state_adjuster.adjust_size( x ); + } + + + stepper_type& stepper( void ) @@ -139,6 +186,8 @@ public: private: stepper_type m_stepper; + size_adjuster< state_type , 1 > m_state_adjuster; + state_type m_xerr; value_type m_atol , m_rtol; bool m_first_step; value_type m_err_old , m_dt_old; diff --git a/boost/numeric/odeint/stepper/rosenbrock4_dense_output.hpp b/boost/numeric/odeint/stepper/rosenbrock4_dense_output.hpp index ca70cbff..863431e5 100644 --- a/boost/numeric/odeint/stepper/rosenbrock4_dense_output.hpp +++ b/boost/numeric/odeint/stepper/rosenbrock4_dense_output.hpp @@ -79,10 +79,6 @@ public: return *this; } - ~rosenbrock4_dense_output( void ) - { - } - template< class StateType > @@ -131,6 +127,15 @@ public: } + template< class StateType > + void adjust_size( const StateType &x ) + { + m_stepper.adjust_size( x ); + m_state_adjuster.adjust_size( x ); + } + + + const state_type& current_state( void ) const { diff --git a/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/concepts.html b/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/concepts.html index f0087739..df10ca1d 100644 --- a/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/concepts.html +++ b/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/concepts.html @@ -3,8 +3,8 @@ Concepts - - + + @@ -13,9 +13,9 @@

-PrevUpHomeNext +PrevUpHomeNext
-
+
@@ -35,7 +35,7 @@

The odeint library defines three concepts for stepping objects.

-
+

Basic stepper @@ -59,12 +59,8 @@ -

-

-

-

@@ -147,25 +143,27 @@

Methods

-
    -
  • -Stepper() - Constructor. -
  • -
  • -Stepper( - container_type &x ) Constructor - that allocates internal memory to store intermediate results of the same - size as x. -
  • -
  • void do_step( DynamicalSystem - &system - , container_type - &x - , time_type - t , - time_type dt - )
  • +
      +
    • + Stepper() + Constructor. +
    • +
    • + Stepper( + container_type &x ) + Constructor that allocates internal memory to store intermediate results + of the same size as x. +
    • +
    • + void do_step( DynamicalSystem + &system + , container_type + &x + , time_type + t , + time_type dt + ) +

    Executes one timestep with the given parameters: @@ -268,44 +266,53 @@ The result of this method is the (approximate) state of the system x(t+dt) and is stored in the variable x (in-place). Note, that the time t is not automatically increased by this method.

    -
    • void do_step( DynamicalSystem - &system - , container_type - &x - , const - container_type &dxdt , time_type t - , time_type - dt )
    +
    • + void do_step( DynamicalSystem + &system + , container_type + &x + , const + container_type &dxdt , time_type t + , time_type + dt ) +

    The same as above but with the additional parameter dxdt that represents the derivative x'(t) = f(x,t) at the time t.

    -
      -
    • -void adjust_size( const container_type &x ) Adjusts - the internal memory to store intermediate results of the same size as - x. This function must - be called whenever the system size changes during the integration. -
    • -
    • -order_type order_step() Returns the order of the algorithm. If - n is the order of a method, then the result - of one iteration with the timestep dt - is accurate up to dt^n. That means the - error made by the time discretization is of order dt^(n+1). -
    • +
        +
      • + void adjust_size( const container_type &x ) + Adjusts the internal memory to store intermediate results of the same + size as x. This function + must be called whenever the system size changes + during the integration. +
      • +
      • + order_type order_step() Returns the order of the algorithm. + If n is the order of a method, then + the result of one iteration with the timestep dt + is accurate up to dt^n. That means the + error made by the time discretization is of order dt^(n+1). +

      Stepper that model this concept

      -
        -
      • stepper_euler
      • -
      • stepper_rk4
      • -
      • stepper_rk78_fehlberg
      • +
          +
        • + stepper_euler +
        • +
        • + stepper_rk4 +
        • +
        • + stepper_rk78_fehlberg +
      -
      +

      Error stepper @@ -330,26 +337,28 @@

      Methods

      -
        -
      • -Error_Stepper() - Constructor. -
      • -
      • -Error_Stepper( - container_type &x ) Constructor - that allocates internal memory to store intermediate results of the same - size as x. -
      • -
      • void do_step( DynamicalSystem - &system - , container_type - &x - , time_type - t , - time_type dt - , container_type - &xerr)
      • +
          +
        • + Error_Stepper() + Constructor. +
        • +
        • + Error_Stepper( + container_type &x ) + Constructor that allocates internal memory to store intermediate results + of the same size as x. +
        • +
        • + void do_step( DynamicalSystem + &system + , container_type + &x + , time_type + t , + time_type dt + , container_type + &xerr) +

        Executes one timestep with the given parameters: @@ -471,45 +480,54 @@ Note, that the time t is not automatically increased by this method.

        -
        • void do_step( DynamicalSystem - &system - , container_type - &x - , const - container_type &dxdt , time_type t - , time_type - dt , - container_type &xerr)
        +
        • + void do_step( DynamicalSystem + &system + , container_type + &x + , const + container_type &dxdt , time_type t + , time_type + dt , + container_type &xerr) +

        The same as above but with the additional parameter dxdt that represents the derivative x'(t) = f(x,t) at the time t.

        -
          -
        • -void adjust_size( const container_type &x ) Adjusts - the internal memory to store intermediate results of the same size as - x. This function must - be called whenever the system size changes during the integration. -
        • -
        • -order_type order_error_step() Returns the order of the result x(t+dt) of the algorithm. -
        • -
        • -order_type order_error() Returns the order of the error estimation - of the algorithm. -
        • +
            +
          • + void adjust_size( const container_type &x ) + Adjusts the internal memory to store intermediate results of the same + size as x. This function + must be called whenever the system size changes + during the integration. +
          • +
          • + order_type order_error_step() Returns the order of the result x(t+dt) of the algorithm. +
          • +
          • + order_type order_error() Returns the order of the error estimation + of the algorithm. +

          Stepper that model this concept

          -
            -
          • stepper_rk5_ck
          • -
          • stepper_rk78_fehlberg
          • -
          • stepper_half_step
          • +
              +
            • + stepper_rk5_ck +
            • +
            • + stepper_rk78_fehlberg +
            • +
            • + stepper_half_step +
          -
          +

          Controlled stepper @@ -529,11 +547,13 @@

          Methods

          -
          • Controlled_Stepper( - time_type abs_err, time_type - rel_err, - time_type factor_x, time_type - factor_dxdt )
          +
          • + Controlled_Stepper( + time_type abs_err, time_type + rel_err, + time_type factor_x, time_type + factor_dxdt ) +

          Constructor that initializes the controlled stepper with several parameters of the error control. The controlled stepper assures that the error done @@ -552,22 +572,20 @@ If the estimated error is less than half of the desired error, an increased stepsize will be suggested.

          -
            -
          • -Controlled_Stepper( - container_type &x, time_type - abs_err, - time_type rel_err, time_type - factor_x, - time_type factor_dxdt - ) Same as above, but with additional - allocation of the internal memory to store intermediate results of the - same size as x. -
          • -
          • controlled_step_result try_step( DynamicalSystem &system, container_type &x, time_type - &t, time_type - &dt - )
          • +
              +
            • + Controlled_Stepper( + container_type &x, time_type abs_err, time_type + rel_err, + time_type factor_x, time_type + factor_dxdt ) + Same as above, but with additional allocation of the internal memory + to store intermediate results of the same size as x. +
            • +
            • + controlled_step_result try_step( + DynamicalSystem &system, container_type &x, time_type &t, time_type &dt ) +

            Tries one timestep with the given parameters @@ -680,43 +698,50 @@ dt now containes the suggested reduced stepsize that should give an error below the desired level.

            -
              -
            • -controlled_step_result try_step( DynamicalSystem &system, container_type &x, const - container_type &dxdt, time_type &t, time_type - &dt - ) Same as above but with the additional - parameter dxdt that that - represents the derivative x'(t) = f(x,t) - at the time t. -
            • -
            • -void adjust_size( const container_type &x ) Adjusts - the internal memory to store intermediate results of the same size as - x. This function must - be called whenever the system size changes during the integration. -
            • -
            • -order_type order_error_step() Returns the order of the result x(t+dt) of the algorithm. -
            • +
                +
              • + controlled_step_result try_step( + DynamicalSystem &system, container_type &x, const container_type + &dxdt, time_type + &t, time_type + &dt + ) Same as above but with the additional + parameter dxdt that that + represents the derivative x'(t) = f(x,t) + at the time t. +
              • +
              • + void adjust_size( const container_type &x ) + Adjusts the internal memory to store intermediate results of the same + size as x. This function + must be called whenever the system size changes + during the integration. +
              • +
              • + order_type order_error_step() Returns the order of the result x(t+dt) of the algorithm. +

              Stepper that model this concept

              -
                -
              • controlled_stepper_standard
              • -
              • controlled_stepper_bs
              • +
                  +
                • + controlled_stepper_standard +
                • +
                • + controlled_stepper_bs +
              -

              + -

              + -

              +

              @@ -733,7 +758,7 @@


              -PrevUpHomeNext +PrevUpHomeNext
              diff --git a/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/extend_odeint.html b/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/extend_odeint.html index f7ff8f38..e3370fad 100644 --- a/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/extend_odeint.html +++ b/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/extend_odeint.html @@ -3,8 +3,8 @@ Extend odeint - - + + @@ -13,9 +13,9 @@

              -PrevUpHomeNext +PrevUpHomeNext
              -
              + -

              + -
              +

              Adapt your own containers @@ -40,7 +40,7 @@ gsl_vector, gsl_matrix, ublas::matrix, blitz::matrix, thrust

              -
              +

              Adapt your own operations @@ -60,7 +60,7 @@
              -PrevUpHomeNext +PrevUpHomeNext
              diff --git a/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/getting_started.html b/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/getting_started.html index 9da67767..337f12bb 100644 --- a/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/getting_started.html +++ b/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/getting_started.html @@ -3,8 +3,8 @@ Getting started - - + + @@ -13,9 +13,9 @@

              -PrevUpHomeNext +PrevUpHomeNext
              -
              +
              @@ -26,7 +26,7 @@
              Short Example
              -
              +
              @@ -52,7 +52,7 @@ algorithms are implemented:

              -

              Table 1.1. Stepper Algorithms

              +

              Table 1.1. Stepper Algorithms

              @@ -331,7 +331,7 @@ sometimes also refered as lattices ODEs.

              -
              +

              Usage, Compilation, Headers @@ -354,7 +354,7 @@

              -
              +

              Short Example @@ -369,9 +369,7 @@ = (x,p):

              -

              -

              - +

              /* The type of container used to hold the state vector */
               typedef std::vector< double > state_type;
              @@ -385,8 +383,6 @@
                   dxdt[1] = -x[0] - gam*x[1];
               }
               
              -

              -

              @@ -403,16 +399,12 @@ start:

              -

              -

              - +

              state_type x(2);
               x[0] = 1.0; // start at x=1.0, p=0.0
               x[1] = 0.0;
               
              -

              -

              @@ -422,15 +414,11 @@ stepper (5th order) and uses adaptive stepsize.

              -

              -

              - +

              size_t steps = integrate( harmonic_oscillator ,
                                         x , 0.0 , 10.0 , 0.1 );
               
              -

              -

              @@ -447,9 +435,7 @@ rhs must then be implemented as a functor having defined the ()-operator:

              -

              -

              - +

              /* The rhs of x' = f(x) defined as a class */
               class harm_osc {
              @@ -466,24 +452,18 @@
                   }
               };
               
              -

              -

              which can be used via

              -

              -

              - +

              harm_osc ho(0.15);
               steps = integrate( ho ,
                                  x , 0.0 , 10.0 , 0.1 );
               
              -

              -

              @@ -493,9 +473,7 @@ to provide a reasonable observer. An example is

              -

              -

              - +

              struct push_back_state_and_time
               {
              @@ -512,8 +490,6 @@
               	}
               };
               
              -

              -

              @@ -521,9 +497,7 @@ pass this container to the integration function:

              -

              -

              - +

              vector<state_type> x_vec;
               vector<double> times;
              @@ -538,8 +512,6 @@
                   cout << times[i] << '\t' << x_vec[i][0] << '\t' << x_vec[i][1] << '\n';
               }
               
              -

              -

              @@ -562,7 +534,7 @@


              -PrevUpHomeNext +PrevUpHomeNext
              diff --git a/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/reference.html b/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/reference.html index 1da49848..af437ee9 100644 --- a/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/reference.html +++ b/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/reference.html @@ -3,8 +3,8 @@ Reference - - + + @@ -12,9 +12,9 @@

              -PrevUpHome +PrevUpHome
              -
              +
              @@ -27,13 +27,13 @@
              Operations
              Resizing
              -
              +
              -

              Table 1.3. Stepper Algorithms

              +

              Table 1.3. Stepper Algorithms

              @@ -201,7 +201,7 @@
              -

              -PrevUpHome +PrevUpHome
              diff --git a/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/tutorial.html b/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/tutorial.html index e4693973..2761c0dd 100644 --- a/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/tutorial.html +++ b/libs/numeric/odeint/doc/html/boost_sandbox_numeric_odeint/tutorial.html @@ -3,8 +3,8 @@ Tutorial - - + + @@ -13,9 +13,9 @@

              -PrevUpHomeNext +PrevUpHomeNext
              -
              +
              @@ -34,7 +34,7 @@ topics
              References
              -
              +
              -
              +

              Define the ODE @@ -71,9 +71,7 @@ to just define a function, e.g:

              -

              -

              - +

              /* The type of container used to hold the state vector */
               typedef std::vector< double > state_type;
              @@ -87,8 +85,6 @@
                   dxdt[1] = -x[0] - gam*x[1];
               }
               
              -

              -

              @@ -102,9 +98,7 @@ parameter structure as above:

              -

              -

              - +

              /* The rhs of x' = f(x) defined as a class */
               class harm_osc {
              @@ -121,8 +115,6 @@
                   }
               };
               
              -

              -

              @@ -130,7 +122,7 @@ which allows for cleaner code.

              -
              +

              Stepper Types @@ -147,7 +139,7 @@ choose:

              -

              Table 1.2. Stepper Algorithms

              +

              Table 1.2. Stepper Algorithms

              @@ -424,7 +416,7 @@ which kind of steppers should be applied.

              -
              +

              Integration with Constant Step Size @@ -442,15 +434,11 @@ ) function from odeint:

              -

              -

              - +

              explicit_rk4< state_type > stepper;
               integrate_const( stepper , harmonic_oscillator , x , 0.0 , 10.0 , 0.01 );
               
              -

              -

              @@ -463,20 +451,16 @@ method which can used directly. So, you write down the above example as

              -

              -

              - +

              const double dt = 0.01;
               for( double t=0.0 ; t<10.0 ; t+= dt )
               	stepper.do_step( harmonic_oscillator , x , t , dt );
               
              -

              -

              -
              +

              Integration with Adaptive Step Size @@ -491,15 +475,11 @@ with 4th order error estimation and coefficients introduced by Cash-Karp.

              -

              -

              - +

              typedef explicit_error_rk54_ck< state_type > error_stepper_type;
               error_stepper_type rk54;
               
              -

              -

              @@ -524,16 +504,12 @@ stepper create by make_controlled_stepper_standard.

              -

              -

              - +

              typedef controlled_error_stepper< error_stepper_type > controlled_stepper_type;
               controlled_stepper_type controlled_stepper;
               integrate_adaptive( controlled_stepper , harmonic_oscillator , x , 0.0 , 10.0 , 0.01 );
               
              -

              -

              @@ -557,7 +533,7 @@ The full cpp file for this example can be found here: ../../examples/harmonic_oscillator.cpp

              -
              +
              -
              +

              Gravitation and energy conservation @@ -630,7 +606,7 @@ f( q​i ).

              -
              +

              Define the system function @@ -640,9 +616,7 @@ space as well as the velocity. Therefore, we use the operators from Boost.Operators:

              -

              -

              - +

              /*the point type */
               template< class T , size_t Dim >
              @@ -672,8 +646,6 @@
               //...
               // more operators
               
              -

              -

              @@ -682,9 +654,7 @@ type we use std::tr1::array

              -

              -

              - +

              // we simulate 5 planets and the sun
               const size_t n = 6;
              @@ -693,8 +663,6 @@
               typedef std::tr1::array< point_type , n > container_type;
               typedef std::tr1::array< double , n > mass_type;
               
              -

              -

              @@ -709,9 +677,7 @@ As system function we have to provide f(p) and f(q):

              -

              -

              - +

              const double gravitational_constant = 2.95912208286e-4;
               
              @@ -728,14 +694,10 @@
               	}
               };
               
              -

              -

              -

              -

              - +

              struct solar_system_momentum
               {
              @@ -762,8 +724,6 @@
               	}
               };
               
              -

              -

              @@ -779,9 +739,7 @@ apply here:

              -

              -

              - +

              typedef symplectic_rkn_sb3a_mclachlan< container_type > stepper_type;
               const double dt = 100.0;
              @@ -792,8 +750,6 @@
               		make_pair( boost::ref( q ) , boost::ref( p ) ) ,
               		0.0 , 200000.0 , dt , streaming_observer( cout ) );
               
              -

              -

              @@ -807,9 +763,7 @@ is also passed, but this is not a problem at all:

              -

              -

              - +

              struct streaming_observer
               {
              @@ -827,8 +781,6 @@
               	}
               };
               
              -

              -

              @@ -836,7 +788,7 @@

              -
              +

              Chaotic systems and Lyapunov exponents @@ -845,7 +797,7 @@ blah blah

              -
              +

              Stiff systems @@ -873,9 +825,7 @@ is needed. Here is the definition of the above example

              -

              -

              - +

              typedef boost::numeric::ublas::vector< double > vector_type;
               typedef boost::numeric::ublas::matrix< double > matrix_type;
              @@ -902,8 +852,6 @@
               	}
               };
               
              -

              -

              @@ -914,16 +862,16 @@ just templatize the operator():

              -

              -

              - +

              -
              struct stiff_system
              +
              typedef boost::numeric::ublas::vector< double > vector_type;
              +typedef boost::numeric::ublas::matrix< double > matrix_type;
              +
              +struct stiff_system
               {
               	template< class State >
               	void operator()( const State &x , State &dxdt , double t )
               	{
              -		...
               	}
               };
               
              @@ -932,12 +880,9 @@
               	template< class State , class Matrix >
               	void operator()( const State &x , Matrix &J , const double &t , State &dfdt )
               	{
              -		...
               	}
               };
               
              -

              -

              @@ -953,34 +898,29 @@ all the other stepper:

              -

              -

              - +

              typedef rosenbrock4< double > stepper_type;
               typedef rosenbrock4_controller< stepper_type > controlled_stepper_type;
              +typedef rosenbrock4_dense_output< controlled_stepper_type > dense_output_type;
               
               vector_type x( 3 , 1.0 );
               
              -size_t num_of_steps = integrate_const( controlled_stepper_type() ,
              +size_t num_of_steps = integrate_const( dense_output_type() ,
               		make_pair( stiff_system() , stiff_system_jacobi() ) ,
               		x , 0.0 , 50.0 , 0.01 ,
               		cout << phoenix::arg_names::arg2 << " " << phoenix::arg_names::arg1[0] << "\n" );
               
              -

              -

              - During the integration approximately XX steps have been done. Comparing to + During the integration approximately 71 steps have been done. Comparing to a classical Runge-Kutta solver this is a very good result. For example the Dormand-Prince 5 method with step size control and dense output yields ca. - XX steps. + 1531 steps.

              -

              -

              - +

              typedef explicit_error_dopri5< vector_type > dopri5_type;
               typedef controlled_error_stepper< dopri5_type > controlled_dopri5_type;
              @@ -991,10 +931,7 @@
               size_t num_of_steps2 = integrate_const( dense_output_dopri5_type() ,
               		stiff_system() , x2 , 0.0 , 50.0 , 0.01 ,
               		cout << phoenix::arg_names::arg2 << " " << phoenix::arg_names::arg1[0] << "\n" );
              -clog << num_of_steps2 << endl;
               
              -

              -

              @@ -1005,7 +942,7 @@ The full example can be found here: ../../examples/stiff_system.cpp

              -
              +

              Odeint in detail @@ -1025,11 +962,14 @@

              dense output

              +

              + Algebras and operations +

              Overview and table

              -
              +
              -
              +

              Complex state types @@ -1073,7 +1013,7 @@ DNLS

              -
              +

              Lattice systems @@ -1082,7 +1022,7 @@ Fermi-Pasta-Ulam system

              -
              +

              Partial differential equations @@ -1097,7 +1037,7 @@ Ginzburg-Landau

              -
              +

              Using boost::units @@ -1106,7 +1046,7 @@ blah blah

              -
              +

              Using Cuda and Thrust @@ -1115,7 +1055,7 @@ blah blah

              -
              +

              Using matrices as state types @@ -1124,11 +1064,11 @@ Expanding resizing

              -

              + -
              +

              Pass by value or by reference @@ -1137,7 +1077,7 @@ blah blah

              -
              +

              Using boost::range @@ -1146,7 +1086,7 @@ blah blah

              -
              +
              -
              -
              +
              @@ -1217,7 +1157,7 @@


              -PrevUpHomeNext +PrevUpHomeNext
              diff --git a/libs/numeric/odeint/doc/html/index.html b/libs/numeric/odeint/doc/html/index.html index 1658d3fb..febd208b 100644 --- a/libs/numeric/odeint/doc/html/index.html +++ b/libs/numeric/odeint/doc/html/index.html @@ -3,15 +3,15 @@ Chapter 1. boost.sandbox.numeric.odeint - - + +

              -
              Next
              -
              +
              Next
              +

              Chapter 1. boost.sandbox.numeric.odeint

              @@ -23,7 +23,7 @@

              -

              +

              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)

              @@ -93,10 +93,10 @@
              - +

              Last revised: March 23, 2011 at 08:22:08 GMT

              Last revised: March 31, 2011 at 16:36:47 GMT


              -
              Next
              +
              Next
              diff --git a/libs/numeric/odeint/doc/tutorial_stiff_systems.qbk b/libs/numeric/odeint/doc/tutorial_stiff_systems.qbk index 176276db..45583479 100644 --- a/libs/numeric/odeint/doc/tutorial_stiff_systems.qbk +++ b/libs/numeric/odeint/doc/tutorial_stiff_systems.qbk @@ -28,7 +28,7 @@ A well know solver for stiff systems is the so called Rosenbrock method. It has [integrate_stiff_system] -During the integration approximately XX steps have been done. Comparing to a classical Runge-Kutta solver this is a very good result. For example the Dormand-Prince 5 method with step size control and dense output yields ca. XX steps. +During the integration approximately 71 steps have been done. Comparing to a classical Runge-Kutta solver this is a very good result. For example the Dormand-Prince 5 method with step size control and dense output yields ca. 1531 steps. [integrate_stiff_system_alternative] diff --git a/libs/numeric/odeint/examples/stiff_system.cpp b/libs/numeric/odeint/examples/stiff_system.cpp index 14f78ddd..646bcf98 100644 --- a/libs/numeric/odeint/examples/stiff_system.cpp +++ b/libs/numeric/odeint/examples/stiff_system.cpp @@ -79,10 +79,11 @@ int main( int argc , char **argv ) //[ integrate_stiff_system typedef rosenbrock4< double > stepper_type; typedef rosenbrock4_controller< stepper_type > controlled_stepper_type; + typedef rosenbrock4_dense_output< controlled_stepper_type > dense_output_type; vector_type x( 3 , 1.0 ); - size_t num_of_steps = integrate_const( controlled_stepper_type() , + size_t num_of_steps = integrate_const( dense_output_type() , make_pair( stiff_system() , stiff_system_jacobi() ) , x , 0.0 , 50.0 , 0.01 , cout << phoenix::arg_names::arg2 << " " << phoenix::arg_names::arg1[0] << "\n" ); @@ -101,8 +102,9 @@ int main( int argc , char **argv ) size_t num_of_steps2 = integrate_const( dense_output_dopri5_type() , stiff_system() , x2 , 0.0 , 50.0 , 0.01 , cout << phoenix::arg_names::arg2 << " " << phoenix::arg_names::arg1[0] << "\n" ); - clog << num_of_steps2 << endl; //] + clog << num_of_steps2 << endl; + return 0; }