Matlab step() from Simulink dataset

I want to execute the step method in Matlab using data from Simulink. Step() requires tf(numerator, denominator). How do I get them from simulink?

k=1;
Tk=0.2;
s=0.15;
W0=tf(k,[Tk^2 2*s*Tk 1]);%ПФ объекта

Woy_ss=ss(Woy); %Объект в формате пространства состояний
p=eig(Woy_ss.A); %Полюса ПФ объекта

%Матрица управляемости
My=ctrb(Woy_ss.A,Woy_ss.B);
rgy=rank(My);

%Матрица наблюдаемости
Mn=obsv(Woy_ss.A,Woy_ss.C);
rgn=rank(Mn);

Q=[1 0
   0 1];
R=[1];
[K,P,E] = lqr(Woy_ss.A,Woy_ss.B,Q,R);
%К-параметры регулятора
%Р-решение алгебраического уравнения Риккати
%Е-корни характеристического полинома замкнутой системы

eig(Woy_ss.A-Woy_ss.B*K); %Проверка корней Характеристического 
% Полинома замкнутой системы

sim('SIM_lab5.slx');
step(ans.y)

Answer 1

[Sulaymon Eshkabilov]

Hi Oleg,

One can obtain tf or state space representation of the system model from a Simulink model by linearizing it the Simulink model - here is the DEMO that shows how to import a linearized ss() or tf() into MATLAB workspace.

An alternative way is manual derivation of a transfer fucntionf rom your Simulink model, and then use it for step response or for some other distrubances with lsim() - see DOC.

Answer 2

[Paul]

Based on the picture of the Simulink model ....

sys = linearize('SIM_lab5');
step(sys)

Answer 3

[Sam Chak]

Hi @Oleg

Some information is unavailable in your code, but I managed to guess and reproduce the result of the computed LQR gain, K = [0.1075, 0.1561], as shown in the Workspace in the image. It seems that the gain matrix K is obtained by computing the LQR using the A and B matrices of the state-space model (sys) that are converted from the W0 transfer function model.

For some unknown reason, the gain matrix K is applied to another transformed state-space model, as shown in your Simulink block diagram. This is likely to result in an undesired step response. If you wish to employ the LQR in your case, the converted state-space model shall be first transformed to the state-space in controllable canonical form (csys), followed by computing the LQR using the A and B matrices of the controllable canonical state-space model.

Despite csys being technically equivalent to sys in terms of output response, the computed LQR gain matrix from sys is different from the LQR gain matrix from csys. If the control system is designed correctly, both results in MATLAB and Simulink should be identical.

%% Parameters of the original transfer function
k       = 1;
Tk      = 0.2;
sig     = 1.5;  % original value is 0.15 (Oleg)
W0      = tf(k, [Tk^2 2*sig*Tk 1])    % ПФ объекта

W0 =
 
           1
  --------------------
  0.04 s^2 + 0.6 s + 1
 
Continuous-time transfer function.

sys     = ss(W0)            % convert W0 tf model to state-space model form

sys =
 
  A = 
          x1     x2
   x1    -15  -6.25
   x2      4      0
 
  B = 
       u1
   x1   2
   x2   0
 
  C = 
          x1     x2
   y1      0  3.125
 
  D = 
       u1
   y1   0
 
Continuous-time state-space model.

p1      = eig(sys)          % check if p1 = p2

p1 = 2×1

-13.0902 -1.9098

%% LQR design for sys
Q       = eye(2);
R       = 1;
[K1,P1,E1]= lqr(sys.A, sys.B, Q, R) % same as K in the image

K1 = 1×2

0.1075 0.1561

P1 = 2×2

0.0538 0.0781 0.0781 0.3809

E1 = 2×1

-1.9839 -13.2312

%% Transform sys to controllable canonical form
Gp      = minreal(W0)       % transfer function equivalent to W0

Gp =
 
        25
  ---------------
  s^2 + 15 s + 25
 
Continuous-time transfer function.

[nm,dn] = tfdata(Gp, 'v');
A       = [0 dn(1); -dn(3) -dn(2)];
B       = [0; nm(3)];
C       = [1 0];
D       = 0*C*B;
csys    = ss(A, B, C, D)

csys =
 
  A = 
        x1   x2
   x1    0    1
   x2  -25  -15
 
  B = 
       u1
   x1   0
   x2  25
 
  C = 
       x1  x2
   y1   1   0
 
  D = 
       u1
   y1   0
 
Continuous-time state-space model.

p2      = eig(csys)         % since p1 = p2, therefore sys = csys

p2 = 2×1

-1.9098 -13.0902

%% LQR design for csys
[K2,P2,E2] = lqr(A, B, Q, R)

K2 = 1×2

0.4142 0.5803

P2 = 2×2

1.0692 0.0166 0.0166 0.0232

E2 = 2×1

-1.2512 -28.2566

%% closed-loop control system
clsys       = ss(A-B*K2, B*(K2(1)+1), C, D);
step(clsys, 10), grid on