Search code examples
pythonc++scipyodeint

Solution for differential equation is completly different in boost::odeint and scipy.integrate


I am trying to port my rapid prototyping from python to C++. I try to test the notation with a simple differential equation, but the results are very different for starting value [2,0]. Python is declining, while the C++ solution is rising strongly.

It worked for a sample found here: How to incorporate time-varying parameters from lookup table into boost::odeint, c++

but it does not work for my example

TransferF::TransferF(const double& deltaT) : dt(deltaT), t(0.0), y(2)
{
    // initial values
    y[0] = 2.0;  //  x1 
    y[1] = 0.0;  //  x2
}


void TransferF::ode(const state_type &x, state_type &y, double t)
{
    y[0] = x[0];
    y[1] = x[1];
    y[2] = (-2*y[1] - y[0] + 1) / (pow(y[0],2));
}

and the same in py:

def modelt(x,t):
    y = x[0]
    dydt = x[1]
    dy2dt2 = (-2*dydt - y + 1)/ (y **2)
    return [dydt,dy2dt2]

x3 = odeint(modelt,[2,0],timev)

I expected the same results for the time series, but pythons solution is falling, C++ is rising.


Solution

  • The C++ code has a subtle inconsistency. The output vector y should only contain the derivatives y', y", not the function y itself:

    void TransferF::ode(const state_type &x, state_type &y, double t)
    {
        y[0] = x[1];
        y[1] = (-2*x[1] - x[0] + 1) / (pow(x[0],2));
    }