Finding Critical Value with Simulink

I have to find several values of k1 in first order inertia and find the critical value when the system is on the stability boundary with k2=2.I have to find them with simulation in order to find the critical value of k1.

In attachemnt you can find my diagram.

[Jon]

You refer to an attachement with your diagram. I don't see any file attached, just a screen shot with part of a simulink diagram.

Also, from at least what I can see of your model, everything is linear, you would be better off doing this analysis using the linear analysis tools, included with the Control System Toolbox. Do you have this toolbox?

Answer 1

[Sam Chak]

Hi @Esin Derin

It's a bit confusing with all the manual switches around. You should also provide the closed-loop system for investigation. Given

With , if my mental calculation is correct, then the closed-loop system should be (please check)

From this point, you can test over a range of and find that should make the system stable.

s = tf('s');
k1 = 0.375;
k2 = 2;
G = k1/(1.5*s^2 + s);
H = k2/(12*s + 1);
Gcl = minreal(feedback(G, H))

Gcl =
 
            0.25 s + 0.02083
  ------------------------------------
  s^3 + 0.75 s^2 + 0.05556 s + 0.04167
 
Continuous-time transfer function.

pole(Gcl)

ans = 
  -0.7500 + 0.0000i
  -0.0000 + 0.2357i
  -0.0000 - 0.2357i

damp(Gcl)       % damping ratio is practically zero, so oscillations are expected.

                                                                       
         Pole              Damping       Frequency      Time Constant  
                                       (rad/seconds)      (seconds)    
                                                                       
 -4.16e-17 + 2.36e-01i     1.77e-16       2.36e-01         2.40e+16    
 -4.16e-17 - 2.36e-01i     1.77e-16       2.36e-01         2.40e+16    
 -7.50e-01                 1.00e+00       7.50e-01         1.33e+00    

step(Gcl, 200)

k1 = 0.00978;               % no overshoot if k1 < 0.00979
G = k1/(1.5*s^2 + s);
Gcl = minreal(feedback(G, H))

Gcl =
 
          0.00652 s + 0.0005433
  -------------------------------------
  s^3 + 0.75 s^2 + 0.05556 s + 0.001087
 
Continuous-time transfer function.

pole(Gcl)

ans = 
  -0.6694 + 0.0000i
  -0.0403 + 0.0008i
  -0.0403 - 0.0008i

damp(Gcl)

                                                                       
         Pole              Damping       Frequency      Time Constant  
                                       (rad/seconds)      (seconds)    
                                                                       
 -4.03e-02 + 7.92e-04i     1.00e+00       4.03e-02         2.48e+01    
 -4.03e-02 - 7.92e-04i     1.00e+00       4.03e-02         2.48e+01    
 -6.69e-01                 1.00e+00       6.69e-01         1.49e+00    

step(Gcl)

Edit: To run this in Simulink, you can probably try Paul's trick as follows:

k1vec   = [0.00978 (0.00978+0.375)/2 0.375]; % choose whatever gain values from 0 < k1 <= 0.375
for ii  = 1:numel(k1vec)
    k1  = k1vec(ii);
    out = sim('Esin'); % use the actual name of the .slx file
    % do whatever processing is needed on out
end