[boost] library proposal: odeint

Subject: [boost] library proposal: odeint
From: Karsten Ahnert (karsten.ahnert_at_[hidden])
Date: 2009-12-13 10:02:34

Hi,

I am new on this list and want to start a discussion if a library for
solving ordinary differential equations (ODEs) would fit into boost.

We are developing a lib - odeint - for such problems. odeint can be
found in the sandbox: https://boost.org/svn/boost/sandbox/odeint, but
the current stage is far away from being ready for release or submission.

ODEs are usually solved iteratively. Therefore, we have introduced two
stepper concepts, one for steppers with- and one for steppers without
error estimation. At the moment some explicit solvers are implemented:
Euler, RK4, Cash Karp, RK78,midpoint and Burlisch-Stoer.

Here are some examples:

typedef std::tr1::array< double , 3 > state_type;

// define the ODE
void lorenz( state_type &x , state_type &dxdt , double t )
{
    const double sigma = 10.0 , R = 28.0 , b = 8.0 / 3.0
    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];
}

state_type x = {{ 1.0 , 0.0 , 0.0 }};
ode_step_euler< state_type > stepper;
double t = 0.0;

// iterate the ODE by hand with constant stepsize
for( size_t oi=0 ; oi<10000 ; ++oi,t+=dt )
    stepper.next_step( lorenz , x , t , 0.01 );

// integrate the ODE with constant stepsize
// from t = 0.0 to t = 100.0 with dt = 0.01
// and write t and x[0] to cout
integrate_const( rk4 , lorenz , 0.0 , 0.01 , x , 100.0 ,
               cout << _1 << tab << _2[0] << "\n" );

// stepsize controller for rk5 cash-karp
controlled_stepper_standard< stepper_rk5_ck< state_type > >
        rk5_controlled( rk5 , eps_abs, eps_rel, 1.0, 1.0);

// solve the ODE with an adaptive stepsize controller
size_t steps = integrate_adaptive( controlled_rk5, lorenz, x, 0.0,
       10.0 ,  0.01, cout << _1 << tab << _2[0] << "\n" );

Best regards,

Karsten