## How to manually simulate linear system without lsim

I have a 3rd order continuous system of the following form identified with System Identification Toolbox:

y = Cx + Du + e

dx/dt = Ax + Bu + Ke

I can simulate it by using lsim as follows:

simulatedOutput = lsim(systemModel, inputVector, timeVector, initialStateVector)

Now I need to get it out of Matlab and convert it to C++.

As a first step, let’s say I re-create what lsim does within Matlab, since Matlab syntax is easy to convert to C++.

I made an attempt that I thought was pretty reasonable, and I don’t understand why this doesn’t work. The output that was simulated manually by running a script looks nothing like the output from lsim.

Fig. 1: output from lsim (expected)

Fig. 2: output from manual simulation (actual – code below)

The following is my attempt to simulate the identified system by using a Matlab script.

The specific CSV file used here: https://github.com/01binary/systemid/blob/main/input.csv

% A weights (3×3 matrix)

A = [ …

0.0356, -3.4131, -14.9525;

-1.0591, 85.8128, 334.3098;

1.5729, -95.1123, -175.5517;

];

% B weights (3×1 vector)

B = [ …

-0.0019;

0.0403;

-0.0225;

];

% C weights (1×3 vector)

C = [ -5316.9, 24.87, 105.92 ];

% D weight (scalar)

D = 0;

% K weights (3×1 vector)

K = [ …

-0.0025;

-0.0582;

0.0984;

];

% Initial state

x0 = [ …

-0.0458;

0.0099;

-0.0139;

];

% Read input file

csv = readmatrix(‘input.csv’);

% Select column vectors

time = csv(:,1);

measurement = csv(:,2);

input = csv(:,3);

% Simulate

x = x0;

output = zeros(length(time), 1);

for i = 1:length(time)

if i > 1

timeStep = time(i) – time(i-1);

else

timeStep = 0;

end

u = input(i);

y = systemModel(x, u, 0, C, D);

x = systemState(x, u, 0, A, B, K) * timeStep;

output(i) = y;

end

function prediction = systemModel(x, u, e, C, D)

prediction = …

sum(C * x) + … % Add contribution of state

D * u + … % Add contribution of input

e; % Add contribution of disturbance

end

function x1 = systemState(x0, u, e, A, B, K)

x1 = …

A * x0 + … % Add contribution of state

B * u + … % Add contribution of input

e * K; % Add contribution of disturbance

end

I am trying to keep everything really simple and minimal… could someone spot what I am missing?

Thank you!I have a 3rd order continuous system of the following form identified with System Identification Toolbox:

y = Cx + Du + e

dx/dt = Ax + Bu + Ke

I can simulate it by using lsim as follows:

simulatedOutput = lsim(systemModel, inputVector, timeVector, initialStateVector)

Now I need to get it out of Matlab and convert it to C++.

As a first step, let’s say I re-create what lsim does within Matlab, since Matlab syntax is easy to convert to C++.

I made an attempt that I thought was pretty reasonable, and I don’t understand why this doesn’t work. The output that was simulated manually by running a script looks nothing like the output from lsim.

Fig. 1: output from lsim (expected)

Fig. 2: output from manual simulation (actual – code below)

The following is my attempt to simulate the identified system by using a Matlab script.

The specific CSV file used here: https://github.com/01binary/systemid/blob/main/input.csv

% A weights (3×3 matrix)

A = [ …

0.0356, -3.4131, -14.9525;

-1.0591, 85.8128, 334.3098;

1.5729, -95.1123, -175.5517;

];

% B weights (3×1 vector)

B = [ …

-0.0019;

0.0403;

-0.0225;

];

% C weights (1×3 vector)

C = [ -5316.9, 24.87, 105.92 ];

% D weight (scalar)

D = 0;

% K weights (3×1 vector)

K = [ …

-0.0025;

-0.0582;

0.0984;

];

% Initial state

x0 = [ …

-0.0458;

0.0099;

-0.0139;

];

% Read input file

csv = readmatrix(‘input.csv’);

% Select column vectors

time = csv(:,1);

measurement = csv(:,2);

input = csv(:,3);

% Simulate

x = x0;

output = zeros(length(time), 1);

for i = 1:length(time)

if i > 1

timeStep = time(i) – time(i-1);

else

timeStep = 0;

end

u = input(i);

y = systemModel(x, u, 0, C, D);

x = systemState(x, u, 0, A, B, K) * timeStep;

output(i) = y;

end

function prediction = systemModel(x, u, e, C, D)

prediction = …

sum(C * x) + … % Add contribution of state

D * u + … % Add contribution of input

e; % Add contribution of disturbance

end

function x1 = systemState(x0, u, e, A, B, K)

x1 = …

A * x0 + … % Add contribution of state

B * u + … % Add contribution of input

e * K; % Add contribution of disturbance

end

I am trying to keep everything really simple and minimal… could someone spot what I am missing?

Thank you! I have a 3rd order continuous system of the following form identified with System Identification Toolbox:

y = Cx + Du + e

dx/dt = Ax + Bu + Ke

I can simulate it by using lsim as follows:

simulatedOutput = lsim(systemModel, inputVector, timeVector, initialStateVector)

Now I need to get it out of Matlab and convert it to C++.

As a first step, let’s say I re-create what lsim does within Matlab, since Matlab syntax is easy to convert to C++.

I made an attempt that I thought was pretty reasonable, and I don’t understand why this doesn’t work. The output that was simulated manually by running a script looks nothing like the output from lsim.

Fig. 1: output from lsim (expected)

Fig. 2: output from manual simulation (actual – code below)

The following is my attempt to simulate the identified system by using a Matlab script.

The specific CSV file used here: https://github.com/01binary/systemid/blob/main/input.csv

% A weights (3×3 matrix)

A = [ …

0.0356, -3.4131, -14.9525;

-1.0591, 85.8128, 334.3098;

1.5729, -95.1123, -175.5517;

];

% B weights (3×1 vector)

B = [ …

-0.0019;

0.0403;

-0.0225;

];

% C weights (1×3 vector)

C = [ -5316.9, 24.87, 105.92 ];

% D weight (scalar)

D = 0;

% K weights (3×1 vector)

K = [ …

-0.0025;

-0.0582;

0.0984;

];

% Initial state

x0 = [ …

-0.0458;

0.0099;

-0.0139;

];

% Read input file

csv = readmatrix(‘input.csv’);

% Select column vectors

time = csv(:,1);

measurement = csv(:,2);

input = csv(:,3);

% Simulate

x = x0;

output = zeros(length(time), 1);

for i = 1:length(time)

if i > 1

timeStep = time(i) – time(i-1);

else

timeStep = 0;

end

u = input(i);

y = systemModel(x, u, 0, C, D);

x = systemState(x, u, 0, A, B, K) * timeStep;

output(i) = y;

end

function prediction = systemModel(x, u, e, C, D)

prediction = …

sum(C * x) + … % Add contribution of state

D * u + … % Add contribution of input

e; % Add contribution of disturbance

end

function x1 = systemState(x0, u, e, A, B, K)

x1 = …

A * x0 + … % Add contribution of state

B * u + … % Add contribution of input

e * K; % Add contribution of disturbance

end

I am trying to keep everything really simple and minimal… could someone spot what I am missing?

Thank you! system MATLAB Answers — New Questions