(*
 * ch4-dfan.m - ducted fan example applied to Chapter 4
 *
 * RMM, 4 Oct 94
 *
 *)

<<Graphics`Graphics`

(* system parameters *)
m = 1.5;	J = 0.0475;
r = 0.25; 	g = 10;
d = 0.1;	l = 0.05;

(* Plant and controller definitions *)
P = r/(J s^2 + m g l);
K = k(s + a) / (s/b + 1);	k = 15;	    a = 25;	b = 300;

L = P K;
S = Together[ 1/(1 + P K) ];
T = Together[ P K / (1 + P K) ];

xyticks = {{1,10,100,1000}, {0.00001, 0.0001, 0.001, 0.01,0.1,1,10}};
bode[G_, opts___] :=
  LogLogPlot[G, {s, 1, 1000},
    Ticks->xyticks, GridLines->xyticks,
    AxesLabel -> {"rad/sec", "dB"},
    opts];
SetOptions[Plot, DisplayFunction->Identity];

(* Nominal performance weight *)
W1 = 20/((s/30) + 1)^2;
Display[
  "perfweight.mps",
  bode[W1, PlotLabel -> "Performance weight"]];
Display[
  "nomperf.mps",
  bode[W1 S, PlotLabel -> "Nominal performance"]];

(* Perturbation #1: 25% change in r *)
W2a = 0.2;
Display[
  "robstab-r.mps",
  bode[W2a T, PlotLabel -> "Robust Stability, 25% change in r"]];
Display[
  "robperf-r.mps",
  bode[Abs[W1 S] + Abs[W2a T],
    PlotLabel -> "Robust Performance, 25% change in r"]];

(* Perturbation #2: 50% uncertainty in l *)
W2b = 0.025 m g / r;
Display[
  "robstab-l.mps",
  bode[W2b P S, PlotLabel->"Robust stability, 50% change in l"]];

(* Perturbation #3: Up to 0.1 sec time lag in motor dynamics *)
W2c = 0.21 s / (0.1 s + 1);
Display[
  "robstab-delay.mps",
  bode[W2c T, PlotLabel->"Robust stability, up to 0.1 sec delay"]];

Display[
  "robperf-delay.mps",
  bode[Abs[W1 S] + Abs[W2c T],
    PlotLabel->"Robust performance, up to 0.1 sec delay"]];
