cannot use `lsim()` function properly

This is my Simulink model

I get into trouble when I use MATLAB code to achieve the same response as Simulink

Here is my .m file

M = 1.5;
m = 0.5;
l = 10;
g = 9.8;

%% state space
A = [0 1 0 0; 0 0 m*g/M 0; 0 0 0 1; 0 0 g/(M*l)*(M+m) 0];
B = [0; 1/M; 0; 1/(M*l)];
C = [0 0 1 0];
D = 0;
sys_ss = ss(A, B, C, D);

% pid
Kp = 99;
Ki = 49;
Kd = 29;
Gc = pid(Kp, Ki, Kd);
% feedback closed-up system
sys_fb = feedback(Gc*sys_ss, 1);
%%
% time
t = 0: 0.01: 10;
r = zeros(size(t));

% start
[y, t] = lsim(sys_fb, r, t, [2; 10; -5; 0]);

Error using DynamicSystem/lsim
Cannot simulate state trajectory for models with singular E matrix. Use the "isproper" command to check if the model is proper and to extract an equivalent explicit representation.

% plot
figure();
plot(t, y);
title('PID控制器 + 状态空间模型闭环系统输出');
xlabel('时间 (s)');
ylabel('输出');

you can directly run this code

when I use lsim()，I meet an error “Cannot simulate response with initial condition for models with singular E matrix. Use the "isproper" command to check if the model is proper and to extract an equivalent explicit representation.”

How does it happen? Please help you, thanks

Answer 1

[Paul]

Hi 嘉伟,

Here is the plant model

M = 1.5;
m = 0.5;
l = 10;
g = 9.8;
%% state space
A = [0 1 0 0; 0 0 m*g/M 0; 0 0 0 1; 0 0 g/(M*l)*(M+m) 0];
B = [0; 1/M; 0; 1/(M*l)];
C = [0 0 1 0];
D = 0;
sys_ss = ss(A, B, C, D);

Here is the PID controller.

% pid
Kp = 99;
Ki = 49;
Kd = 29;
Gc = pid(Kp, Ki, Kd)

Gc =
 
             1          
  Kp + Ki * --- + Kd * s
             s          

  with Kp = 99, Ki = 49, Kd = 29
 
Continuous-time PID controller in parallel form.

We see that its transfer function is improper, i.e., the order of the numerator is higher than the order of the denominator

tf(Gc)

ans =
 
  29 s^2 + 99 s + 49
  ------------------
          s
 
Continuous-time transfer function.

When converted to state space, it's in descriptor form: E*xdot = A*x + B*u; y = C*x + D*u

ss(Gc)

ans =
 
  A = 
       x1  x2  x3
   x1   0   0   0
   x2   0   1   0
   x3   0   0   1
 
  B = 
       u1
   x1   1
   x2   0
   x3  -1
 
  C = 
       x1  x2  x3
   y1  49  29   0
 
  D = 
       u1
   y1  99
 
  E = 
       x1  x2  x3
   x1   1   0   0
   x2   0   0   1
   x3   0   0   0
 
Continuous-time state-space model.

For these inputs, feedback returns a state space model that is also in descriptor form, i.e., has non-zero E matrix

% feedback closed-up system
sys_fb = feedback(Gc*sys_ss, 1);
sys_fb.E

ans = 7×7

1 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 1

That's what cause the error you see.

However, that form of Gc is not what's used in the Simulink model, assuming the PID block is using the default parameters. From the PID Block doc page, we see that the derivative term in the PID control isn't a pure derivative. Rather, it has the form

Kd*(N*s/(s + N), with N = 100 as the default value

The same can be adhieved in Matlab using the fourth argument to pid

Gc = pid(Kp, Ki, Kd, 0.01)

Gc =
 
             1            s    
  Kp + Ki * --- + Kd * --------
             s          Tf*s+1 

  with Kp = 99, Ki = 49, Kd = 29, Tf = 0.01
 
Continuous-time PIDF controller in parallel form.

Now, Gc is proper

tf(Gc)

ans =
 
  2999 s^2 + 9949 s + 4900
  ------------------------
        s^2 + 100 s
 
Continuous-time transfer function.

And the closed loop is not in descriptor form

% feedback closed-up system
sys_fb = feedback(Gc*sys_ss, 1);
sys_fb.E

ans =

     []

But now there is another problem. sys_fb has six states

size(sys_fb)

State-space model with 1 outputs, 1 inputs, and 6 states.

But the call to lsim only specifies four initial conditions, which will result in an error

%%
% time
t = 0: 0.01: 10;
r = zeros(size(t));
% start
try
[y, t] = lsim(sys_fb, r, t, [2; 10; -5; 0]);
catch
    disp('lsim error')
end

lsim error

There is also another complication. The first input to the feedback command was Gc*sys_ss. But if you look at the doc page for series, that's not really what we want, which is sys_ss*Gc. Those two forms are the same from an input/output perspective, but are different with respect to the relationship between the output and the states.

So, assuming that those initial conditions apply to the states of the plant, and making the input to the feedback command sys_ss*Gc, and assuming that the initial conditions in the PID block in Simulink use the default values of zero, I think what we want is

sys_fb = feedback(sys_ss*Gc, 1);
[y, t] = lsim(sys_fb, r, t, [2; 10; -5; 0; 0; 0;]);

figure
plot(t,y),grid

[嘉伟]

Hello @Paul, I sincerely appreciate your response.

As a Matlab beginner, I have spent several hours on this issue and asked for help on various forums but received no replies.

Just as I was about to sleep (perhaps at a different time than yours due to the time difference between Beijing, China and your location), I saw your answer which was very patient and detailed. I really thank you for your help.Wishing you all the best and happiness!