%% Effect of Disk Shape on Disk Margins
% Varying the eccentricity parameter E amounts to changing the shape of the uncertainty disk \left(1+\delta \right)\;from which the gain and phase margins are derived. Disks of different shapes give more or less conservative estimates of the margins. Thus, you can examine the margins resulting from different E values to get the least conservative estimate for your system.
% To see how the disk margins change with the eccentricity parameter for a particular open-loop transfer function L, vary E between large negative and large positive values. Plot the resulting gain and phase margins.

L = tf(25,[1 10 10 10]);
% Generate vector of E values
E = [linspace(0,1,10) logspace(0,2,10)];
E = unique([E -E]);
% Compute disk margins
for ct=1:numel(E)
    DM(ct) = diskmargin(ss(L),E(ct));
end
% Extract and plot gain and phase margin values against E
GM = cat(1,DM.GainMargin);
PM = cat(1,DM.PhaseMargin);

subplot(211) 
plot(E,GM(:,1),E,GM(:,2))
title('GM'), grid

subplot(212)
plot(E,PM(:,1),E,PM(:,2))
title('PM'), grid

% Consider the gain-margin plot. For this system, large negative values of E give a very conservative estimate for the upper margin, but a much less conservative estimate of the lower margin. The reverse is true for large positive values of E. The true margins of the system must encompass all these estimates. Therefore, you can conclude that the loop remains stable for variations in loop gain within the complete range between the lowest lower margin and the highest upper margin, or factors of about 0.08–1.9.
% For phase margins, the least conservative estimate of about ±25° occurs closer to E = 0.
%
% To see how variation in E values corresponds to different disk shapes, examine the Nyquist plot of the system and the disks corresponding to a few of the E values. You can find the disk for a particular value of E by setting \Delta =\alpha e^{i\theta } and varying \theta between 0 and \pi. The disk margin \alpha is returned by diskmargin in the DiskMargin field of DM.

figure
nyquist(L,'b')
hold on

Evalues = [-3,0,3];
disk = [];

for k = 1:length(Evalues)
    E = Evalues(k);
    DM = diskmargin(L,E);
    theta = 0:0.01:2*pi;
    delta = DM.DiskMargin * exp(1i*theta);
    ExcludeRegion = (delta*(1+E)/2-1) ./ (1+delta*(1-E)/2);
    disk(:,k) = ExcludeRegion';
end

plot(disk)
axis equal
legend('L','E = -3','E = 0','E = 3')

%
% For visibility, the plot shows only three E values, –3, 0, and 3. The plot shows how the shape of the disk, and thus the resulting disk margins, vary as E changes.
%
% Copyright 2018 The MathWorks, Inc.