octave_packages/control-2.3.52/MDSSystem.m

   1 %% -*- texinfo -*-
   2 %% Robust control of a mass-damper-spring system.
   3 %% Type @code{which MDSSystem} to locate,
   4 %% @code{edit MDSSystem} to open and simply
   5 %% @code{MDSSystem} to run the example file.
   6
   7 % ===============================================================================
   8 % Robust Control of a Mass-Damper-Spring System     Lukas Reichlin    August 2011
   9 % ===============================================================================
  10 % Reference: Gu, D.W., Petkov, P.Hr. and Konstantinov, M.M.
  11 %            Robust Control Design with Matlab, Springer 2005
  12 % ===============================================================================
  13
  14 % Tabula Rasa
  15 clear all, close all, clc
  16
  17 % ===============================================================================
  18 % System Model
  19 % ===============================================================================
  20 %                +---------------+
  21 %                | d_m   0    0  |
  22 %          +-----|  0   d_c   0  |<----+
  23 %      u_m |     |  0    0   d_k |     | y_m
  24 %      u_c |     +---------------+     | y_c
  25 %      u_k |                           | y_k
  26 %          |     +---------------+     |
  27 %          +---->|               |-----+
  28 %                |     G_nom     |
  29 %        u ----->|               |-----> y
  30 %                +---------------+
  31
  32 % Nominal Values
  33 m_nom = 3;   % mass
  34 c_nom = 1;   % damping coefficient
  35 k_nom = 2;   % spring stiffness
  36
  37 % Perturbations
  38 p_m = 0.4;   % 40% uncertainty in the mass
  39 p_c = 0.2;   % 20% uncertainty in the damping coefficient
  40 p_k = 0.3;   % 30% uncertainty in the spring stiffness
  41
  42 % State-Space Representation
  43 A =   [            0,            1
  44         -k_nom/m_nom, -c_nom/m_nom ];
  45
  46 B1 =  [            0,            0,            0
  47                 -p_m,   -p_c/m_nom,   -p_k/m_nom ];
  48
  49 B2 =  [            0
  50              1/m_nom ];
  51
  52 C1 =  [ -k_nom/m_nom, -c_nom/m_nom
  53                    0,        c_nom
  54                k_nom,            0 ];
  55
  56 C2 =  [            1,            0 ];
  57
  58 D11 = [         -p_m,   -p_c/m_nom,   -p_k/m_nom
  59                    0,            0,            0
  60                    0,            0,            0 ];
  61
  62 D12 = [      1/m_nom
  63                    0
  64                    0 ];
  65
  66 D21 = [            0,            0,            0 ];
  67
  68 D22 = [            0 ];
  69
  70 inname = {'u_m', 'u_c', 'u_k', 'u'};   % input names
  71 outname = {'y_m', 'y_c', 'y_k', 'y'};  % output names
  72
  73 G_nom = ss (A, [B1, B2], [C1; C2], [D11, D12; D21, D22], ...
  74             'inputname', inname, 'outputname', outname);
  75
  76 G = G_nom(4, 4);                       % extract output y and input u
  77
  78
  79 % ===============================================================================
  80 % Frequency Analysis of Uncertain System
  81 % ===============================================================================
  82
  83 % Uncertainties: -1 <= delta_m, delta_c, delta_k <= 1
  84 [delta_m, delta_c, delta_k] = ndgrid ([-1, 0, 1], [-1, 0, 1], [-1, 0, 1]);
  85
  86 % Bode Plots of Perturbed Plants
  87 w = logspace (-1, 1, 100);             % frequency vector
  88 figure (1)
  89
  90 for k = 1 : numel (delta_m)
  91   Delta = diag ([delta_m(k), delta_c(k), delta_k(k)]);
  92   G_per = lft (Delta, G_nom);
  93   bode (G_per, w)
  94   subplot (2, 1, 1)
  95   hold on
  96   subplot (2, 1, 2)
  97   hold on
  98 endfor
  99
 100
 101 % ===============================================================================
 102 % Mixed Sensitivity H-infinity Controller Design (S over KS Method)
 103 % ===============================================================================
 104 %                                    +-------+
 105 %             +--------------------->|  W_p  |----------> e_p
 106 %             |                      +-------+
 107 %             |                      +-------+
 108 %             |                +---->|  W_u  |----------> e_u
 109 %             |                |     +-------+
 110 %             |                |    +---------+
 111 %             |                |  ->|         |->
 112 %  r   +    e |   +-------+  u |    |  G_nom  |
 113 % ----->(+)---+-->|   K   |----+--->|         |----+----> y
 114 %        ^ -      +-------+         +---------+    |
 115 %        |                                         |
 116 %        +-----------------------------------------+
 117
 118 % Weighting Functions
 119 s = tf ('s');                          % transfer function variable
 120 W_p = 0.95 * (s^2 + 1.8*s + 10) / (s^2 + 8.0*s + 0.01);  % performance weighting
 121 W_u = 10^-2;                           % control weighting
 122
 123 % Synthesis
 124 K_mix = mixsyn (G, W_p, W_u);          % mixed-sensitivity H-infinity synthesis
 125
 126 % Interconnections
 127 L_mix = G * K_mix;                     % open loop
 128 T_mix = feedback (L_mix);              % closed loop
 129
 130 % Plotting
 131 figure (2)
 132 bode (K_mix)                           % bode plot
 133
 134 figure (3)
 135 step (T_mix, 10)                       % step response for 10 seconds
 136
 137
 138 % ===============================================================================
 139 % H-infinity Loop-Shaping Design (Normalized Coprime Factor Perturbations)
 140 % ===============================================================================
 141
 142 % Settings
 143 W1 = 8 * (2*s + 1) / (0.9*s);          % precompensator
 144 W2 = 1;                                % postcompensator
 145 factor = 1.1;                          % suboptimal controller
 146
 147 % Synthesis
 148 K_ncf = ncfsyn (G, W1, W2, factor);    % positive feedback controller
 149
 150 % Interconnections
 151 K_ncf = -K_ncf;                        % negative feedback controller
 152 L_ncf = G * K_ncf;                     % open loop
 153 T_ncf = feedback (L_ncf);              % closed loop
 154
 155 % Plotting
 156 figure (4)
 157 bode (K_ncf)                           % bode plot
 158
 159 figure (5)
 160 step (T_ncf, 10)                       % step response for 10 seconds
 161
 162 % ===============================================================================