Search code examples
c++matlabodeodeint

can't get same results Matlab generates


Take a look at the following transfer function:

enter image description here

With Matlab Simulink:

enter image description here

The result is

enter image description here

In State-space representation, the system can be modeled as follows:

enter image description here

In Matlab, we can model the system in the state-space representation:

enter image description here

enter image description here

which yields the following plot:

enter image description here

which is exactly the result generated by using transfer function. I'm trying to generate same results with odeint but failed. This is the code

#include <iostream>
#include <Eigen/Dense>
#include <boost/numeric/odeint.hpp>
#include <iomanip>
#include <fstream>


using namespace std;
using namespace boost::numeric::odeint;

typedef std::vector< double > state_type;

void equations(const state_type &y, state_type &dy, double x)
{
    Eigen::MatrixXd A(3, 3), B(3,1);

    /*
        x = y[0]
       dx = y[1] = dy[0]
      ddx = y[2] = dy[1]
     dddx =        dy[2]

    */
    const double r(1);

    A <<   0,   1,  0,
           0,   0,  1,
         -24, -26, -9;

    B << 0, 0, 1;      
    //#####################( ODE Equations )################################
    Eigen::MatrixXd X(3, 1), dX(3,1);
    X << y[0], y[1], y[2];
    dX = A*X + B*r;

    dy[0] = dX(0,0);
    dy[1] = dX(1,0);
    dy[2] = dX(2,0);
}


int main(int argc, char **argv)
{
    const double dt = 0.01;
    runge_kutta_dopri5 < state_type > stepper;

    state_type y(3); 
    // initial values
    y[0] = 0.0; //  x1 
    y[1] = 0.0; //  x2
    y[2] = 0.0; //  x3

    ofstream data("file.txt");

    for (double t(0.0); t <= 5.0; t += dt){
        data << t << " " << 2*y[0] << " " << 7*y[1] << " " << 1*y[2] << std::endl;
        stepper.do_step(equations, y, t, dt);
    }

    return 0;

}

And this is the result for all state vector

enter image description here

None of the preceding variables match the results generated by Matlab. Any suggestions to fix this code?


Solution

  • Look at the expression you have for y. When you multiply a 1x3 matrix with a 3x1 matrix, the result should be a 1x1 matrix, where the value of the single element is the dot product of the two matrices. What you're currently doing is element-wise multiplication when you write to data instead of calculating the dot product.