Finding the modulus margin for three different open loop systems through matlab?

I have three different open loop systems in Matlab with all their Nyquist plots. Someone can easily find the Gm and Pm margins with the margin command, but I can not find anything regarding the modulus margin. Does somebody know how to find this value through Matlab?

Answer 1

[Paul]

Hi Yarno,

I'm not aware of any function that will compute the modulus margin (which is too bad). But one way to get what should be reasonable approximation would be something like this assuming you know something about shape of L and can specify an upper bound of frequency over which to search. There may be a better way to do this.

Define a loop transfer function

L = tf(25,[1 10 10 10]); % loop transfer function

Verify closed loop system is stable.

pole(feedback(L,1))

ans = 
  -9.3303 + 0.0000i
  -0.3349 + 1.9076i
  -0.3349 - 1.9076i

Compute the modulus margin and the frequency at which it occurs and the value of L(s) at that frequency. Could use more sophisticated methods, in principal.

w = linspace(0,1000,1e6);
[modulusmargin,index] = min(abs(freqresp(1 + L,w)));
modulusmargin

modulusmargin = 0.4021

wc = w(index)

wc = 2.0110

Lc = freqresp(L,wc);

The additive perturbation is

-1 - Lc

ans = -0.2888 + 0.2798i

Plot the Nyquist plot and show the modulus margin

figure
nyquist(L)
hold on;
h1 = plot([-1 real(Lc)],[0 imag(Lc)],'-gx');

The concept of Disk Margins might also be of interest. Let's compute the disk margin and add that point to the plot.

DM = diskmargin(L);
Lc = freqresp(L,DM.Frequency);
h2 = plot([-1 real(Lc)],[0 imag(Lc)],'-ro');
legend([h1 h2],{'Modulus Margin' 'Disk Margin'})

Zoom in to see the difference

copyobj(gca,figure)

axis([-1 0 -1 1])

There probably is other functionality in the Robust Control Toolbox to compute the modulus margin. For example, robstab gives a result very close to that computed above. I'm not that familiar with the RCT.

Lp = L + ucomplex('pert',0);
[stabmarg,wcu,info] = robstab(feedback(Lp,1))

stabmarg = struct with fields:
           LowerBound: 0.4015
           UpperBound: 0.4023
    CriticalFrequency: 2.0211

wcu = struct with fields:
    pert: -0.2954 + 0.2730i

info = struct with fields:
                Model: 1
            Frequency: 2.0211
               Bounds: [0.4015 0.4023]
    WorstPerturbation: [1×1 struct]
          Sensitivity: [1×1 struct]

[Yarno Basten]

Thank you so much Paul. This really helped!