Strange behaviour of step response (unstable system)

Hi!

The system defned in the following code should be stable (verified also in simulink). The step response looks correct in the beginning 25sec) but becomes unstable for simulation up to 50sec. Whats happening? Nummerical problems?

s = tf('s');

ks=2;
a=1;
b=2;
 
G=ks/s/(s-a)/(s+b);
 
kp=0.4;
Tv=2*a*b/(b-a)
 
Gr=kp*(1+Tv*s)
 
Gcl=Gr*G/(1+Gr*G);
t=0:0.0001:50;
figure;
step(Gcl,t);

Answer 1

[Paul]

When you do transfer function math like this, matlab doesn't attempt to deal with poles and zeros that should cancel. Then you can run into problems when the pole/zero cancellations that should be exact but aren't because of rounding errors. You can see from your expression for Gcl that Gcl should have a third order denominator. But Gcl as computed is:

>> Gcl

Gcl =
 
     3.2 s^4 + 4 s^3 - 5.6 s^2 - 1.6 s
  ---------------------------------------
  s^6 + 2 s^5 + 0.2 s^4 - 1.6 s^2 - 1.6 s
 
Continuous-time transfer function.

Now look at the poles and zeros:

>> [p,z]=pzmap((Gcl))

p =

   0.0000 + 0.0000i
  -2.0000 + 0.0000i
   1.0000 + 0.0000i
  -0.1107 + 1.0076i
  -0.1107 - 1.0076i
  -0.7785 + 0.0000i


z =

         0
   -2.0000
    1.0000
   -0.2500
   

Note that the first three poles seem to be equivalent to the first three zeros. Mathematically, they should be. But are they numerically:

>> [z(1:3)-p(1:3)].'

ans =

   1.0e-15 *

         0         0   -0.8882
         

So we see that third pole at s=1 doesn't perfectly cancel with the zero. As a result, your step response will eventually go unstable. If you're confident that the those three poles and zeros should cancel based on your problem, you can do:

>> minreal(Gcl)

ans =
 
        3.2 s + 0.8
  -----------------------
  s^3 + s^2 + 1.2 s + 0.8
 
Continuous-time transfer function.

Which is probably the answer you're looking for. Or, you could have gotten that answer directly using:

>> feedback(Gr*G,1)

ans =
 
        3.2 s + 0.8
  -----------------------
  s^3 + s^2 + 1.2 s + 0.8
 
Continuous-time transfer function.

Answer 2

[Benedikt]

Thank you all for this detailed analysis.

Very helpful!