mvpldemo.m

控制系统计算机辅助设计——MATLAB语言与应用(源代码）

% MVPLDEMO    Graphical plots demo for MFD Toolbox.

echo off
%       J-M Boyle  20th August 1987
% Copyright (c) 1987,1993 by GEC Engineering Research Centre and 
% Cambridge Control Ltd
% History:
%       Made matlab 4.0 compatible, 12.7.93, JMM.
%       MRN0016, MRN0017, MRN0023, MRN0026, MRN0036

if exist('genrdata.mat') ~= 2, % If data file does not exist ...
  disp(' Generating data for demonstration')
  genrprep;                % ... then generate it.
end
load genrdata
clf,subplot(1,1,1),hold off % Made V4.0 compatible, 12.7.93, JMM.
%clg,subplot(111),hold off  % Original line
clc
echo on

% Plotting Functions available in the 
% Multivariable Frequency Domain Toolbox.
% =======================================
%
% This demo goes through some of the plotting functions 
% available in the Multivariable Frequency Domain Toolbox.

% The demo uses the Multivariable Frequency Response (MVFR) 
% matrix created during the 'Generating and Displaying a MVFR matrix'
% demo. Remember that the MVFR matrix represents an open loop unstable 
% chemical reactor.


% To refresh your memory, here is that MVFR matrix. Note that it has 
% been calculated at 9 frequency points.

pause % strike any key to continue
clc

more on  % Added 12.7.93, JMM.
fdisp(w, mvfr)

pause % strike any key to continue.
more off % Added 12.7.93, JMM.
clc

% A Bode Plot of the Eigenvalues.
% -------------------------------

% We shall start with Bode Plots. In this example, 
% we'll plot the eigenvalues.

mvfr_eig = feig(w, mvfr);   % Calculate the eigenvalues.

        % Strike any key to display the Bode plot.
        % Strike a further key to continue.
pause

plotbode(w, mvfr_eig)
axes_handles = get(gcf,'Children'); axes(axes_handles(2)); 
title('Bode Plot')
% Note the trickery needed to get the title at the top of the
% Magnitude plot, if you are doing it from a file.
% From the keyboard you can just use `title'.
pause
clc

% That plot didn't look too smooth! Perhaps it is not all that 
% surprising since we were only using 9 frequency points.

% The function FINSERT can be used to insert component matrices, 
% calculated at specified frequencies, into an MVFR matrix.

% It can be a good idea to work with a few frequency points to start 
% with (in order to speed up calculations). Then, when you know where 
% you want more points, you can use the function FINSERT.

pause % Strike any key to continue.
format compact
clc

% Let's use the function FINSERT to add 32 points. We shall add 
% 20 points between 0.01 and 0.316 rad/sec and another 12 points
% between 0.25 and 1.35 rad/sec.

% First, we generate a vector WADD, containing the frequency 
% points at which we wish to insert MVFR matrices.

wadd1 = logspace(-2,-0.5,20); wadd2 = [25:10:135]/100;
wadd = [wadd1,wadd2];

[mvfr, w] = finsert(w, mvfr, wadd, a,b,c,d);

% Here is the resulting frequency vector. Note that the new points 
% have been inserted in the correct positions. We won't show the 
% MVFR matrix as it's as it's rather long now !

pause % strike any key to display the frequency vector.

w
pause % Strike any key to continue.
format
clc

% We can now have a look at the new Bode plot 
% with the new frequency points.

      % Strike any key to display the Bode plot.
      % Strike a further key to continue.
pause

mvfr_eig = feig(w, mvfr);
plotbode(w, mvfr_eig)
axes_handles = get(gcf,'Children'); axes(axes_handles(2)); 
title('Bode Plot {With Extra Points}')
pause

% The Bode plot is now much smoother. However, a discontinuity 
% exists between w = 2 & 3 rad/s. This type of numerical problem 
% can be removed by sorting the eigenvalues into continuous curves 
% using the function CSORT.

% Let's sort the eigenvalues and plot the results.  We can use 
% SUBPLOT to place the magnitudes and the phases on the same screen.

mvfr_eig = csort(mvfr_eig);

pause % Strike any key to start the Bode plot
plotbode(w, mvfr_eig)
axes_handles = get(gcf,'Children'); axes(axes_handles(2)); 
title('Bode Plot {with sorted eigenvalues}')
pause

% Sorting the eigenvalues has produced continuous curves.

pause % Strike any key to continue.
clc, clf % Now V4.0 compatible, 12.7.93, JMM.

% Nyquist Plots.
% --------------

% We can produce Nyquist plots.
% For example, we can plot the individual elements of the MVFR matrix
% all on the same plot.

      % Strike any key to display the Nyquist plot of the MVFR matrix.
      % Strike a further key to clear the plot.
pause

plotnyq(mvfr)
title(' Nyquist Plot of MVFR Matrix [Incorrect] ')
% Note that adding the title is easy with Nyquist plots.
pause
clc

% The plot is a mess because we forgot to use FGET. FGET is used 
% to put the elements of the MVFR matrix we want to plot into 
% individual columns.
% It is recommended that you use FGET in the following way:

% >>  'plot type'( FGET(w, mvfr, index) )
%  where:
%     plot type is the type of plot required,
%     mvfr      is the MVFR matrix whose elements are to be plotted,
%     w         is the associated frequency vector.
%     index     is a two column matrix whose rows specify the elements
%               to be selected. { Default: Select all the elements }.
%
% Example: To pick out and produce Nyquist plot the [1,1] element of
%          the MVFR matrix:
% >> plotnyq( FGET( w, mvfr, [1,1] ))

% Strike any key to display the correct Nyquist plots.
% Strike a further key to clear the plots.

pause
clc,clf
plotnyq( fget(w, mvfr)), title(' Nyquist Plot of MVFR Matrix ')
pause
clc,clf

% Alternatively, we can plot the Nyquist array elements on 
% individual plots using the following commands.

subplot(2,2,1),
plotnyq(fget(w, mvfr, [1,1])), title('[1,1] element'),
subplot(2,2,2),
plotnyq(fget(w, mvfr, [1,2])), title('[1,2] element'),
subplot(2,2,3)
plotnyq(fget(w, mvfr, [2,1])), title('[2,1] element'),
subplot(2,2,4)
plotnyq(fget(w, mvfr, [2,2])), title('[2,2] element')
axis([-3,3,0,1.5])
pause
clc

% Gershgorin Circles.
% -------------------

% Gershgorin circles can be used to establish the region
% within which a system's eigenvalues lie.

% These circles are usually plotted row or column wise on 
% the diagonal elements of the frequency response matrix.

% In this example, we plot the row Gershgorin circles on 
% the diagonal elements. 

pause    % Strike any key to start the plotting process.
hold off

subplot(2,2,1)
plotnyq(fget(w, mvfr, [1,1])), title('[1,1] element')
hold on, plotnyq(frgersh(w, mvfr, 1),'-w'); hold off

subplot(2,2,4)
axis([-3,3,0,1.5])
hold on, plotnyq(frgersh(w, mvfr, 2),'-w'); 
pause
hold off

clf,clc

% Nichols Charts.
% ---------------

% We have included a function for plotting Nichols charts.

      % Strike any key to display the Nichols chart of all
      % the MVFR elements simultaneously.
      % Strike a further key to clear the plot.

pause

plotnic(fget(w, mvfr)), title('Nichols chart for MVFR matrix')
pause
clc, clf
% M and N circles.
% ----------------
% We shall now show how to plot M and N-circles. Such circles 
% are usually plotted on a Nichols chart or a Nyquist Diagram.
% In this demo, we shall draw them on a Nichols chart.

% M and N-circles are used to indicate the closed loop gain and phase,
% respectively, of a system from the open loop frequency response plot.

% The function MCIRC(m) calculates M-circles for a given vector 'm' of gain
% values in dB.

% The function NCIRC(n) calculates N-circles for a given vector 'n' of
% negative phases in degrees.

% To make it easier to see the M and N-circles, we shall only 
% plot the Nichols chart of the [2,1] element of the MVFR matrix. 
% We shall then place M and N circles on the resulting chart.

pause % Strike any key to continue.
clc


% We shall plot M circles for 0, 1 and 3 dB.
% and N-circles for -140, -120, and -90 degrees

m = [0,1,3]; mc = mcirc(m);
n = [-140, -120, -90]; nc = ncirc(n);

% The MVFR element will be white, the M-circles red, and
% the N-circles green.

pause % Strike any key to continue.

plotnic(fget(w, mvfr, [2,1]),'-w'), hold on
plotnic(mc,'--r')           % The M-circles
plotnic(nc,'--g')           % The N-circles
title('Nichols plot with M and N-circles')
pause
clc

% Frequency points can be marked on Nichols (or Nyquist) plots 
% using the MKNIC (or MKNYQ) function.

      % Strike any key to display Nichols chart with frequencies marked.
      % Strike a further key to clear the Nichols plot.
pause

mknic(fget(w ,mvfr, [2,1]), w,1)

% If the values marked using MKNYQ overlap on your screen, you can 
% modify the parameters at the beginning of the MKPLOT.M file.

pause % Strike any key to continue.
clc
more on
% Multivariable Plotting Functions
% -------------------------------

% To simplify the process of plotting and marking an
% element of an MVFR matrix we provide MultiVariable
% versions of the plotting functions:
%           MVBODE, MVNYQ, MVNIC
% These functions perform an FGET and then plot the elements,
% and in the case of MVNYQ and MVNIC also mark the plot with
% the frequency points
%
% For example to reproduce a plot similar to the one just shown
% you would type :
%
% mvnic(w,mvfr,[2,1])
% plotnic([nc,mc]), title('Nichols plot with M and N-circles')
%
% The slight difference would be in the intervals of the marked
% frequency points. MVNYQ and MVNIC use the default setting for
% MKNYQ and MKNIC of 4 plots per screen whereas the previous plot used
% the 1 plot per screen setting and so showed more frequency points.

pause % Strike any key to return to the Demonstrations Menu.

echo off
more off
subplot(1,1,1),clf,hold off