ODE response vs CCF response of mass-damper-spring system

Could someone please explain why is that the (zero-input) response obtained using ode45 solver is different than the CCF (controllable canonical form) response as shown below? However, if the C matrix is changed to [1 0], the responses match.

clear all
syms x(t)
x0 = [0;2];
m = 3;
b = 0.6;
k = 1.5;
Fa = 0;
eqn = m*diff(x,2) + b*diff(x,1) + k*x == Fa;
V = odeToVectorField(eqn);
M = matlabFunction(V,'vars',{'t','Y'});
ySol = ode45(M,[0 50],x0);
time = 0:0.1:50;
yValues_pos = deval(ySol,time,1);
s = tf('s');
G = 1/(m*s^2+b*s+k);

% CCF
Accf = [0,      1;
        -k/m,   -b/m];
Bccf = [0;
        1];
Cccf = [1/m, 0];
Dccf = 0;
for i = 1:length(time)
    yccf(i) = Cccf*expm(Accf*time(i))*x0;
end
figure()
hold on
plot(time,yccf,'b')
plot(time,yValues_pos,'r--')
legend('ccf','ODE45')

clear all
syms x(t)
x0 = [0;2];
m = 3;
b = 0.6;
k = 1.5;
Fa = 0;
eqn = m*diff(x,2) + b*diff(x,1) + k*x == Fa;
V = odeToVectorField(eqn);
M = matlabFunction(V,'vars',{'t','Y'});
ySol = ode45(M,[0 50],x0);
time = 0:0.1:50;
yValues_pos = deval(ySol,time,1);
s = tf('s');
G = 1/(m*s^2+b*s+k);

% CCF
Accf = [0,      1;
        -k/m,   -b/m];
Bccf = [0;
        1];
Cccf = [1, 0];
Dccf = 0;
for i = 1:length(time)
    yccf(i) = Cccf*expm(Accf*time(i))*x0;
end
figure()
hold on
plot(time,yccf,'b')
plot(time,yValues_pos,'r--')
legend('ccf','ODE45')

Answer 1

[Sam Chak]

Hi @J AI

Because the output matrix really means for this output equation

where is the output vector and is the state vector (come from your ODE).

Thus, if and , then it implies that the measured output is

.

In fact, the division by m only occurs at the dynamical level.

.

Solving this ODE returns the solution for the and , the states of the system.

So, if you want to compare the ODE result with the result using the state-space model, naturally you want to make for only or for both and .

[J AI]

My question - to be more precise - is that why the CCF has 1/m instead of 1 in the C matrix? Why does CCF scale the output by 1/m?

[Sam Chak]

As I explained previously, the division by m only occurs at the dynamical level; now as shown symbolically below. It is unnecessary to scale the output by , or .

syms x(t) m b k
x0  = [0;2];
% m = 3;
% b = 0.6;
% k = 1.5;
Fa  = 0;
eqn = m*diff(x,2) + b*diff(x,1) + k*x == Fa;
V   = odeToVectorField(eqn)

V =