LFRT extension for the Skew Mu Toolbox (SMT)

• Introduction
• List of functions
• µ-analysis with frequency sweeping
• Skew µ-analysis
• Comments
• Installation

Introduction

Skew µ-analysis is an extension of classical µ-analysis in which uncertainties are split into two subsets. Uncertainties of the first subset are treated as for classical µ-analysis, i.e., we look for the worst case with minimum deviation from the origin. Uncertainties of the second subset are just bounded, usually in the unit ball. For more details, see report 1 and report 2.

The Skew Mu Toolbox written by G. Ferreres and J.M. Biannic offers several tools for skew µ-analysis but also for µ-analysis with frequency sweeping and in some cases both. The objective of this extension of the LFR Toolbox is to propose simple interfacing functions with predefined default options. For example, in most case it suffices to invoke

• [...] = skewmu(sys,sknames) or
• [...] = usweepmu(sys,sknames) or
• [...] = lsweepmu(sys)

It is also possible to call directly the functions of the Skew Mu Toolbox without the above interfacing functions. For this purpose the function lfr2smt retrieves form a dynamic LFR-object the input arguments to be used with SMT functions:

• [sys_M,blk] = lfr2smt(sys,sknames);

List of functions.

• skewmu: skew-µ-analysis on a frequency gridding.
• usweepmu: upper bound of µ based on continuous frequency sweeping. Skew-µ-analysis in some cases.
• lsweepmu: lower bound of µ based on continuous frequency sweeping.
• lfr2smt: converts an LFR-objetc to a system in SMT format.
• skewblk: sub-function of the above function, converts the block structure of an LFR-objetc to the SMT format.

µ-analysis with frequency sweeping.

For illustrating the interest of using µ-analysis based on continuous frequency sweeping, a system with numerous flexible modes will be considered. This system has 60 states and a 30x30 Δ-matrix. This is a very challenging problem, classical µ-analysis would lead to useless results.

For loading the LFR-object corresponding to this system, type load demoskmu (the uncertainties are already normalized otherwise you should apply the function normalizelfr before µ-analysis).

The upper-bound analysis is performed by invoking the function usweepmu and the lower-bound (it is more precisely a lower bound of the peak values of µ) by invoking the function lsweepmu.

% Load system sys (load demoskmu_v6 for Matlab 6)
   load demoskmu
   size(sys)

   LFR-object with 0 output(s), 0 input(s) and 60 state(s).
   Uncertainty blocks (globally (30 x 30)):
    Name  Dims  Type   Real/Cplx   Full/Scal   Bounds
    x1    10x10 LTI       r           s        [-1,1]
    x2    10x10 LTI       r           s        [-1,1]
    x3    8x8   LTI       r           s        [-1,1]
    x4    2x2   LTI       r           s        [-1,1]

% Upper bound
   [u_mumax,u_muval,u_pulsint] = usweepmu(sys);
   plot_bar(u_pulsint,u_muval,1,'g-');
   hold on

% Lower bound
   [l_mumax,l_pulsmax,l_paramax] = lsweepmu(sys);
   plot(l_pulsmax(1),l_mumax(1),'r*');

The figures below show the results: the green line is the upper bound and the red dot is the lower bound. The second figure is a zoom of the first one, it shows that classical µ-analysis based on frequency gridding has no chance to identify the peak.

It is now possible to check that the worst case corresponding to the parameters values l_paramax (obtained from the lower-bound of µ) corresponds to a system at the limit of stability (below, the function uplft evaluates the LFR-object sys at xi's equal to l_paramax).

   sys_worst_case = uplft(sys,{'x1','x2','x3','x4'},l_paramax{1});
   eig(sys_worst_case)

   ans =

  -0.1774 +30.0199i
  -0.1774 -30.0199i
  -1.0057 +29.0008i
  -1.0057 -29.0008i
  -2.4737 +27.9070i
  -2.4737 -27.9070i
  -0.9305 +26.9988i
  -0.9305 -26.9988i
  -1.9288 +25.9416i
  -1.9288 -25.9416i
  -1.6400 +24.9580i
  -1.6400 -24.9580i
  -2.2912 +23.9009i
  -2.2912 -23.9009i
  -1.6750 +22.9483i
  -1.6750 -22.9483i
  -1.8648 +21.9291i
  -1.8648 -21.9291i
  -1.5264 +20.9517i
  -1.5264 -20.9517i
  -0.0446 +20.0064i
  -0.0446 -20.0064i
  -1.6534 +18.9335i
  -1.6534 -18.9335i
  -1.4699 +17.9447i
  -1.4699 -17.9447i
  -0.6965 +16.9899i
  -0.6965 -16.9899i
  -0.1498 +16.0029i
  -0.1498 -16.0029i
  -0.0664 +15.0029i
  -0.0664 -15.0029i
  -0.3848 +13.9972i
  -0.3848 -13.9972i
  -0.4360 +12.9948i
  -0.4360 -12.9948i
   0.0000 +12.0017i   %<- pole crossing the imaginary axis
   0.0000 -12.0017i   %<- at the frequency l_pulsmax
  -0.7018 +10.9789i
  -0.7018 -10.9789i
  -0.1609 + 9.9998i
  -0.1609 - 9.9998i
  -0.1867 + 8.9989i
  -0.1867 - 8.9989i
  -0.0488 + 0.9988i
  -0.0488 - 0.9988i
  -0.7711 + 7.9633i
  -0.7711 - 7.9633i
  -0.0476 + 7.0003i
  -0.0476 - 7.0003i
  -0.0182 + 1.9999i
  -0.0182 - 1.9999i
  -0.2021 + 2.9932i
  -0.2021 - 2.9932i
  -0.1107 + 4.9990i
  -0.1107 - 4.9990i
  -0.2059 + 3.9948i
  -0.2059 - 3.9948i
  -0.4349 + 5.9845i
  -0.4349 - 5.9845i

Skew µ-analysis.

In this section is illustrated the function skewmu which is based on frequency gridding. For classical µ-problems (sknames = []), it is the standard function of the Mu or LMI toolboxes (or of the Robust Control Toolbox version 3 if you use Matlab 7) which are invoked. In the skew-µ case, specific algorithms of the SMT toolbox are used (the upper bound algorithm is slower (LMI-based) and the lower bound algorithm (power algorithm) is less efficient).

For illustrating the use of the function skewmu we shall consider an LFR-object with complex uncertainties for better convergence of the lower bounds (if all uncertainties are real use preferably lsweepmu but there are some restrictions).

% Generation of an academic example
   lfrs x1 x2 [-2 -2] [2 3] [0.1 0.2]
   lfrs c1 c2 complex
   a = [-1+x1*x2 -20*(1+c2) c1*x2 c1*c2;...
        20+c1*c2 -1 x2*c1 x2*c2;...
        x1*x2*c1 5*c1 -0.1+c2*c1 -2;...
        1+c1 3+c2 2+x2^3 -0.1+c2*x1];
   b = ones(4,1); c = ones(1,4); d=0;
   sys = eval(lfr(a,b,c,d,'c'));
   size(sys)

   LFR-object with 1 output(s), 1 input(s) and 4 state(s).
   Uncertainty blocks (globally (26 x 26)):
    Name  Dims  Type   Real/Cplx   Full/Scal   Bounds
    c1    8x8   LTI       c           s        |c1| < 1
    c2    7x7   LTI       c           s        |c2| < 1
    x1    3x3   LTI       r           s        [-2,2], nominal=0.1
    x2    8x8   LTI       r           s        [-2,3], nominal=0.2

% before using µ-analysis, the parameters must be normalized
   sysn = normalizelfr(sys);

% Order reduction is recommended
   sysn = minlfr(sysn);
   size(sysn)

   LFR-object with 1 output(s), 1 input(s) and 4 state(s).
   Uncertainty blocks (globally (19 x 19)):
    Name  Dims  Type   Real/Cplx   Full/Scal   Bounds
    c1    5x5   LTI       c           s        |c1| < 1
    c2    5x5   LTI       c           s        |c2| < 1
    x1    3x3   LTI       r           s        [-1,1]
    x2    6x6   LTI       r           s        [-1,1]

First, is performed a standard µ-analysis. In this case, it is a function of the Mu (or LMI with non-default options) standard Matlab toolbox which is used.

Then, a skew-µ analysis is performed considering that x1 can take any value in the unit ball. For the skew-µ upperbound, an heuristic algorithm is used (it doesn't converge at 2 or 3 frequencies for the considered example but with option 'u0' it is possible to switch to an LMI-based skew-µ upperbound but CPU time increases considerably).

% Standard mu
   [bounds1,params1,deltas1,puls1] = skewmu(sysn);
   semilogx(puls1,bounds1(:,1),'r-',puls1,bounds1(:,2),'m.');

% Skew-mu
   [bounds2,params2,deltas2,puls2] = skewmu(sysn,{'x1'});
   hold on
   semilogx(puls2,bounds2(:,1),'b-',puls2,bounds2(:,2),'c.');

The next figure shows the results: the red line corresponds to the standard upper bound, the magenta dotted line to the standard lower bound, the blue line to the skew upper bound and the cyan dotted line to the skew lower bound. As expected, the skew-µ values are larger that the µ values (robust stability with -1 < x1 < 1).

% Worst case: normalized system
   names = {'c1','c2','x1','x2'};
   worstparam = params2{45}

   worstparam =

      0.0626+0.0138i  -0.0032-0.0640i  1.0000  0.0641 % <- x1 = 1

   eig(uplft(sysn,names,worstparam))

   ans =

     -1.2819 -19.9756i
     -0.1028 - 2.0908i
     -0.1029 + 1.9594i
      0.0000 +19.9750i % <- worst case pole crossing
                       %    the imginary axis

% Worst case: original system
   actualparam = actualval(sysn,names,worstparam)
   eig(uplft(sys,names,actualparam))

   ans =

     -0.0000 +19.9750i
     -1.2819 -19.9756i
     -0.1029 + 1.9594i
     -0.1028 - 2.0908i

Miscellaneous comments.

• Never use µ-analysis in case of a single real uncertain parameter !
• The "skew" and "sweep" approaches were introduced separately but the function usweepmu can be used for both simultaneously if invoked with options 'u0' or 'u2'.
• In this tutorial, robust LHP-stability was considered, however most functions of the Skew Mu Toolbox allow us to consider robust D-stability.
  • With function skewmu you can define the input argument points as a set of complex number delimiting any kind of domain.
  • With function usweepmu it is possible to define a truncated sector as a stability domain (input argument domain).
• Two functions of the SMT with very limited field of application are not interfaced with the LFR toolbox, these functions are mu_zd (upper bound, non-repeated real uncertainties) and mu_dailey (lower bound, non-repeated real uncertainties).

Installation

• You should have Matlab version 6 or higher and the Mu and LMI toolboxes.
• This extension of the LFR Toolbox LFRT V2 is an interface with the Skew Mu Toolbox, so, both toolboxes should be installed.
• Download and unzip this file.
• This will generate a directory smt_interface which should be placed in the main directory LFRTv20 (for example) of the LFR Toolbox.
• As a final point, you must update your Matlab path as follows :

  path('.../toolbox/LFRTv20/smt_interface',path);