Robust Matlab

Robust Control Toolbox™
Getting Started Guide
R2013b
Gary Balas
Richard Chiang
Andy Packard
Michael Safonov
How to Contact MathWorks
www.mathworks.com Web
comp.soft-sys.matlab Newsgroup
www.mathworks.com/contact_TS.html Technical Support
[email protected] Product enhancement suggestions

[email protected] Bug reports
[email protected] Documentation error reports
[email protected] Order status, license renewals, passcodes
[email protected] Sales, pricing, and general information
508-647-7000 (Phone)
508-647-7001 (Fax)
The MathWorks, Inc.

3 Apple Hill Drive
Natick, MA 01760-2098
For contact information about worldwide offices, see the MathWorks Web site.
Robust Control Toolbox™ Getting Started Guide
© COPYRIGHT 2005–2013 by The MathWorks, Inc.
The software described in this document is furnished under a license agreement. The software may be used
or copied only under the terms of the license agreement. No part of this manual may be photocopied or
reproduced in any form without prior written consent from The MathWorks, Inc.
FEDERAL ACQUISITION: This provision applies to all acquisitions of the Program and Documentation
by, for, or through the federal government of the United States. By accepting delivery of the Program
or Documentation, the government hereby agrees that this software or documentation qualifies as
commercial computer software or commercial computer software documentation as such terms are used
or defined in FAR 12.212, DFARS Part 227.72, and DFARS 252.227-7014. Accordingly, the terms and
conditions of this Agreement and only those rights specified in this Agreement, shall pertain to and govern
the use, modification, reproduction, release, performance, display, and disclosure of the Program and
Documentation by the federal government (or other entity acquiring for or through the federal government)
and shall supersede any conflicting contractual terms or conditions. If this License fails to meet the
government’s needs or is inconsistent in any respect with federal procurement law, the government agrees
to return the Program and Documentation, unused, to The MathWorks, Inc.
Trademarks
MATLAB and Simulink are registered trademarks of The MathWorks, Inc. See
www.mathworks.com/trademarks for a list of additional trademarks. Other product or brand
names may be trademarks or registered trademarks of their respective holders.
Patents
MathWorks products are protected by one or more U.S. patents. Please see
www.mathworks.com/patents for more information.
Revision History
September 2005 First printing New for Version 3.0.2 (Release 14SP3)
March 2006 Online only Revised for Version 3.1 (Release 2006a)
September 2006 Online only Revised for Version 3.1.1 (Release 2006b)
September 2007 Online only Revised for Version 3.3 (Release 2007b)
March 2008 Online only Revised for Version 3.3.1 (Release 2008a)
October 2008 Online only Revised for Version 3.3.2 (Release 2008b)
April 2011 Online only Revised for Version 3.6 (Release 2011a)
Contents
Introduction
1
Robust Control Toolbox Product Description . . . . . . . . 1-2
Key Features . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1-2
Product Requirements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1-3
Modeling Uncertainty . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1-4
System with Uncertain Parameters ................. 1-5
Worst-Case Performance . . . . . . . . . . . . . . . . . . . . . . . . . . . 1-9
Worst-Case Performance of Uncertain System . . . . . . . . 1-10
Loop-Shaping Controller Design . . . . . . . . . . . . . . . . . . . . 1-12
Model Reduction and Approximation . . . . . . . . . . . . . . . . 1-16
Reduce Order of Synthesized Controller . . . . . . . . . . . . . 1-17
LMI Solvers . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1-21
Extends Control System Toolbox Capabilities . . . . . . . . 1-22
Acknowledgments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1-23
Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1-25
v
Multivariable Loop Shaping
2
Tradeoff Between Performance and Robustness . . . . . . 2-2
Norms and Singular Values . . . . . . . . . . . . . . . . . . . . . . . . . 2-4

Properties of Singular Values . . . . . . . . . . . . . . . . . . . . . . . . 2-4
Typical Loop Shapes, S and T Design . . . . . . . . . . . . . . . . 2-6

Robustness in Terms of Singular Values . . . . . . . . . . . . . . . 2-7
Guaranteed Gain/Phase Margins in MIMO Systems . . . . . 2-12
Using LOOPSYN to Do H-Infinity Loop Shaping . . . . . . 2-15
Loop-Shaping Control Design of Aircraft Model . . . . . . 2-16

Design Specifications . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2-18
MATLAB Commands for a LOOPSYN Design . . . . . . . . . . 2-18
Fine-Tuning the LOOPSYN Target Loop Shape Gd to

Meet Design Goals . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2-21
Mixed-Sensitivity Loop Shaping . . . . . . . . . . . . . . . . . . . . 2-22
Mixed-Sensitivity Loop-Shaping Controller Design . . . 2-24
Loop-Shaping Commands . . . . . . . . . . . . . . . . . . . . . . . . . . 2-26
Model Reduction for Robust Control

3
Why Reduce Model Order? ......................... 3-2
Hankel Singular Values . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3-3
Model Reduction Techniques . . . . . . . . . . . . . . . . . . . . . . . 3-5
vi Contents
Approximate Plant Model by Additive Error
Methods . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3-7
Approximate Plant Model by Multiplicative Error

Method . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3-9
Using Modal Algorithms . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3-12

Rigid Body Dynamics . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3-12
Reducing Large-Scale Models . . . . . . . . . . . . . . . . . . . . . . . 3-16
Normalized Coprime Factor Reduction . . . . . . . . . . . . . . 3-17
Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3-19
Robustness Analysis
4
Create Models of Uncertain Systems . . . . . . . . . . . . . . . . 4-2
Creating Uncertain Models of Dynamic Systems . . . . . . . . 4-2
Creating Uncertain Parameters . . . . . . . . . . . . . . . . . . . . . . 4-3
Quantifying Unmodeled Dynamics . . . . . . . . . . . . . . . . . . . 4-6
Robust Controller Design .......................... 4-9
MIMO Robustness Analysis . . . . . . . . . . . . . . . . . . . . . . . . . 4-14

Adding Independent Input Uncertainty to Each
Channel . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 4-15
Closed-Loop Robustness Analysis . . . . . . . . . . . . . . . . . . . . 4-19
Nominal Stability Margins . . . . . . . . . . . . . . . . . . . . . . . . . . 4-21
Robustness of Stability Model Uncertainty . . . . . . . . . . . . . 4-22
Worst-Case Gain Analysis . . . . . . . . . . . . . . . . . . . . . . . . . . 4-23
Summary of Robustness Analysis Tools . . . . . . . . . . . . . . 4-27
vii
H-Infinity and Mu Synthesis
5
Interpretation of H-Infinity Norm . . . . . . . . . . . . . . . . . . . 5-2
Norms of Signals and Systems . . . . . . . . . . . . . . . . . . . . . . . 5-2
Using Weighted Norms to Characterize Performance . . . . 5-3
H-Infinity Performance . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5-9

Performance as Generalized Disturbance Rejection . . . . . . 5-9
Robustness in the H-Infinity Framework . . . . . . . . . . . . . . 5-15
Active Suspension Control Design . . . . . . . . . . . . . . . . . . . 5-17
Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5-38
Tuning Fixed Control Architectures

6
What Is a Fixed-Structure Control System? . . . . . . . . . . 6-3
Choosing an Approach to Fixed Control Structure

Tuning . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-4
Difference Between Fixed-Structure Tuning and

Traditional H-Infinity Synthesis . . . . . . . . . . . . . . . . . . 6-5
Bibliography . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-5
Structure of Control System for Tuning With

looptune . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-6
Set Up Your Control System for Tuning with

looptune . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-8
Set Up Your Control System for looptunein MATLAB . . . . 6-8
Set Up Your Control System for looptune in Simulink . . . . 6-8
viii Contents
Performance and Robustness Specifications for
looptune . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-10
Tune MIMO Control System for Specified

Bandwidth . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-12
What Is hinfstruct? . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-19
Formulating Design Requirements as H-Infinity

Constraints . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-20
Structured H-Infinity Synthesis Workflow . . . . . . . . . . . 6-21
Build Tunable Closed-Loop Model for Tuning with

hinfstruct . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-22
Constructing the Closed-Loop System Using Control
System Toolbox Commands . . . . . . . . . . . . . . . . . . . . . . . 6-22
Constructing the Closed-Loop System Using Simulink
Control Design Commands . . . . . . . . . . . . . . . . . . . . . . . . 6-26
Tune the Controller Parameters . . . . . . . . . . . . . . . . . . . . 6-29
Interpret the Outputs of hinfstruct . . . . . . . . . . . . . . . . . . 6-30

Output Model is Tuned Version of Input Model . . . . . . . . . 6-30
Interpreting gamma . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-30
Validate the Controller Design . . . . . . . . . . . . . . . . . . . . . . 6-31

Validating the Design in MATLAB . . . . . . . . . . . . . . . . . . . 6-31
Validating the Design in Simulink . . . . . . . . . . . . . . . . . . . 6-32
Set Up Your Control System for Tuning with systune . . 6-35

Set Up Your Control System for systune in MATLAB . . . . 6-35
Set Up Your Control System for systune in Simulink . . . . 6-35
Specifying Design Requirements for systune . . . . . . . . . 6-37
Tune Controller Against Set of Plant Models . . . . . . . . . 6-39
ix
Supported Blocks for Tuning in Simulink . . . . . . . . . . . . 6-40
Tuning Other Blocks . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-40
Speed Up Tuning with Parallel Computing Toolbox

Software . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-42
Tuning Control Systems with SYSTUNE . . . . . . . . . . . . . 6-44
Tuning Feedback Loops with LOOPTUNE . . . . . . . . . . . 6-54
Tuning Control Systems in Simulink . . . . . . . . . . . . . . . . 6-62
Building Tunable Models . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-72
Using Design Requirement Objects . . . . . . . . . . . . . . . . . . 6-81
Validating Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-102
Tuning Multi-Loop Control Systems . . . . . . . . . . . . . . . . . 6-112
Using Parallel Computing to Accelerate Tuning . . . . . . 6-124
PID Tuning for Setpoint Tracking vs. Disturbance

Rejection . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-129
Decoupling Controller for a Distillation Column . . . . . 6-141
Tuning of a Digital Motion Control System . . . . . . . . . . . 6-154
Control of a Linear Electric Actuator . . . . . . . . . . . . . . . . 6-175
Multi-Loop PID Control of a Robot Arm . . . . . . . . . . . . . 6-187
Active Vibration Control in Three-Story Building . . . . 6-205
Tuning of a Two-Loop Autopilot . . . . . . . . . . . . . . . . . . . . 6-220
x Contents
Multi-Loop Control of a Helicopter . . . . . . . . . . . . . . . . . . 6-238
Fixed-Structure Autopilot for a Passenger Jet . . . . . . . 6-248
Fault-Tolerant Control of a Passenger Jet . . . . . . . . . . . 6-262
Fixed-Structure H-infinity Synthesis with

HINFSTRUCT . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-274
MIMO Control of Diesel Engine . . . . . . . . . . . . . . . . . . . . . 6-286
Digital Control of Power Stage Voltage . . . . . . . . . . . . . . 6-304
Gain-Scheduled Controllers
7
Gain-Scheduled Control Systems . . . . . . . . . . . . . . . . . . . . 7-2
Plant Models for Gain-Scheduled Control . . . . . . . . . . . . 7-4

Gain Scheduling for Linear Parameter-Varying Plants . . . 7-4
Gain Scheduling for Nonlinear Plants . . . . . . . . . . . . . . . . . 7-5
Parametric Gain Surfaces . . . . . . . . . . . . . . . . . . . . . . . . . . 7-8
Tuning Gain-Scheduled Controllers . . . . . . . . . . . . . . . . . 7-13
Validating Gain-Scheduled Controllers . . . . . . . . . . . . . . 7-14
Improving Gain-Scheduled Tuning Results . . . . . . . . . . 7-15

Normalize the Scheduling Variables . . . . . . . . . . . . . . . . . . 7-15
Changing Requirements With Operating Condition . . . . . . 7-16
Tunable Gain Surface With Two Scheduling

Variables . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 7-18
xi
Gain-Scheduled PID Controller . . . . . . . . . . . . . . . . . . . . . 7-22
Tuning of Gain-Scheduled Three-Loop Autopilot . . . . . 7-24
Index
xii Contents
1
Introduction
• “Robust Control Toolbox Product Description” on page 1-2

• “Product Requirements” on page 1-3
• “Modeling Uncertainty” on page 1-4
• “System with Uncertain Parameters” on page 1-5
• “Worst-Case Performance” on page 1-9
• “Worst-Case Performance of Uncertain System” on page 1-10
• “Loop-Shaping Controller Design” on page 1-12
• “Model Reduction and Approximation” on page 1-16
• “Reduce Order of Synthesized Controller” on page 1-17
• “LMI Solvers” on page 1-21
• “Extends Control System Toolbox Capabilities” on page 1-22
• “Acknowledgments” on page 1-23
• “Bibliography” on page 1-25
1 Introduction
Robust Control Toolbox Product Description

Design robust controllers for uncertain plants
Robust Control Toolbox™ provides functions, algorithms, and blocks for
analyzing and tuning control systems for performance and robustness. You
can create uncertain models by combining nominal dynamics with uncertain
elements, such as uncertain parameters or unmodeled dynamics. You can
analyze the impact of plant model uncertainty on control system performance
and identify worst-case combinations of uncertain elements. H-infinity and
mu-synthesis techniques let you design controllers that maximize robust
stability and performance.
The toolbox automatically tunes both SISO and MIMO controllers. These can
include decentralized, fixed-structure controllers with multiple tunable blocks
spanning multiple feedback loops. The toolbox lets you tune one controller
against a set of plant models. You can also tune gain-scheduled controllers.
You can specify multiple tuning objectives, such as reference tracking,
disturbance rejection, stability margins, and closed-loop pole locations.
Key Features
• Modeling of systems with uncertain parameters or neglected dynamics
• Worst-case analysis of stability margins and sensitivity to disturbances
• Automatic tuning of centralized, decentralized, and multiloop controllers
• Automatic tuning of gain-scheduled controllers
• Robustness analysis and controller tuning in Simulink®
• H-infinity and mu-synthesis algorithms
• General-purpose LMI solvers
1-2
Product Requirements
Product Requirements
Robust Control Toolbox software requires that you have installed Control
System Toolbox™ software.
1-3
1 Introduction
Modeling Uncertainty
At the heart of robust control is the concept of an uncertain LTI system.
Model uncertainty arises when system gains or other parameters are not
precisely known, or can vary over a given range. Examples of real parameter
uncertainties include uncertain pole and zero locations and uncertain gains.
You can also have unstructured uncertainties, by which is meant complex
parameter variations satisfying given magnitude bounds.
With Robust Control Toolbox software you can create uncertain LTI models
as MATLAB® objects specifically designed for robust control applications. You
can build models of complex systems by combining models of subsystems
using addition, multiplication, and division, as well as with Control System
Toolbox commands like feedback and lft.
For information about LTI model types, see “Linear System Representation”.
1-4
System with Uncertain Parameters

For instance, consider the two-cart "ACC Benchmark" system [13] consisting
of two frictionless carts connected by a spring shown as follows.
ACC Benchmark Problem
The system has the block diagram model shown below, where the individual
carts have the respective transfer functions.
1
G1 ( s ) =
m1 s2
1
G2 ( s ) = .
m2 s2
The parameters m1, m2, and k are uncertain, equal to one plus or minus 20%:
m1 = 1 – 0.2
m2 = 1 – 0.2
k = 1 – 0.2
1-5
1 Introduction
"ACC Benchmark" Two-Cart System Block Diagram y1 = P(s) u1
The upper dashed-line block has transfer function matrix F(s):
⎡ 0 ⎤ ⎡1⎤
F ( s) = ⎢ ⎥ [1 −1] + ⎢ ⎥ ⎡⎣0 G2 ( s ) ⎤⎦ .
⎣G1 ( s ) ⎦ ⎣ −1⎦
This code builds the uncertain system model P shown above:
% Create the uncertain real parameters m1, m2, & k

m1 = ureal('m1',1,'percent',20);
m2 = ureal('m2',1,'percent',20);
k = ureal('k',1,'percent',20);
s = zpk('s'); % create the Laplace variable s

G1 = ss(1/s^2)/m1; % Cart 1
G2 = ss(1/s^2)/m2; % Cart 2
% Now build F and P

F = [0;G1]*[1 -1]+[1;-1]*[0,G2];
P = lft(F,k) % close the loop with the spring k
1-6
The variable P is a SISO uncertain state-space (USS) object with four states
and three uncertain parameters, m1, m2, and k. You can recover the nominal
plant with the command
zpk(P.nominal)
which returns
Zero/pole/gain:
1
--------------
s^2 (s^2 + 2)
If the uncertain model P(s) has LTI negative feedback controller
100 ( s + 1)
3
C ( s) =
( 0.001s + 1)3
then you can form the controller and the closed-loop system y1 = T(s) u1 and
view the closed-loop system’s step response on the time interval from t=0 to
t=0.1 for a Monte Carlo random sample of five combinations of the three
uncertain parameters k, m1, and m2 using this code:
C=100*ss((s+1)/(.001*s+1))^3 % LTI controller

T=feedback(P*C,1); % closed-loop uncertain system
step(usample(T,5),.1);
The resulting plot is shown below.
1-7
1 Introduction
Monte Carlo Sampling of Uncertain System’s Step Response
1-8
Worst-Case Performance
Worst-Case Performance
To be robust, your control system should meet your stability and performance
requirements for all possible values of uncertain parameters. Monte Carlo
parameter sampling via usample can be used for this purpose as shown in
Monte Carlo Sampling of Uncertain System’s Step Response on page 1-8, but
Monte Carlo methods are inherently hit or miss. With Monte Carlo methods,
you might need to take an impossibly large number of samples before you hit
upon or near a worst-case parameter combination.
Robust Control Toolbox software gives you a powerful assortment of

robustness analysis commands that let you directly calculate upper and lower
bounds on worst-case performance without random sampling.
Worst-Case Robustness Analysis Commands

loopmargin Comprehensive analysis of feedback loop
loopsens Sensitivity functions of feedback loop
ncfmargin Normalized coprime stability margin of feedback
loop
robustperf Robust performance of uncertain systems
robuststab Stability margins of uncertain systems
wcgain Worst-case gain of an uncertain system
wcmargin Worst-case gain/phase margins for feedback loop
wcsens Worst-case sensitivity functions of feedback loop
1-9
1 Introduction
Worst-Case Performance of Uncertain System

Returning to the “System with Uncertain Parameters” on page 1-5, the closed
loop system is:
T=feedback(P*C,1); % Closed-loop uncertain system
This uncertain state-space model T has three uncertain parameters, k, m1,

and m2, each equal to 1±20% uncertain variation. To analyze whether the
closed-loop system T is robustly stable for all combinations of values for these
three parameters, you can execute the commands:
[StabilityMargin,Udestab,REPORT] = robuststab(T);
REPORT
This displays the REPORT:
Uncertain System is robustly stable to modeled uncertainty.

-- It can tolerate up to 311% of modeled uncertainty.
-- A destabilizing combination of 500% the modeled uncertainty exists,
causing an instability at 44.3 rad/s.
The report tells you that the control system is robust for all parameter
variations in the ±20% range, and that the smallest destabilizing combination
of real variations in the values k, m1, and m2 has sizes somewhere between
311% and 500% greater than ±20%, i.e., between ±62.2% and ±100%. The
value Udestab returns an estimate of the 500% destabilizing parameter
variation combination:
Udestab =
k: 1.2174e-005
m1: 1.2174e-005
m2: 2.0000.
1-10
Worst-Case Performance of Uncertain System
Uncertain System Closed-Loop Bode Plots
You have a comfortable safety margin of between 311% to 500% larger

than the anticipated ±20% parameter variations before the closed loop
goes unstable. But how much can closed-loop performance deteriorate for
parameter variations constrained to lie strictly within the anticipated ±20%
range? The following code computes worst-case peak gain of T, and estimates
the frequency and parameter values at which the peak gain occurs:
[PeakGain,Uwc] = wcgain(T);
Twc=usubs(T,Uwc);
% Worst case closed-loop system T
Trand=usample(T,4);
% 4 random samples of uncertain system T
bodemag(Twc,'r',Trand,'b-.',{.5,50}); % Do bode plot
legend('T_{wc} - worst-case',...
'T_{rand} - random samples',3);
The resulting plot is shown in Uncertain System Closed-Loop Bode Plots

on page 1-11.
1-11
1 Introduction
Loop-Shaping Controller Design

One of the most powerful yet simple controller synthesis tools is loopsyn.
Given an LTI plant, you specify the shape of the open-loop systems frequency
response plot that you want, then loopsyn computes a stabilizing controller
that best approximates your specified loop shape.
For example, consider the 2-by-2 NASA HiMAT aircraft model (Safonov,
Laub, and Hartmann [8]) depicted in the following figure. The control
variables are elevon and canard actuators (δe and δc). The output variables
are angle of attack (α) and attitude angle (θ). The model has six states:
⎡ x1 ⎤ ⎡  ⎤
⎢x ⎥ ⎢  ⎥
⎢ 2⎥ ⎢ ⎥
⎢ x3 ⎥ ⎢  ⎥
x=⎢ ⎥=⎢ ⎥
⎢ x4 ⎥ ⎢  ⎥
⎢ x5 ⎥ ⎢ xe ⎥
⎢ ⎥ ⎢ ⎥
⎢⎣ x6 ⎥⎦ ⎢⎣ x ⎦⎥
where xe and xδ are elevator and canard actuator states.
1-12
Aircraft Configuration and Vertical Plane Geometry
You can enter the state-space matrices for this model with the following code:
% NASA HiMAT model G(s)

ag =[ -2.2567e-02 -3.6617e+01 -1.8897e+01 -3.2090e+01 3.2509e+00 -7.6257e-01;
9.2572e-05 -1.8997e+00 9.8312e-01 -7.2562e-04 -1.7080e-01 -4.9652e-03;
1.2338e-02 1.1720e+01 -2.6316e+00 8.7582e-04 -3.1604e+01 2.2396e+01;
0 0 1.0000e+00 0 0 0;
0 0 0 0 -3.0000e+01 0;
0 0 0 0 0 -3.0000e+01];
bg = [ 0 0;
0 0;
0 0;
0 0;
30 0;
0 30];
cg = [ 0 1 0 0 0 0;
1-13
1 Introduction
0 0 0 1 0 0];
dg = [ 0 0;
0 0];
G=ss(ag,bg,cg,dg);
To design a controller to shape the frequency response (sigma) plot so that the
system has approximately a bandwidth of 10 rad/s, you can set as your target
desired loop shape Gd(s)=10/s, then use loopsyn(G,Gd) to find a loop-shaping
controller for G that optimally matches the desired loop shape Gd by typing:
s=zpk('s'); w0=10; Gd=w0/(s+.001);

[K,CL,GAM]=loopsyn(G,Gd); % Design a loop-shaping controller K
% Plot the results

sigma(G*K,'r',Gd,'k-.',Gd/GAM,'k:',Gd*GAM,'k:',{.1,30})
figure ;T=feedback(G*K,eye(2));
sigma(T,ss(GAM),'k:',{.1,30});grid
The value of γ= GAM returned is an indicator of the accuracy to which

the optimal loop shape matches your desired loop shape and is an upper
bound on the resonant peak magnitude of the closed-loop transfer function
T=feedback(G*K,eye(2)). In this case, γ = 1.6024 = 4 dB — see the next
figure.
1-14
MIMO Robust Loop Shaping with loopsyn(G,Gd)
The achieved loop shape matches the desired target Gd to within about γ dB.
1-15
1 Introduction
Model Reduction and Approximation

Complex models are not always required for good control. Unfortunately,
however, optimization methods (including methods based on H∞, H2, and
µ-synthesis optimal control theory) generally tend to produce controllers
with at least as many states as the plant model. For this reason, Robust
Control Toolbox software offers you an assortment of model-order reduction
commands that help you to find less complex low-order approximations to
plant and controller models.
Model Reduction Commands

reduce Main interface to model approximation algorithms
balancmr Balanced truncation model reduction
bstmr Balanced stochastic truncation model reduction
hankelmr Optimal Hankel norm model approximations
modreal State-space modal truncation/realization
ncfmr Balanced normalized coprime factor model reduction
schurmr Schur balanced truncation model reduction
slowfast State-space slow-fast decomposition
stabsep State-space stable/antistable decomposition
imp2ss Impulse response to state-space approximation
Among the most important types of model reduction methods are minimize
bounds methods on additive, multiplicative, and normalized coprime factor
(NCF) model error. You can access all three of these methods using the
command reduce.
1-16
Reduce Order of Synthesized Controller

For instance, the NASA HiMAT model considered in “Loop-Shaping Controller
Design” on page 1-12 has eight states, and the optimal loop-shaping controller
turns out to have 16 states. Using model reduction, you can remove at least
some of the states without appreciably affecting stability or closed-loop
performance. For controller order reduction, the NCF model reduction is
particularly useful, and it works equally well with controllers that have poles
anywhere in the complex plane.
For the NASA HiMAT design in the last section, you can type
hankelsv(K,'ncf','log');
which displays a logarithmic plot of the NCF Hankel singular values — see
the following figure.
1-17
1 Introduction
Hankel Singular Values of Coprime Factorization of K
Theory says that, without danger of inducing instability, you can confidently
discard at least those controller states that have NCF Hankel singular values
that are much smaller than ncfmargin(G,K).
Compute ncfmargin(G,K) and add it to your Hankel singular values plot.
hankelsv(K,'ncf','log');v=axis;
hold on; plot(v(1:2), ncfmargin(G,K)*[1 1],'--'); hold off
1-18
Five of the 16 NCF Hankel Singular Values of HiMAT Controller K Are Small
Compared to ncfmargin(G,K)
In this case, you can safely discard 5 of the 16 states of K and compute an
11-state reduced controller by typing:
K1=reduce(K,11,'errortype','ncf');
The result is plotted in the following figure.
sigma(G*K1,'b',G*K,'r--',{.1,30});
1-19
1 Introduction
HiMAT with 11-State Controller K1 vs. Original 16-State Controller K
The picture above shows that low-frequency gain is decreased considerably for
inputs in one vector direction. Although this does not affect stability, it affects
performance. If you wanted to better preserve low-frequency performance,
you would discard fewer than five of the 16 states of K.
1-20
LMI Solvers
LMI Solvers
At the core of many emergent robust control analysis and synthesis routines
are powerful general-purpose functions for solving a class of convex nonlinear
programming problems known as linear matrix inequalities. The LMI
capabilities are invoked by Robust Control Toolbox software functions
that evaluate worst-case performance, as well as functions like hinfsyn
and h2hinfsyn. Some of the main functions that help you access the LMI
capabilities of the toolbox are shown in the following table.
Specification of LMIs
lmiedit GUI for LMI specification
setlmis Initialize the LMI description
lmivar Define a new matrix variable
lmiterm Specify the term content of an LMI
newlmi Attach an identifying tag to new LMIs
getlmis Get the internal description of the LMI system
LMI Solvers
feasp Test feasibility of a system of LMIs
gevp Minimize generalized eigenvalue with LMI constraints
mincx Minimize a linear objective with LMI constraints
dec2mat Convert output of the solvers to values of matrix
variables
Evaluation of LMIs/Validation of Results

evallmi Evaluate for given values of the decision variables
showlmi Return the left and right sides of an evaluated LMI
1-21
1 Introduction
Extends Control System Toolbox Capabilities

Robust Control Toolbox software is designed to work with Control System
Toolbox software. Robust Control Toolbox software extends the capabilities
of Control System Toolbox software and leverages the LTI and plotting
capabilities of Control System Toolbox software. The major analysis and
synthesis commands in Robust Control Toolbox software accept LTI object
inputs, e.g., LTI state-space systems produced by commands such as:
G=tf(1,[1 2 3])
G=ss([-1 0; 0 -1], [1;1],[1 1],3)
The uncertain system (USS) objects in Robust Control Toolbox software

generalize the Control System Toolbox LTI SS objects and help ease the
task of analyzing and plotting uncertain systems. You can do many of
the same algebraic operations on uncertain systems that are possible for
LTI objects (multiply, add, invert), and Robust Control Toolbox software
provides USS uncertain system extensions of Control System Toolbox software
interconnection and plotting functions like feedback, lft, and bode.
1-22
Acknowledgments
Acknowledgments
Professor Andy Packard is with the Faculty of Mechanical Engineering
at the University of California, Berkeley. His research interests include
robustness issues in control analysis and design, linear algebra and numerical
algorithms in control problems, applications of system theory to aerospace
problems, flight control, and control of fluid flow.
Professor Gary Balas is with the Faculty of Aerospace Engineering &

Mechanics at the University of Minnesota and is president of MUSYN Inc.
His research interests include aerospace control systems, both experimental
and theoretical.
Dr. Michael Safonov is with the Faculty of Electrical Engineering at the

University of Southern California. His research interests include control
and decision theory.
Dr. Richard Chiang is employed by Boeing Satellite Systems, El Segundo,

CA. He is a Boeing Technical Fellow and has been working in the aerospace
industry over 25 years. In his career, Richard has designed 3 flight control
laws, 12 spacecraft attitude control laws, and 3 large space structure vibration
controllers, using modern robust control theory and the tools he built in this
toolbox. His research interests include robust control theory, model reduction,
and in-flight system identification. Working in industry instead of academia,
Richard serves a unique role in our team, bridging the gap between theory
and reality.
The linear matrix inequality (LMI) portion of Robust Control Toolbox software
was developed by these two authors:
Dr. Pascal Gahinet is employed by MathWorks. His research interests

include robust control theory, linear matrix inequalities, numerical linear
algebra, and numerical software for control.
Professor Arkadi Nemirovski is with the Faculty of Industrial Engineering

and Management at Technion, Haifa, Israel. His research interests include
convex optimization, complexity theory, and nonparametric statistics.
1-23
1 Introduction
The structured H∞ synthesis (hinfstruct) portion of Robust Control Toolbox

software was developed by the following author in collaboration with Pascal
Gahinet:
Professor Pierre Apkarian is with ONERA (The French Aerospace Lab)

and the Institut de Mathématiques at Paul Sabatier University, Toulouse,
France. His research interests include robust control, LMIs, mathematical
programming, and nonsmooth optimization techniques for control.
1-24
Bibliography
Bibliography
[1] Boyd, S.P., El Ghaoui, L., Feron, E., and Balakrishnan, V., Linear Matrix
Inequalities in Systems and Control Theory, Philadelphia, PA, SIAM, 1994.
[2] Dorato, P. (editor), Robust Control, New York, IEEE Press, 1987.
[3] Dorato, P., and Yedavalli, R.K. (editors), Recent Advances in Robust
Control, New York, IEEE Press, 1990.
[4] Doyle, J.C., and Stein, G., “Multivariable Feedback Design: Concepts
for a Classical/Modern Synthesis,” IEEE Trans. on Automat. Contr., 1981,
AC-26(1), pp. 4-16.
[5] El Ghaoui, L., and Niculescu, S., Recent Advances in LMI Theory for
Control, Philadelphia, PA, SIAM, 2000.
[6] Lehtomaki, N.A., Sandell, Jr., N.R., and Athans, M., “Robustness Results
in Linear-Quadratic Gaussian Based Multivariable Control Designs,” IEEE
Trans. on Automat. Contr., Vol. AC-26, No. 1, Feb. 1981, pp. 75-92.
[7] Safonov, M.G., Stability and Robustness of Multivariable Feedback

Systems, Cambridge, MA, MIT Press, 1980.
[8] Safonov, M.G., Laub, A.J., and Hartmann, G., “Feedback Properties of
Multivariable Systems: The Role and Use of Return Difference Matrix,” IEEE
Trans. of Automat. Contr., 1981, AC-26(1), pp. 47-65.
[9] Safonov, M.G., Chiang, R.Y., and Flashner, H., “H∞ Control Synthesis
for a Large Space Structure,” Proc. of American Contr. Conf., Atlanta, GA,
June 15-17, 1988.
[10] Safonov, M.G., and Chiang, R.Y., “CACSD Using the State-Space L∞
Theory — A Design Example,” IEEE Trans. on Automatic Control, 1988,
AC-33(5), pp. 477-479.
[11] Sanchez-Pena, R.S., and Sznaier, M., Robust Systems Theory and
Applications, New York, Wiley, 1998.
1-25
1 Introduction
[12] Skogestad, S., and Postlethwaite, I., Multivariable Feedback Control,

New York, Wiley, 1996.
[13] Wie, B., and Bernstein, D.S., “A Benchmark Problem for Robust
Controller Design,” Proc. American Control Conf., San Diego, CA, May 23-25,
1990; also Boston, MA, June 26-28, 1991.
[14] Zhou, K., Doyle, J.C., and Glover, K., Robust and Optimal Control,
Englewood Cliffs, NJ, Prentice Hall, 1996.
1-26
2
Multivariable Loop Shaping
• “Tradeoff Between Performance and Robustness” on page 2-2

• “Norms and Singular Values” on page 2-4
• “Typical Loop Shapes, S and T Design” on page 2-6
• “Using LOOPSYN to Do H-Infinity Loop Shaping” on page 2-15
• “Loop-Shaping Control Design of Aircraft Model” on page 2-16
• “Fine-Tuning the LOOPSYN Target Loop Shape Gd to Meet Design Goals”
on page 2-21
• “Mixed-Sensitivity Loop Shaping” on page 2-22
• “Mixed-Sensitivity Loop-Shaping Controller Design” on page 2-24
• “Loop-Shaping Commands” on page 2-26
2 Multivariable Loop Shaping
Tradeoff Between Performance and Robustness

When the plant modeling uncertainty is not too big, you can design high-gain,
high-performance feedback controllers. High loop gains significantly larger
than 1 in magnitude can attenuate the effects of plant model uncertainty
and reduce the overall sensitivity of the system to plant noise. But if your
plant model uncertainty is so large that you do not even know the sign of
your plant gain, then you cannot use large feedback gains without the risk
that the system will become unstable. Thus, plant model uncertainty can
be a fundamental limiting factor in determining what can be achieved with
feedback.
Multiplicative Uncertainty: Given an approximate model of the plant G0

of a plant G, the multiplicative uncertainty ΔM of the model G0 is defined
as Δ M = G0−1 ( G − G0 )
or, equivalently,
G = ( I + Δ M ) G0 .
Plant model uncertainty arises from many sources. There might be

small unmodeled time delays or stray electrical capacitance. Imprecisely
understood actuator time constants or, in mechanical systems, high-frequency
torsional bending modes and similar effects can be responsible for plant
model uncertainty. These types of uncertainty are relatively small at lower
frequencies and typically increase at higher frequencies.
In the case of single-input/single-output (SISO) plants, the frequency at

which there are uncertain variations in your plant of size |ΔM|=2 marks
a critical threshold beyond which there is insufficient information about
the plant to reliably design a feedback controller. With such a 200% model
uncertainty, the model provides no indication of the phase angle of the true
plant, which means that the only way you can reliably stabilize your plant is
to ensure that the loop gain is less than 1. Allowing for an additional factor
of 2 margin for error, your control system bandwidth is essentially limited
2-2
Tradeoff Between Performance and Robustness
to the frequency range over which your multiplicative plant uncertainty ΔM

has gain magnitude |ΔM|<1.
2-3
Norms and Singular Values

For MIMO systems the transfer functions are matrices, and relevant
measures of gain are determined by singular values, H∞, and H2 norms, which
are defined as follows:
H2 and H• Norms The H2-norm is the energy of the impulse response

of plant G. The H∞-norm is the peak gain of G across all frequencies and all
input directions.
Another important concept is the notion of singular values.
Singular Values: The singular values of a rank r matrix A ∈ C m×n , denoted

σi, are the nonnegative square roots of the eigenvalues of A* A ordered such
that σ1 ≥ σ2 ≥ ... ≥σp > 0, p ≤ min{m, n}.
If r < p then there are p – r zero singular values, i.e., σr+1 = σr+2 = ... =σp = 0.
The greatest singular value σ1 is sometimes denoted
 ( A ) = 1.
When A is a square n-by-n matrix, then the nth singular value (i.e., the least
singular value) is denoted
 ( A)  n.
Properties of Singular Values

Some useful properties of singular values are:
2-4
Norms and Singular Values
Ax
 ( A ) = max x∈C h
x
Ax
 ( A ) = min x∈C h
x
These properties are especially important because they establish that the
greatest and least singular values of a matrix A are the maximal and minimal
"gains" of the matrix as the input vector x varies over all possible directions.
For stable continuous-time LTI systems G(s), the H2-norm and the H∞-norms
are defined terms of the frequency-dependent singular values of G(jω):
H2-norm:
p
⎡ 1 ⎤ ∞
( )
2
G 2 ⎢⎣ 2 ⎥⎦ ∫−∞ ∑  i ( G ( j ) ) d
i=1
H∞-norm:
G 2
sup  ( G ( j ) )

where sup denotes the least upper bound.
2-5
Typical Loop Shapes, S and T Design

Consider the multivariable feedback control system shown in the following
figure. In order to quantify the multivariable stability margins and
performance of such systems, you can use the singular values of the closed-loop
transfer function matrices from r to each of the three outputs e, u, and y, viz.
def −1
S ( s) = ( I + L ( s))
def −1
R ( s) = K ( s) ( I + L ( s))
def −1
T ( s) = L ( s) ( I + L ( s)) = I − S ( s)
where the L(s) is the loop transfer function matrix
L ( s) = G ( s) K ( s). (2-1)
Block Diagram of the Multivariable Feedback Control System
The two matrices S(s) and T(s) are known as the sensitivity function and
complementary sensitivity function, respectively. The matrix R(s) has no
common name. The singular value Bode plots of each of the three transfer
function matrices S(s), R(s), and T(s) play an important role in robust
multivariable control system design. The singular values of the loop transfer
function matrix L(s) are important because L(s) determines the matrices
S(s) and T(s).
2-6
Robustness in Terms of Singular Values

The singular values of S(jω) determine the disturbance attenuation, because
S(s) is in fact the closed-loop transfer function from disturbance d to plant
output y — see Block Diagram of the Multivariable Feedback Control System
on page 2-6. Thus a disturbance attenuation performance specification can
be written as
 ( S ( j ) ) ≤ W1−1 ( j ) (2-2)
where W1−1 ( j ) is the desired disturbance attenuation factor. Allowing
W1 ( j ) to depend on frequency ω enables you to specify a different

attenuation factor for each frequency ω.
The singular value Bode plots of R(s) and of T(s) are used to measure
the stability margins of multivariable feedback designs in the face of
additive plant perturbations ΔA and multiplicative plant perturbations ΔM,
respectively. See the following figure.
Consider how the singular value Bode plot of complementary sensitivity T(s)
determines the stability margin for multiplicative perturbations ΔM. The
multiplicative stability margin is, by definition, the "size" of the smallest
stable ΔM(s) that destabilizes the system in the figure below when ΔA = 0.
2-7
Additive/Multiplicative Uncertainty
Taking  ( Δ M ( j ) ) to be the definition of the "size" of ΔM(jω), you have the

following useful characterization of "multiplicative" stability robustness:
Multiplicative Robustness: The size of the smallest destabilizing

multiplicative uncertainty ΔM(s) is:
1
 ( Δ M ( j ) ) = .
 ( T ( j ) )
The smaller is  ( T ( j ) ) , the greater will be the size of the smallest

destabilizing multiplicative perturbation, and hence the greater will be the
stability margin.
A similar result is available for relating the stability margin in the face of
additive plant perturbations ΔA(s) to R(s) if you take  ( Δ A ( j ) ) to be the

definition of the "size" of ΔA(jω) at frequency ω.
2-8
Additive Robustness: The size of the smallest destabilizing additive

uncertainty ΔA is:
1
 ( Δ A ( j ) ) = .
 ( R ( j ) )
As a consequence of robustness theorems 1 and 2, it is common to specify the

stability margins of control systems via singular value inequalities such as
 ( R { j}) ≤ W2−1 ( j ) (2-3)
 ( T { j}) ≤ W3−1 ( j ) (2-4)
where |W2(jω)| and |W3(jω)| are the respective sizes of the largest
anticipated additive and multiplicative plant perturbations.
It is common practice to lump the effects of all plant uncertainty into a

single fictitious multiplicative perturbation ΔM, so that the control design
requirements can be written
1
≥ W1 ( j ) ;  i ( T [ j ]) ≤ W3−1 ( j )
 i ( S ( j ) )
as shown in Singular Value Specifications on L, S, and T on page 2-12.
It is interesting to note that in the upper half of the figure (above the 0 dB
line),
1
 ( L ( j ) ) ≈
 ( S ( j ) )
while in the lower half of Singular Value Specifications on L, S, and T on page

2-12 (below the 0 dB line),
2-9
 ( L ( j ) ) ≈  ( T ( j ) ) .
This results from the fact that
def −1 −1
S ( s) = ( I + L ( s)) ≈ L ( s)
if  ( L ( s ) ) 1 , and
def −1
T ( s) = L ( s) ( I + L ( s)) ≈ L ( s)
if  ( L ( s ) ) 1.
2-10
2-11
Singular Value Specifications on L, S, and T
Thus, it is not uncommon to see specifications on disturbance attenuation

and multiplicative stability margin expressed directly in terms of forbidden
regions for the Bode plots of σi(L(jω)) as "singular value loop shaping"
requirements, either as specified upper/lower bounds or as a target desired
loop shape — see the preceding figure.
Guaranteed Gain/Phase Margins in MIMO Systems

For those who are more comfortable with classical single-loop concepts, there
are the important connections between the multiplicative stability margins
predicted by  ( T ) and those predicted by classical M-circles, as found on the
Nichols chart. Indeed in the single-input/single-output case,
2-12
L ( j )
 ( T ( j ) ) =
1 + L ( j )
which is precisely the quantity you obtain from Nichols chart M-circles. Thus,
T ∞ is a multiloop generalization of the closed-loop resonant peak magnitude

which, as classical control experts will recognize, is closely related to the
damping ratio of the dominant closed-loop poles. Also, it turns out that you
can relate T ∞ , S ∞ to the classical gain margin GM and phase margin θM in

each feedback loop of the multivariable feedback system of Block Diagram of
the Multivariable Feedback Control System on page 2-6 via the formulas:
1
GM ≥ 1 +
T∞
1
GM ≥ 1 +
1
1−
S∞
⎛ 1 ⎞
 M ≥ 2 sin −1 ⎜ ⎟
⎜2 T ⎟
⎝ ∞ ⎠
⎛ 1 ⎞
 M ≥ 2 sin −1 ⎜ ⎟.
⎜2 T ⎟
⎝ ∞ ⎠
(See [6].) These formulas are valid provided S ∞ and T ∞ are larger than 1,
as is normally the case. The margins apply even when the gain perturbations
or phase perturbations occur simultaneously in several feedback channels.
The infinity norms of S and T also yield gain reduction tolerances. The gain
reduction tolerance gm is defined to be the minimal amount by which the gains
in each loop would have to be decreased in order to destabilize the system.
Upper bounds on gm are as follows:
2-13
1
gM ≤ 1 −
T∞
1
gM ≤ .
1
1+
S∞
2-14
Using LOOPSYN to Do H-Infinity Loop Shaping
Using LOOPSYN to Do H-Infinity Loop Shaping

The command loopsyn lets you design a stabilizing feedback controller to
optimally shape the open loop frequency response of a MIMO feedback control
system to match as closely as possible a desired loop shape Gd — see the
preceding figure. The basic syntax of the loopsyn loop-shaping controller
synthesis command is:
K = loopsyn(G,Gd)
Here G is the LTI transfer function matrix of a MIMO plant model, Gd is

the target desired loop shape for the loop transfer function L=G*K, and K is
the optimal loop-shaping controller. The LTI controller K has the property
that it shapes the loop L=G*K so that it matches the frequency response of Gd
as closely as possible, subject to the constraint that the compensator must
stabilize the plant model G.
2-15
Loop-Shaping Control Design of Aircraft Model

To see how the loopsyn command works in practice to address robustness
and performance tradeoffs, consider again the NASA HiMAT aircraft model
taken from the paper of Safonov, Laub, and Hartmann [8]. The longitudinal
dynamics of the HiMAT aircraft trimmed at 25000 ft and 0.9 Mach are
unstable and have two right-half-plane phugoid modes. The linear model
has state-space realization G(s) = C(Is – A)–1B with six states, with the first
four states representing angle of attack (α) and attitude angle (θ) and their
rates of change, and the last two representing elevon and canard control
actuator dynamics — see Aircraft Configuration and Vertical Plane Geometry
on page 2-17.
ag =[
-2.2567e-02 -3.6617e+01 -1.8897e+01 -3.2090e+01 3.2509e+00 -7.6257e-01;
9.2572e-05 -1.8997e+00 9.8312e-01 -7.2562e-04 -1.7080e-01 -4.9652e-03;
1.2338e-02 1.1720e+01 -2.6316e+00 8.7582e-04 -3.1604e+01 2.2396e+01;
0 0 1.0000e+00 0 0 0;
0 0 0 0 -3.0000e+01 0;
0 0 0 0 0 -3.0000e+01];
bg = [0 0;
0 0;
0 0;
0 0;
30 0;
0 30];
cg = [0 1 0 0 0 0;
0 0 0 1 0 0];
dg = [0 0;
0 0];
G=ss(ag,bg,cg,dg);
The control variables are elevon and canard actuators (δe and δc). The output
variables are angle of attack (α) and attitude angle (θ).
2-16
Aircraft Configuration and Vertical Plane Geometry
This model is good at frequencies below 100 rad/s with less than 30%
variation between the true aircraft and the model in this frequency range.
However as noted in [8], it does not reliably capture very high-frequency
behaviors, because it was derived by treating the aircraft as a rigid body and
neglecting lightly damped fuselage bending modes that occur at somewhere
between 100 and 300 rad/s. These unmodeled bending modes might cause as
much as 20 dB deviation (i.e., 1000%) between the frequency response of
the model and the actual aircraft for frequency ω > 100 rad/s. Other effects
like control actuator time delays and fuel sloshing also contribute to model
inaccuracy at even higher frequencies, but the dominant unmodeled effects
are the fuselage bending modes. You can think of these unmodeled bending
modes as multiplicative uncertainty of size 20 dB, and design your controller
using loopsyn, by making sure that the loop has gain less than –20 dB at, and
beyond, the frequency ω > 100 rad/s.
2-17
Design Specifications
The singular value design specifications are
• Robustness Spec.: –20 dB/decade roll-off slope and –20 dB loop gain
at 100 rad/s
• Performance Spec.: Minimize the sensitivity function as much as
possible.
Both specs can be accommodated by taking as the desired loop shape
Gd(s)=8/s
MATLAB Commands for a LOOPSYN Design

%% Enter the desired loop shape Gd
s=zpk('s'); % Laplace variable s
Gd=8/s; % desired loop shape
%% Compute the optimal loop shaping controller K
[K,CL,GAM]=loopsyn(G,Gd);
%% Compute the loop L, sensitivity S and
%% complementary sensitivity T:
L=G*K;
I=eye(size(L));
S=feedback(I,L); % S=inv(I+L);
T=I-S;
%% Plot the results:
% step response plots
step(T);title('\alpha and \theta command step responses');
% frequency response plots
figure;
sigma(I+L,'--',T,':',L,'r--',Gd,'k-.',Gd/GAM,'k:',...
Gd*GAM,'k:',{.1,100});grid
legend('1/\sigma(S) performance',...
'\sigma(T) robustness',...
'\sigma(L) loopshape',...
'\sigma(Gd) desired loop',...
'\sigma(Gd) \pm GAM, dB');
2-18
The plots of the resulting step and frequency response for the NASA HiMAT
loopsyn loop-shaping controller design are shown in the following figure. The
number ±GAM, dB (i.e., 20log10(GAM)) tells you the accuracy with which
your loopsyn control design matches the target desired loop:
 ( GK ) , db ≥ G d , db − GAM, db ( < c )
 ( GK ) , db ≥ G d , db + GAM, db ( > c ).
HiMAT Closed Loop Step Responses
2-19
LOOPSYN Design Results for NASA HiMAT
2-20
Fine-Tuning the LOOPSYN Target Loop Shape Gd to Meet Design Goals
Fine-Tuning the LOOPSYN Target Loop Shape Gd to Meet

Design Goals
If your first attempt at loopsyn design does not achieve everything you
wanted, you will need to readjust your target desired loop shape Gd. Here are
some basic design tradeoffs to consider:
• Stability Robustness. Your target loop Gd should have low gain (as small
as possible) at high frequencies where typically your plant model is so poor
that its phase angle is completely inaccurate, with errors approaching
±180° or more.
• Performance. Your Gd loop should have high gain (as great as possible)
at frequencies where your model is good, in order to ensure good control
accuracy and good disturbance attenuation.
• Crossover and Roll-Off. Your desired loop shape Gd should have its 0 dB
crossover frequency (denoted ωc) between the above two frequency ranges,
and below the crossover frequency ωc it should roll off with a negative slope
of between –20 and –40 dB/decade, which helps to keep phase lag to less
than –180° inside the control loop bandwidth (0 < ω < ωc).
Other considerations that might affect your choice of Gd are the

right-half-plane poles and zeros of the plant G, which impose ffundamental
limits on your 0 dB crossover frequency ωc [12]. For instance, your 0 dB
crossover ωc must be greater than the magnitude of any plant right-half-plane
poles and less than the magnitude of any right-half-plane zeros.
max pi < c < min zi .

Re( pi )>0 Re( zi )>0
If you do not take care to choose a target loop shape Gd that conforms to
these fundamental constraints, then loopsyn will still compute the optimal
loop-shaping controller K for your Gd, but you should expect that the optimal
loop L=G*K will have a poor fit to the target loop shape Gd, and consequently it
might be impossible to meet your performance goals.
2-21
Mixed-Sensitivity Loop Shaping

A popular alternative approach to loopsyn loop shaping is H∞ mixed-sensitivity
loop shaping, which is implemented by the Robust Control Toolbox software
command:
K=mixsyn(G,W1,[],W3)
With mixsyn controller synthesis, your performance and stability robustness

specifications equations (2-2) and (2-4) are combined into a single infinity
norm specification of the form
Ty1u1 ≤1
∞
where (see MIXSYN H∞ Mixed-Sensitivity Loop Shaping Ty1 u1 on page 2-23):
def ⎡ W S⎤
1
Ty1u1 = ⎢ .
⎣W3 T ⎥⎦
The term Ty1u1 is called a mixed-sensitivity cost function because it

∞
penalizes both sensitivity S(s) and complementary sensitivity T(s). Loop
shaping is achieved when you choose W1 to have the target loop shape for
frequencies ω < ωc, and you choose 1/W3 to be the target for ω > ωc. In choosing
design specifications W1 and W3 for a mixsyn controller design, you need to
ensure that your 0 dB crossover frequency for the Bode plot of W1 is below the
0 dB crossover frequency of 1/W3, as shown in Singular Value Specifications
on L, S, and T on page 2-12, so that there is a gap for the desired loop shape
Gd to pass between the performance bound W1 and your robustness bound
W3−1 . Otherwise, your performance and robustness requirements will not

be achievable.
2-22
Mixed-Sensitivity Loop Shaping
MIXSYN H• Mixed-Sensitivity Loop Shaping Ty1 u1
2-23
Mixed-Sensitivity Loop-Shaping Controller Design

To do a mixsyn H∞ mixed-sensitivity synthesis design on the HiMAT model,
start with the plant model G created in “Mixed-Sensitivity Loop-Shaping
Controller Design” on page 2-24 and type the following commands:
% Set up the performance and robustness bounds W1 & W3

s=zpk('s'); % Laplace variable s
MS=2;AS=.03;WS=5;
W1=(s/MS+WS)/(s+AS*WS);
MT=2;AT=.05;WT=20;
W3=(s+WT/MT)/(AT*s+WT);
% Compute the H-infinity mixed-sensitivity optimal sontroller K1
[K1,CL1,GAM1]=mixsyn(G,W1,[],W3);
% Next compute and plot the closed-loop system.
% Compute the loop L1, sensitivity S1, and comp sensitivity T1:
L1=G*K1;
I=eye(size(L1));
S1=feedback(I,L1); % S=inv(I+L1);
T1=I-S1;
% Plot the results:
% step response plots
step(T1,1.5);
title('\alpha and \theta command step responses');
% frequency response plots
figure;
sigma(I+L1,'--',T1,':',L1,'r--',... W1/GAM1,'k--',GAM1/W3,'k-.',{.1,100});grid
legend('1/\sigma(S) performance',...
'\sigma(T) robustness',...
'\sigma(L) loopshape',...
'\sigma(W1) performance bound',...
'\sigma(1/W3) robustness bound');
The resulting mixsyn singular value plots for the NASA HiMAT model are
shown below.
2-24
Mixed-Sensitivity Loop-Shaping Controller Design
MIXSYN Design Results for NASA HiMAT
2-25
Loop-Shaping Commands
Robust Control Toolbox software gives you several choices for shaping the
frequency response properties of multiinput/multioutput (MIMO) feedback
control loops. Some of the main commands that you are likely to use for
loop-shaping design, and associated utility functions, are listed below:
MIMO Loop-Shaping Commands

loopsyn H∞ loop-shaping controller synthesis
ltrsyn LQG loop-transfer recovery
mixsyn H∞ mixed-sensitivity controller synthesis
ncfsyn Glover-McFarlane H∞ normalized coprime factor loop-
shaping controller synthesis
MIMO Loop-Shaping Utility Functions

augw Augmented plant for weighted H2 and H∞ mixed-
sensitivity control synthesis
makeweight Weights for H∞ mixed sensitivity (mixsyn, augw)
sigma Singular value plots of LTI feedback loops
2-26
3
Model Reduction for Robust

Control
• “Why Reduce Model Order?” on page 3-2

• “Hankel Singular Values” on page 3-3
• “Model Reduction Techniques” on page 3-5
• “Approximate Plant Model by Additive Error Methods” on page 3-7
• “Approximate Plant Model by Multiplicative Error Method” on page 3-9
• “Using Modal Algorithms” on page 3-12
• “Reducing Large-Scale Models” on page 3-16
• “Normalized Coprime Factor Reduction” on page 3-17
3 Model Reduction for Robust Control
Why Reduce Model Order?

In the design of robust controllers for complicated systems, model reduction
fits several goals:
1 To simplify the best available model in light of the purpose for which the
model is to be used—namely, to design a control system to meet certain
specifications.
2 To speed up the simulation process in the design validation stage, using a

smaller size model with most of the important system dynamics preserved.
3 Finally, if a modern control method such as LQG or H∞ is used for which

the complexity of the control law is not explicitly constrained, the order of
the resultant controller is likely to be considerably greater than is truly
needed. A good model reduction algorithm applied to the control law can
sometimes significantly reduce control law complexity with little change in
control system performance.
Model reduction routines in this toolbox can be put into two categories:
• Additive error method — The reduced-order model has an additive error

bounded by an error criterion.
• Multiplicative error method — The reduced-order model has a
multiplicative or relative error bounded by an error criterion.
The error is measured in terms of peak gain across frequency (H∞ norm), and
the error bounds are a function of the neglected Hankel singular values.
3-2
Hankel Singular Values
Hankel Singular Values

In control theory, eigenvalues define a system stability, whereas Hankel
singular values define the “energy” of each state in the system. Keeping
larger energy states of a system preserves most of its characteristics in terms
of stability, frequency, and time responses. Model reduction techniques
presented here are all based on the Hankel singular values of a system.
They can achieve a reduced-order model that preserves the majority of the
system characteristics.
Mathematically, given a stable state-space system (A,B,C,D), its Hankel

singular values are defined as [1]
 H = i ( PQ )
where P and Q are controllability and observability grammians satisfying
AP + PAT = − BBT
AT Q + QA = −C T C.
For example,
rng(1234,'twister');
G = rss(30,4,3);
hankelsv(G)
returns a Hankel singular value plot as follows:
3-3
which shows that system G has most of its “energy” stored in states 1 through
15 or so. Later, you will see how to use model reduction routines to keep a
15-state reduced model that preserves most of its dynamic response.
3-4
Model Reduction Techniques
Model Reduction Techniques

Robust Control Toolbox software offers several algorithms for model
approximation and order reduction. These algorithms let you control the
absolute or relative approximation error, and are all based on the Hankel
singular values of the system.
Robust control theory quantifies a system uncertainty as either additive or

multiplicative types. These model reduction routines are also categorized
into two groups: additive error and multiplicative error types. In other
words, some model reduction routines produce a reduced-order model Gred
of the original model G with a bound on the error G − Gred ∞ , the peak gain
across frequency. Others produce a reduced-order model with a bound on
the relative error G −1 ( G − Gred ) .

∞
These theoretical bounds are based on the “tails” of the Hankel singular
values of the model, i.e.,
Additive Error Bound

n
G − Gred ∞
≤ 2∑ i
k+1 (3-1)
where σi are denoted the ith Hankel singular value of the original system G.
Multiplicative (Relative) Error Bound
( 1+ )
n
G −1 ( G − Gred ) ≤ ∏ ⎛⎜ 1 + 2 i 2
i +  i ⎞⎟ − 1
k+1 ⎝ ⎠
∞ (3-2)
where σi are denoted the ith Hankel singular value of the phase matrix of the
model G (see the bstmr reference page).
3-5
Top-Level Model Reduction Command
Method Description
reduce Main interface to model approximation algorithms
Normalized Coprime Balanced Model Reduction Command
Method Description
ncfmr Normalized coprime balanced truncation
Additive Error Model Reduction Commands
Method Description
balancmr Square-root balanced model truncation
schurmr Schur balanced model truncation
hankelmr Hankel minimum degree approximation
Multiplicative Error Model Reduction Command
Method Description
bstmr Balanced stochastic truncation
Additional Model Reduction Tools
Method Description
modreal Modal realization and truncation
slowfast Slow and fast state decomposition
stabsep Stable and antistable state projection
3-6
Approximate Plant Model by Additive Error Methods
Approximate Plant Model by Additive Error Methods

Given a system in LTI form, the following commands reduce the system to
any desired order you specify. The judgment call is based on its Hankel
singular values as shown in the previous paragraph.
G = rss(30,4,3);
% balanced truncation to models with sizes 12:16
[g1,info1] = balancmr(G,12:16); % or use reduce
% Schur balanced truncation by specifying `MaxError'
[g2,info2] = schurmr(G,'MaxError',[1,0.8,0.5,0.2]);
sigma(G,'b-',g1,'r--',g2,'g-.')
shows a comparison plot of the original model G and reduced models g1 and g2.
3-7
To determine whether the theoretical error bound is satisfied,
norm(G-g1(:,:,1),'inf') % 1.2556
info1.ErrorBound(1) % 6.2433
or plot the model error vs. error bound via the following commands:
[sv,w] = sigma(G-g1(:,:,1));
loglog(w,sv,w,info1.ErrorBound(1)*ones(size(w)))
xlabel('rad/sec');ylabel('SV');
title('Error Bound and Model Error')
3-8
Approximate Plant Model by Multiplicative Error Method

In most cases, multiplicative error model reduction method bstmr tends to
bound the relative error between the original and reduced-order models
across the frequency range of interest, hence producing a more accurate
reduced-order model than the additive error methods. This characteristic is
obvious in system models with low damped poles.
The following commands illustrate the significance of a multiplicative error

model reduction method as compared to any additive error type. Clearly, the
phase-matching algorithm using bstmr provides a better fit in the Bode plot.
G = rss(30,1,1);
[gr,infor] = reduce(G,'algo','balance','order',7);
[gs,infos] = reduce(G,'algo','bst','order',7);
figure(1);bode(G,'b-',gr,'r--');
title('Additive Error Method')
figure(2);bode(G,'b-',gs,'r--');
title('Relative Error Method')
3-9
3-10
Therefore, for some systems with low damped poles/zeros, the balanced
stochastic method (bstmr) produces a better reduced-order model fit in those
frequency ranges to make multiplicative error small. Whereas additive error
methods such as balancmr, schurmr, or hankelmr only care about minimizing
the overall “absolute” peak error, they can produce a reduced-order model
missing those low damped poles/zeros frequency regions.
3-11
Using Modal Algorithms
Rigid Body Dynamics

In many cases, a model’s jω-axis poles are important to keep after model
reduction, e.g., rigid body dynamics of a flexible structure plant or integrators
of a controller. A unique routine, modreal, serves the purpose nicely.
modreal puts a system into its modal form, with eigenvalues appearing on
the diagonal of its A-matrix. Real eigenvalues appear in 1-by-1 blocks, and
complex eigenvalues appear in 2-by-2 real blocks. All the blocks are ordered
in ascending order, based on their eigenvalue magnitudes, by default, or
descending order, based on their real parts. Therefore, specifying the number
of jω-axis poles splits the model into two systems with one containing only
jω-axis dynamics, the other containing the non-jω axis dynamics.
G = rss(30,1,1);
[Gjw,G2] = modreal(G,1); % only one rigid body dynamics
G2.d = Gjw.d; Gjw.d = 0; % put DC gain of G into G2
subplot(211);sigma(Gjw);ylabel('Rigid Body')
subplot(212);sigma(G2);ylabel('Nonrigid Body')
3-12
Further model reduction can be done on G2 without any numerical difficulty.

After G2 is further reduced to Gred, the final approximation of the model is
simply Gjw+Gred.
This process of splitting jω-axis poles has been built in and automated in
all the model reduction routines (balancmr, schurmr, hankelmr, bstmr,
hankelsv) so that users need not worry about splitting the model.
The following single command creates a size 8 reduced-order model from

its original 30-state model:
G = rss(30,1,1);
[gr,info] = reduce(G); % choose a size of 8 at prompt
bode(G,'b-',gr,'r--')
3-13
Without specifying the size of the reduced-order model, a Hankel singular

value plot is shown below.
3-14
The default algorithm balancmr of reduce has done a great job of

approximating a 30-state model with just eight states. Again, the rigid body
dynamics are preserved for further controller design.
3-15
Reducing Large-Scale Models

For some really large size problems (states > 200), modreal turns out
to be the only way to start the model reduction process. Because of the
size and numerical properties associated with those large size, and low
damped dynamics, most Hankel based routines can fail to produce a good
reduced-order model.
modreal puts the large size dynamics into the modal form, then truncates the
dynamic model to an intermediate stage model with a comfortable size of 50
or so states. From this point on, those more sophisticated Hankel singular
value based routines can further reduce this intermediate stage model, in a
much more accurate fashion, to a smaller size for final controller design.
For a typical 240-state flexible spacecraft model in the spacecraft industry,

applying modreal and bstmr (or any other additive routines) in sequence can
reduce the original 240-state plant dynamics to a seven-state three-axis model
including rigid body dynamics. Any modern robust control design technique
mentioned in this toolbox can then be easily applied to this smaller size plant
for a controller design.
3-16
Normalized Coprime Factor Reduction
Normalized Coprime Factor Reduction

A special model reduction routine ncfmr produces a reduced-order model
by truncating a balanced coprime set of a given model. It can directly
simplify a modern controller with integrators to a smaller size by balanced
truncation of the normalized coprime factors. It does not need modreal for
pre-/postprocessing as the other routines do. However, any integrators in the
model will not be preserved.
rng(89,'twister');
K= rss(30,4,3);
[Kred,info2] = ncfmr(K);
Again, without specifying the size of the reduced-order model, any model
reduction routine presented here will plot a Hankel singular value bar chart
and prompt you for a reduced model size. In this case, enter 15.
Then, plot the singular values of the original and reduced-order models.
sigma(K,Kred)
legend('Original (30-state)','Kred (15-state)')
3-17
If integral control is important, previously mentioned methods (except ncfmr)

can nicely preserve the original integrator(s) in the model.
3-18
Bibliography
Bibliography
[1] Glover, K., “All Optimal Hankel Norm Approximation of Linear
Multivariable Systems, and Their L∝ - Error Bounds,” Int. J. Control, Vol. 39,
No. 6, 1984, pp. 1145-1193.
[2] Zhou, K., Doyle, J.C., and Glover, K., Robust and Optimal Control,
Englewood Cliffs, NJ, Prentice Hall, 1996.
[3] Safonov, M.G., and Chiang, R.Y., “A Schur Method for Balanced Model
Reduction,” IEEE Trans. on Automat. Contr., Vol. 34, No. 7, July 1989,
pp. 729-733.
[4] Safonov, M.G., Chiang, R.Y., and Limebeer, D.J.N., “Optimal Hankel
Model Reduction for Nonminimal Systems,” IEEE Trans. on Automat. Contr.,
Vol. 35, No. 4, April 1990, pp. 496-502.
[5] Safonov, M.G., and Chiang, R.Y., “Model Reduction for Robust Control:
A Schur Relative Error Method,” International J. of Adaptive Control and
Signal Processing, Vol. 2, 1988, pp. 259-272.
[6] Obinata, G., and Anderson, B.D.O., Model Reduction for Control System
Design, London, Springer-Verlag, 2001.
3-19
3-20
4
Robustness Analysis
• “Create Models of Uncertain Systems” on page 4-2

• “Robust Controller Design” on page 4-9
• “MIMO Robustness Analysis” on page 4-14
• “Summary of Robustness Analysis Tools” on page 4-27
4 Robustness Analysis
Create Models of Uncertain Systems

Dealing with and understanding the effects of uncertainty are important tasks
for the control engineer. Reducing the effects of some forms of uncertainty
(initial conditions, low-frequency disturbances) without catastrophically
increasing the effects of other dominant forms (sensor noise, model
uncertainty) is the primary job of the feedback control system.
Closed-loop stability is the way to deal with the (always present) uncertainty
in initial conditions or arbitrarily small disturbances.
High-gain feedback in low-frequency ranges is a way to deal with the effects

of unknown biases and disturbances acting on the process output. In this
case, you are forced to use roll-off filters in high-frequency ranges to deal with
high-frequency sensor noise in a feedback system.
Finally, notions such as gain and phase margins (and their generalizations)
help quantify the sensitivity of stability and performance in the face of model
uncertainty, which is the imprecise knowledge of how the control input
directly affects the feedback variables.
Robust Control Toolbox software has built-in features allowing you to specify
model uncertainty simply and naturally. The primary building blocks, called
uncertain elements (or uncertain Control Design Blocks) are uncertain real
parameters and uncertain linear, time-invariant objects. These can be used to
create coarse and simple or detailed and complex descriptions of the model
uncertainty present within your process models.
Once formulated, high-level system robustness tools can help you analyze the
potential degradation of stability and performance of the closed-loop system
brought on by the system model uncertainty.
Creating Uncertain Models of Dynamic Systems

The two dominant forms of model uncertainty are as follows:
• Uncertainty in parameters of the underlying differential equation models
4-2
• Frequency-domain uncertainty, which often quantifies model uncertainty

by describing absolute or relative uncertainty in the process’s frequency
response
Using these two basic building blocks, along with conventional system
creation commands (such as ss and tf), you can easily create uncertain
system models.
Creating Uncertain Parameters

An uncertain parameter has a name (used to identify it within an uncertain
system with many uncertain parameters) and a nominal value. Being
uncertain, it also has variability, described in one of the following ways:
• An additive deviation from the nominal

• A range about the nominal
• A percentage deviation from the nominal
Create a real parameter, with name ’bw’, nominal value 5, and a percentage
uncertainty of 10%.
bw = ureal('bw',5,'Percentage',10)
This creates a ureal object. View its properties using the get command.
Uncertain Real Parameter: Name bw, NominalValue 5, variability = [-10 10]%

get(bw)
Name: 'bw'
NominalValue: 5
Mode: 'Percentage'
Range: [4.5000 5.5000]
PlusMinus: [-0.5000 0.5000]
Percentage: [-10 10]
AutoSimplify: 'basic'
Note that the range of variation (Range property) and the additive deviation
from nominal (the PlusMinus property) are consistent with the Percentage
property value.
4-3
You can create state-space and transfer function models with uncertain
real coefficients using ureal objects. The result is an uncertain state-space
object, or uss. As an example, use the uncertain real parameter bw to model a
first-order system whose bandwidth is between 4.5 and 5.5 rad/s.
H = tf(1,[1/bw 1])
USS: 1 State, 1 Output, 1 Input, Continuous System
bw: real, nominal = 5, variability = [-10 10]%, 1 occurrence
Note that the result H is an uncertain system, called a uss object. The nominal
value of H is a state-space object. Verify that the pole is at –5.
pole(H.NominalValue)
ans =
-5
Next, use bode and step to examine the behavior of H.
bode(H,{1e-1 1e2});
4-4
step(H)
4-5
While there are variations in the bandwidth and time constant of H, the
high-frequency rolls off at 20 dB/decade regardless of the value of bw. You can
capture the more complicated uncertain behavior that typically occurs at high
frequencies using the ultidyn uncertain element, which is described next.
Quantifying Unmodeled Dynamics

An informal way to describe the difference between the model of a process and
the actual process behavior is in terms of bandwidth. It is common to hear
“The model is good out to 8 radians/second.” The precise meaning is not clear,
but it is reasonable to believe that for frequencies lower than, say, 5 rad/s,
the model is accurate, and for frequencies beyond, say, 30 rad/s, the model is
not necessarily representative of the process behavior. In the frequency range
between 5 and 30, the guaranteed accuracy of the model degrades.
The uncertain linear, time-invariant dynamics object ultidyn can be used

to model this type of knowledge. An ultidyn object represents an unknown
linear system whose only known attribute is a uniform magnitude bound
on its frequency response. When coupled with a nominal model and a
4-6
frequency-shaping filter, ultidyn objects can be used to capture uncertainty

associated with the model dynamics.
Suppose that the behavior of the system modeled by H significantly deviates

from its first-order behavior beyond 9 rad/s, for example, about 5% potential
relative error at low frequency, increasing to 1000% at high frequency where
H rolls off. In order to model frequency domain uncertainty as described above
using ultidyn objects, follow these steps:
1 Create the nominal system Gnom, using tf, ss, or zpk. Gnom itself might
already have parameter uncertainty. In this case Gnom is H, the first-order
system with an uncertain time constant.
2 Create a filter W, called the “weight,” whose magnitude represents the

relative uncertainty at each frequency. The utility makeweight is useful for
creating first-order weights with specific low- and high-frequency gains,
and specified gain crossover frequency.
3 Create an ultidyn object Delta with magnitude bound equal to 1.
The uncertain model G is formed by G = Gnom*(1+W*Delta).
If the magnitude of W represents an absolute (rather than relative)

uncertainty, use the formula G = Gnom + W*Delta instead.
The following commands carry out these steps:
Gnom = H;
W = makeweight(.05,9,10);
Delta = ultidyn('Delta',[1 1]);
G = Gnom*(1+W*Delta)
USS: 2 States, 1 Output, 1 Input, Continuous System
Delta: 1x1 LTI, max. gain = 1, 1 occurrence
bw: real, nominal = 5, variability = [-10 10]%, 1 occurrence
Note that the result G is also an uncertain system, with dependence on both
Delta and bw. You can use bode to make a Bode plot of 20 random samples of
G's behavior over the frequency range [0.1 100] rad/s.
bode(G,{1e-1 1e2},25)
4-7
In the next section, you design and compare two feedback controllers for G.
4-8
Robust Controller Design

In this tutorial, design a feedback controller for G, the uncertain model
created in “Create Models of Uncertain Systems” on page 4-2. The goals of
this design are the usual ones: good steady-state tracking and disturbance
rejection properties. Because the plant model is nominally a first-order lag,
choose a PI control architecture. Given the desired closed-loop damping ratio
ξ and natural frequency ωn, the design equations for KI and KP (based on the
nominal open-loop time constant of 0.2) are
n2 2n
KI = , KP = − 1.
5 5
Follow these steps to design the controller:
1 In order to study how the uncertain behavior of G affects the achievable

closed-loop bandwidth, design two controllers, both achieving ξ=0.707, with
different ωn: 3 and 7.5 respectively.
xi = 0.707;
wn = 3;
K1 = tf([(2*xi*wn/5-1) wn*wn/5],[1 0]);
wn = 7.5;
K2 = tf([(2*xi*wn/5-1) wn*wn/5],[1 0]);
Note that the nominal closed-loop bandwidth achieved by K2 is in a region

where G has significant model uncertainty. It will not be surprising if
the model variations lead to significant degradations in the closed-loop
performance.
2 Form the closed-loop systems using feedback.
T1 = feedback(G*K1,1);
T2 = feedback(G*K2,1);
3 Plot the step responses of 20 samples of each closed-loop system.
tfinal = 3;
stepplot(T1,'b',T2,'r',tfinal)
4-9
The step responses for T2 exhibit a faster rise time because K2 sets a higher
closed loop bandwidth. However, the model variations have a greater effect.
You can use robuststab to check the robustness of stability to the model
variations.
[stabmarg1,destabu1,report1] = robuststab(T1);
stabmarg1
stabmarg1 =
ubound: 4.0241
lbound: 4.0241
destabfreq: 3.4959
[stabmarg2,destabu2,report2] = robuststab(T2);
stabmarg2
stabmarg2 =
4-10
ubound: 1.2545
lbound: 1.2544
destabfreq: 10.5249
The stabmarg variable gives lower and upper bounds on the stability margin.
A stability margin greater than 1 means the system is stable for all values
of the modeled uncertainty. A stability margin less than 1 means there are
allowable values of the uncertain elements that make the system unstable.
The report variable briefly summarizes the analysis.
report1
report1 =
-- A destabilizing combination of 402% the modeled uncertainty exists, causing an instability at
3.5 rad/s.
report2
report2 =
-- A destabilizing combination of 125% the modeled uncertainty exists, causing an instability at
10.5 rad/s.
While both systems are stable for all variations, their performance is clearly
affected to different degrees. To determine how the uncertainty affects
closed-loop performance, you can use wcgain to compute the worst-case
effect of the uncertainty on the peak magnitude of the closed-loop sensitivity
(S=1/(1+GK)) function. This peak gain is typically correlated with the amount
of overshoot in a step response.
To do this, form the closed-loop sensitivity functions and call wcgain.
S1 = feedback(1,G*K1);
S2 = feedback(1,G*K2);
[maxgain1,wcu1] = wcgain(S1);
maxgain1
maxgain1 =
lbound: 1.8684
ubound: 1.9025
critfreq: 3.5152
4-11
[maxgain2,wcu2] = wcgain(S2);
maxgain2
maxgain2 =
lbound: 4.6031
ubound: 4.6671
critfreq: 11.0231
The maxgain variable gives lower and upper bounds on the worst-case peak
gain of the sensitivity transfer function, as well as the specific frequency
where the maximum gain occurs. The wcu variable contains specific values of
the uncertain elements that achieve this worst-case behavior.
You can use usubs to substitute these worst-case values for uncertain
elements, and compare the nominal and worst-case behavior. Use bodemag
and step to make the comparison.
bodemag(S1.NominalValue,'b',usubs(S1,wcu1),'b');
hold on
bodemag(S2.NominalValue,'r',usubs(S2,wcu2),'r');
hold off
4-12
Clearly, while K2 achieves better nominal sensitivity than K1, the nominal
closed-loop bandwidth extends too far into the frequency range where the
process uncertainty is very large. Hence the worst-case performance of K2 is
inferior to K1 for this particular uncertain model.
The next section explores these robustness analysis tools further on a

multiinput, multioutput system.
4-13
MIMO Robustness Analysis

The previous sections focused on simple uncertainty models of single-input
and single-output systems, predominantly from a transfer function
perspective. You can also create uncertain state-space models made up of
uncertain state-space matrices. Moreover, all the analysis tools covered thus
far can be applied to these systems as well.
Consider, for example, a two-input, two-output, two-state system whose

model has parametric uncertainty in the state-space matrices. First create
an uncertain parameter p. Using the parameter, make uncertain A and C
matrices. The B matrix happens to be not-uncertain, although you will add
frequency domain input uncertainty to the model in “Adding Independent
Input Uncertainty to Each Channel” on page 4-15.
p = ureal('p',10,'Percentage',10);
A = [0 p;-p 0];
B = eye(2);
C = [1 p;-p 1];
H = ss(A,B,C,[0 0;0 0]);
You can view the properties of the uncertain system H using the get command.
get(H)
a: [2x2 umat]
b: [2x2 double]
c: [2x2 umat]
d: [2x2 double]
e: []
StateName: {2x1 cell}
StateUnit: {2x1 cell}
NominalValue: [2x2 ss]
Uncertainty: [1x1 struct]
InternalDelay: [0x1 double]
InputDelay: [2x1 double]
OutputDelay: [2x1 double]
Ts: 0
TimeUnit: 'seconds'
InputName: {2x1 cell}
4-14
InputUnit: {2x1 cell}

InputGroup: [1x1 struct]
OutputName: {2x1 cell}
OutputUnit: {2x1 cell}
OutputGroup: [1x1 struct]
Name: ''
Notes: {}
UserData: []
The properties a, b, c, d, and StateName behave in exactly the same manner

as ss objects. The properties InputName, OutputName, InputGroup, and
OutputGroup behave in exactly the same manner as all the system objects
(ss, zpk, tf, and frd). The NominalValue is an ss object.
Adding Independent Input Uncertainty

to Each Channel
The model for H does not include actuator dynamics. Said differently, the
actuator models are unity-gain for all frequencies.
Nevertheless, the behavior of the actuator for channel 1 is modestly uncertain

(say 10%) at low frequencies, and the high-frequency behavior beyond 20
rad/s is not accurately modeled. Similar statements hold for the actuator in
channel 2, with larger modest uncertainty at low frequency (say 20%) but
accuracy out to 45 rad/s.
Use ultidyn objects Delta1 and Delta2 along with shaping filters W1 and W2
to add this form of frequency domain uncertainty to the model.
W1 = makeweight(.1,20,50);
W2 = makeweight(.2,45,50);
Delta1 = ultidyn('Delta1',[1 1]);
Delta2 = ultidyn('Delta2',[1 1]);
G = H*blkdiag(1+W1*Delta1,1+W2*Delta2)
G =
Uncertain continuous-time state-space model with 2 outputs, 2 inputs, 4 s

The model uncertainty consists of the following blocks:
Delta1: Uncertain 1x1 LTI, peak gain = 1, 1 occurrences
4-15
Delta2: Uncertain 1x1 LTI, peak gain = 1, 1 occurrences

p: Uncertain real, nominal = 10, variability = [-10,10]%, 2 occurrences
Type "G.NominalValue" to see the nominal value, "get(G)" to see all propert
Note that G is a two-input, two-output uncertain system, with dependence on

three uncertain elements, Delta1, Delta2, and p. It has four states, two from
H and one each from the shaping filters W1 and W2, which are embedded in G.
You can plot a 2-second step response of several samples of G. The 10%
uncertainty in the natural frequency is obvious.
stepplot(G,2)
4-16
You can create a Bode plot of samples of G. The high-frequency uncertainty

in the model is also obvious. For clarity, start the Bode plot beyond the
resonance.
bodeplot(G,{13 100})
4-17
4-18
Closed-Loop Robustness Analysis

Load the controller and verify that it is two-input and two-output.
load mimoKexample
size(K)
State-space model with 2 outputs, 2 inputs, and 9 states.
You can use the command loopsens to form all the standard plant/controller
feedback configurations, including sensitivity and complementary sensitivity
at both the input and output. Because G is uncertain, all the closed-loop
systems are uncertain as well.
F = loopsens(G,K)
F =
Si: [2x2 uss]

Ti: [2x2 uss]
Li: [2x2 uss]
So: [2x2 uss]
To: [2x2 uss]
Lo: [2x2 uss]
PSi: [2x2 uss]
CSo: [2x2 uss]
Poles: [13x1 double]
Stable: 1
F is a structure with many fields. The poles of the nominal closed-loop system
are in F.Poles, and F.Stable is 1 if the nominal closed-loop system is stable.
In the remaining 10 fields, S stands for sensitivity, T for complementary
sensitivity, and L for open-loop gain. The suffixes i and o refer to the input and
output of the plant (G). Finally, P and C refer to the “plant” and “controller.”
Hence Ti is mathematically the same as
K(I + GK)–1G
while Lo is G*K, and CSo is mathematically the same as
4-19
K(I + GK)–1
You can examine the transmission of disturbances at the plant input to the
plant output using bodemag on F.PSi. Graph some samples along with the
nominal.
bodemag(F.PSi.NominalValue,'r+',F.PSi,'b-',{1e-1 100})
4-20
Nominal Stability Margins

You can use loopmargin to investigate loop-at-a-time gain and phase margins,
loop-at-a-time disk margins, and simultaneous multivariable margins. They
are computed for the nominal system and do not reflect the uncertainty
models within G.
Explore the simultaneous margins individually at the plant input,

individually at the plant output, and simultaneously at both input and output.
[I,DI,SimI,O,DO,SimO,Sim] = loopmargin(G,K);
The third output argument is the simultaneous gain and phase variations
allowed in all input channels to the plant.
SimI
SimI =
GainMargin: [0.1180 8.4769]

PhaseMargin: [-76.5441 76.5441]
Frequency: 6.2287
This information implies that the gain at the plant input can vary in both
channels independently by factors between (approximately) 1/8 and 8, as well
as phase variations up to 76 degrees.
The sixth output argument is the simultaneous gain and phase variations
allowed in all output channels to the plant.
SimO
SimO =
GainMargin: [0.1193 8.3836]

PhaseMargin: [-76.3957 76.3957]
Frequency: 18.3522
Note that the simultaneous margins at the plant output are similar to those
at the input. This is not always the case in multiloop feedback systems.
4-21
The last output argument is the simultaneous gain and phase variations
allowed in all input and output channels to the plant. As expected, when
you consider all such variations simultaneously, the margins are somewhat
smaller than those at the input or output alone.
Sim
Sim =
GainMargin: [0.5671 1.7635]

PhaseMargin: [-30.8882 30.8882]
Frequency: 18.3522
Nevertheless, these numbers indicate a generally robust closed-loop system,

able to tolerate significant gain (more than +/-50% in each channel) and 30
degree phase variations simultaneously in all input and output channels
of the plant.
Robustness of Stability Model Uncertainty

With loopmargin, you determined various margins of the nominal, multiloop
system. These margins are computed only for the nominal system, and do
not reflect the uncertainty explicitly modeled by the ureal and ultidyn
objects. When you work with detailed, complex uncertain system models,
the conventional margins computed by loopmargin might not always be
indicative of the actual stability margins associated with the uncertain
elements. You can use robuststab to check the stability margin of the system
to these specific modeled variations.
In this example, use robuststab to compute the stability margin of the

closed-loop system represented by Delta1, Delta2, and p.
Use any of the closed-loop systems within F = loopsens(G,K). All of them,

F.Si, F.To, etc., have the same internal dynamics, and hence the stability
properties are the same.
[stabmarg,desgtabu,report] = robuststab(F.So);
stabmarg
stabmarg =
4-22
ubound: 2.2175
lbound: 2.2175
destabfreq: 13.7576
report
report =

-- A destabilizing combination of 222% the modeled uncertainty exists,causing an instability at
13.8 rad/s.
This analysis confirms what the loopmargin analysis suggested. The

closed-loop system is quite robust, in terms of stability, to the variations
modeled by the uncertain parameters Delta1, Delta2, and p. In fact, the
system can tolerate more than twice the modeled uncertainty without losing
closed-loop stability.
The next section studies the effects of these variations on the closed-loop
output sensitivity function.
Worst-Case Gain Analysis

You can plot the Bode magnitude of the nominal output sensitivity function.
It clearly shows decent disturbance rejection in all channels at low frequency.
bodemag(F.So.NominalValue,{1e-1 100})
4-23
You can compute the peak value of the maximum singular value of the
frequency response matrix using norm.
[PeakNom,freq] = norm(F.So.NominalValue,'inf')
PeakNom =
1.1317
freq =
7.0483
The peak is about 1.13, occurring at a frequency of 36 rad/s.
4-24
What is the maximum output sensitivity gain that is achieved when the
uncertain elements Delta1, Delta2, and p vary over their ranges? You can
use wcgain to answer this.
[maxgain,wcu] = wcgain(F.So);
maxgain
maxgain =
lbound: 2.1017
ubound: 2.1835
critfreq: 8.5546
The analysis indicates that the worst-case gain is somewhere between 2.1 and
2.2. The frequency where the peak is achieved is about 8.5.
You can replace the values for Delta1, Delta2, and p that achieve the gain
of 2.1, using usubs. Make the substitution in the output complementary
sensitivity, and do a step response.
step(F.To.NominalValue,usubs(F.To,wcu),5)
4-25
The perturbed response, which is the worst combination of uncertain values

in terms of output sensitivity amplification, does not show significant
degradation of the command response. The settling time is increased by about
50%, from 2 to 4, and the off-diagonal coupling is increased by about a factor
of about 2, but is still quite small.
4-26
Summary of Robustness Analysis Tools
Summary of Robustness Analysis Tools

Function Description
ureal Create uncertain real parameter.
ultidyn Create uncertain, linear, time-invariant dynamics.
uss Create uncertain state-space object from uncertain
state-space matrices.
ufrd Create uncertain frequency response object.
loopsens Compute all relevant open and closed-loop
quantities for a MIMO feedback connection.
loopmargin Compute loop-at-a-time as well as MIMO gain and
phase margins for a multiloop system, including
the simultaneous gain/phase margins.
robustperf Robustness performance of uncertain systems.
robuststab Compute the robust stability margin of a nominally
stable uncertain system.
wcgain Compute the worst-case gain of a nominally stable
uncertain system.
wcmargin Compute worst-case (over uncertainty)
loop-at-a-time disk-based gain and phase margins.
wcsens Compute worst-case (over uncertainty) sensitivity
of plant-controller feedback loop.
4-27
4-28
5
H-Infinity and Mu
Synthesis
• “Interpretation of H-Infinity Norm” on page 5-2

• “H-Infinity Performance” on page 5-9
• “Active Suspension Control Design” on page 5-17
5 H-Infinity and Mu Synthesis
Interpretation of H-Infinity Norm
Norms of Signals and Systems

There are several ways of defining norms of a scalar signal e(t) in the
time domain. We will often use the 2-norm, (L2-norm), for mathematical
convenience, which is defined as
1
⎛ ∞ ⎞2
:= ⎜ ∫ e ( t ) dt ⎟ .
2
e2
⎝ −∞ ⎠
If this integral is finite, then the signal e is square integrable, denoted as e

L2. For vector-valued signals
⎡ e1 ( t ) ⎤
⎢ ⎥
e (t)
e (t) = ⎢ 2 ⎥ ,
⎢  ⎥
⎢ ⎥
⎢⎣ en ( t ) ⎥⎦
the 2-norm is defined as
1
⎛ ∞ 2 ⎞2
e2 := ⎜ ∫ e ( t ) dt ⎟
⎝ −∞ 2 ⎠
1
⎛ ∞ ⎞2
= ⎜ ∫ eT ( t ) e ( t ) dt ⎟ .
⎝ −∞ ⎠
In µ-tools the dynamic systems we deal with are exclusively linear, with
state-space model
⎡ x ⎤ ⎡ A B⎤ ⎡ x ⎤
⎢ e ⎥ = ⎢C D ⎥⎦ ⎢⎣ d ⎥⎦
,
⎣ ⎦ ⎣
or, in the transfer function form,
5-2
e(s) = T(s)d(s), T(s):= C(sI – A)–1B + D
Two mathematically convenient measures of the transfer matrix T(s) in the

frequency domain are the matrix H2 and H∞ norms,
1
⎡ 1 ∞ 2 ⎤2
T ( j ) d ⎥
⎣ 2 ∫−∞
T 2 := ⎢
F ⎦
T ∞ := max  ⎡⎣T ( j ) ⎤⎦ ,
∈R
where the Frobenius norm (see the MATLAB norm command) of a complex
matrix M is
M F
:= Trace M * M .( )
Both of these transfer function norms have input/output time-domain
interpretations. If, starting from initial condition x(0) = 0, two signals d and
e are related by
⎡ x ⎤ ⎡ A B⎤ ⎡ x ⎤
⎢ e ⎥ = ⎢C D ⎥⎦ ⎢⎣ d ⎥⎦
,
⎣ ⎦ ⎣
then
• For d, a unit intensity, white noise process, the steady-state variance of e

is T 2.
• The L2 (or RMS) gain from d → e,
e2
max
d ≠0 d 2
is equal to T ∞. This is discussed in greater detail in the next section.
Using Weighted Norms to Characterize Performance

Any performance criterion must also account for
5-3
• Relative magnitude of outside influences

• Frequency dependence of signals
• Relative importance of the magnitudes of regulated variables
So, if the performance objective is in the form of a matrix norm, it should

actually be a weighted norm
WLTWR
where the weighting function matrices WL and WR are frequency dependent,

to account for bandwidth constraints and spectral content of exogenous
signals. The most natural (mathematical) manner to characterize acceptable
performance is in terms of the MIMO · ∞ (H∞) norm. For this reason, this
section now discusses some interpretations of the H∞ norm.
Unweighted MIMO System
Suppose T is a MIMO stable linear system, with transfer function matrix T(s).
For a given driving signal d ( t ) , define e as the output, as shown below.
Note that it is more traditional to write the diagram in Unweighted MIMO

System: Vectors from Left to Right on page 5-4 with the arrows going from
left to right as in Weighted MIMO System on page 5-6.
Unweighted MIMO System: Vectors from Left to Right
The two diagrams shown above represent the exact same system. We prefer to
write these block diagrams with the arrows going right to left to be consistent
with matrix and operator composition.
Assume that the dimensions of T are ne × nd. Let β > 0 be defined as
5-4
 := T ∞
:= max  ⎡⎣T ( j ) ⎤⎦ .
∈R
Now consider a response, starting from initial condition equal to 0. In that

case, Parseval’s theorem gives that
1
⎡ ∞ T ⎤2
⎢⎣ ∫0 e ( t ) e ( t ) dt ⎥⎦

e 2
= ≤ .
d 1
⎡ ∞ T ⎤2
⎢⎣ ∫0 d ( t ) d ( t ) dt ⎥⎦
2 
Moreover, there are specific disturbances d that result in the ratio e 2 d

2
arbitrarily close to β. Because of this, T ∞ is referred to as the L2 (or RMS)
gain of the system.
As you would expect, a sinusoidal, steady-state interpretation of T ∞ is also

possible: For any frequency  ∈ R , any vector of amplitudes a ∈ Rnd , and
any vector of phases  ∈ Rnd , with a 2 ≤ 1, define a time signal
⎡ a sin  t +  ⎤
⎢ 1 ( 1) ⎥
d (t) = ⎢
  ⎥.
⎢ ⎥
⎣ (
⎢ and sin  t + nd ) ⎥
⎦
Applying this input to the system T results in a steady-state response ess of

the form
⎡ b sin  t +  ⎤
⎢ 1 ( 1) ⎥
ess ( t ) = ⎢  ⎥.
⎢ ⎥
⎣ (
⎢bne sin  t + ne ) ⎥
⎦
5-5
The vector b ∈ Rne will satisfy b 2 ≤ β. Moreover, β, as defined in Weighted

MIMO System on page 5-6, is the smallest number such that this is true
for every a 2 ≤ 1,  , and ϕ.
Note that in this interpretation, the vectors of the sinusoidal magnitude

responses are unweighted, and measured in Euclidean norm. If realistic
multivariable performance objectives are to be represented by a single
MIMO · ∞ objective on a closed-loop transfer function, additional scalings
are necessary. Because many different objectives are being lumped into one
matrix and the associated cost is the norm of the matrix, it is important to
use frequency-dependent weighting functions, so that different requirements
can be meaningfully combined into a single cost function. Diagonal weights
are most easily interpreted.
Consider the diagram of Weighted MIMO System on page 5-6, along with
Unweighted MIMO System: Vectors from Left to Right on page 5-4.
Assume that WL and WR are diagonal, stable transfer function matrices, with
diagonal entries denoted Li and Ri.
⎡ L1 0  0 ⎤ ⎡ R1 0  0 ⎤
⎢0 L2  0 ⎥⎥ ⎢0 R2  0 ⎥⎥
WL = ⎢ , WR = ⎢ .
⎢     ⎥ ⎢     ⎥
⎢ ⎥ ⎢ ⎥
⎢⎣ 0 0  Lne ⎦⎥ ⎢⎣ 0 0  Rnd ⎦⎥
Weighted MIMO System
Bounds on the quantity WLTWR ∞ will imply bounds about the sinusoidal
( )
steady-state behavior of the signals d and e = Td in the diagram of
Unweighted MIMO System: Vectors from Left to Right on page 5-4.
Specifically, for sinusoidal signal d , the steady-state relationship between
5-6
( )
e = Td , d and WLTWR ∞ is as follows. The steady-state solution ess ,
denoted as
⎡ e sin  t +  ⎤
⎢ 1 ( 1) ⎥
ess ( t ) = ⎢  ⎥
⎢ ⎥
⎣ (
⎢ ene sin  t + nd ) ⎥
⎦ (5-1)
satisfies
ne
2
∑ WL ( j ) ei
i
≤1
i=1
for all sinusoidal input signals d of the form
⎡ d sin  t +  ⎤
⎢ 1 ( 1) ⎥
d (t) = ⎢
  ⎥
⎢ ⎥
⎣ (
⎢ dne sin  t + nd ) ⎥
⎦ (5-2)
satisfying
2
nd di
∑ 2
≤1
i=1 WRi ( j )
if and only if WLTWR ∞ ≤ 1.
This approximately (very approximately — the next statement is not actually

correct) implies that WLTWR ∞ ≤ 1 if and only if for every fixed frequency  ,
and all sinusoidal disturbances d of the form Equation 5-2 satisfying
di ≤ WRi ( j )
5-7
the steady-state error components will satisfy
1
ei ≤ .
WLi ( j )
This shows how one could pick performance weights to reflect the desired
frequency-dependent performance objective. Use WR to represent the relative
magnitude of sinusoids disturbances that might be present, and use 1/WL to
represent the desired upper bound on the subsequent errors that are produced.
Remember, however, that the weighted H∞ norm does not actually

give element-by-element bounds on the components of e based on
element-by-element bounds on the components of d . The precise bound it
gives is in terms of Euclidean norms of the components of e and d (weighted
appropriately by WL(j  ) and WR(j  )).
5-8
H-Infinity Performance
Performance as Generalized Disturbance Rejection

The modern approach to characterizing closed-loop performance objectives
is to measure the size of certain closed-loop transfer function matrices using
various matrix norms. Matrix norms provide a measure of how large output
signals can get for certain classes of input signals. Optimizing these types
of performance objectives over the set of stabilizing controllers is the main
thrust of recent optimal control theory, such as L1, H2, H∞, and optimal
control. Hence, it is important to understand how many types of control
objectives can be posed as a minimization of closed-loop transfer functions.
Consider a tracking problem, with disturbance rejection, measurement noise,

and control input signal limitations, as shown in Generalized and Weighted
Performance Block Diagram on page 5-11. K is some controller to be designed
and G is the system you want to control.
Typical Closed-Loop Performance Objective
A reasonable, though not precise, design objective would be to design K

to keep tracking errors and control input signal small for all reasonable
reference commands, sensor noises, and external force disturbances.
Hence, a natural performance objective is the closed-loop gain from

exogenous influences (reference commands, sensor noise, and external force
disturbances) to regulated variables (tracking errors and control input signal).
Specifically, let T denote the closed-loop mapping from the outside influences
to the regulated variables:
5-9
You can assess performance by measuring the gain from outside influences
to regulated variables. In other words, good performance is associated with
T being small. Because the closed-loop system is a multiinput, multioutput
(MIMO) dynamic system, there are two different aspects to the gain of T:
• Spatial (vector disturbances and vector errors)

• Temporal (dynamic relationship between input/output signals)
Hence the performance criterion must account for
• Relative magnitude of outside influences

• Frequency dependence of signals
• Relative importance of the magnitudes of regulated variables
So if the performance objective is in the form of a matrix norm, it should

actually be a weighted norm
WLTWR
where the weighting function matrices WL and WR are frequency dependent, to

account for bandwidth constraints and spectral content of exogenous signals.
A natural (mathematical) manner to characterize acceptable performance
is in terms of the MIMO · ∞ (H∞) norm. See “Interpretation of H-Infinity
Norm” on page 5-2 for an interpretation of the H∞ norm and signals.
Interconnection with Typical MIMO Performance Objectives

The closed-loop performance objectives are formulated as weighted closed-loop
transfer functions that are to be made small through feedback. A generic
example, which includes many relevant terms, is shown in block diagram
form in Generalized and Weighted Performance Block Diagram on page 5-11.
In the diagram, G denotes the plant model and K is the feedback controller.
5-10
Wmodel
e1
d2
Wact ~ ~
d2 - e2
~ Wdist Wperf1 e2
d1 d1 ~ ~
Wcmd e1 e3
K G Wperf2 e3
Hsens
~
d3 d3
Wsnois
Generalized and Weighted Performance Block Diagram
The blocks in this figure might be scalar (SISO) and/or multivariable (MIMO),
depending on the specific example. The mathematical objective of H∞ control
is to make the closed-loop MIMO transfer function Ted satisfy Ted ∞ < 1. The
weighting functions are used to scale the input/output transfer functions such
that when T < 1, the relationship between d and e is suitable.
ed ∞
Performance requirements on the closed-loop system are transformed into

the H∞ framework with the help of weighting or scaling functions. Weights
are selected to account for the relative magnitude of signals, their frequency
dependence, and their relative importance. This is captured in the figure
above, where the weights or scalings [Wcmd, Wdist,Wsnois] are used to transform
and scale the normalized input signals [d1,d2,d3] into physical units defined
as [d1, d2, d3]. Similarly weights or scalings [Wact, Wperf1,Wperf2] transform
and scale physical units into normalized output signals [e1, e2, e3]. An
interpretation of the signals, weighting functions, and models follows.
5-11
Signal Meaning
d1 Normalized reference command
Typical reference command in physical units
d1
d2 Normalized exogenous disturbances
Typical exogenous disturbances in physical units
d2
d3 Normalized sensor noise
Typical sensor noise in physical units
d3
e1 Weighted control signals
e1 Actual control signals in physical units

e2 Weighted tracking errors
e2 Actual tracking errors in physical units

e3 Weighted plant errors
e3 Actual plant errors in physical units
Wcmd
Wcmd is included in H∞ control problems that require tracking of a reference

command. Wcmd shapes the normalized reference command signals
(magnitude and frequency) into the actual (or typical) reference signals that
you expect to occur. It describes the magnitude and the frequency dependence
of the reference commands generated by the normalized reference signal.
Normally Wcmd is flat at low frequency and rolls off at high frequency. For
example, in a flight control problem, fighter pilots generate stick input
reference commands up to a bandwidth of about 2 Hz. Suppose that the stick
has a maximum travel of three inches. Pilot commands could be modeled as
normalized signals passed through a first-order filter:
3
Wcmd = .
1
s+1
2 ⋅ 2
5-12
Wmodel
Wmodel represents a desired ideal model for the closed-looped system and is
often included in problem formulations with tracking requirements. Inclusion
of an ideal model for tracking is often called a model matching problem, i.e.,
the objective of the closed-loop system is to match the defined model. For
good command tracking response, you might want the closed-loop system to
respond like a well-damped second-order system. The ideal model would
then be
2
Wmodel =
s2 + 2 +  2
for specific desired natural frequency ω and desired damping ratio ζ. Unit
conversions might be necessary to ensure exact correlation between the ideal
model and the closed-loop system. In the fighter pilot example, suppose that
roll-rate is being commanded and 10º/second response is desired for each inch
of stick motion. Then, in these units, the appropriate model is:
2
Wmodel = 10 .
s2 + 2 +  2
Wdist
Wdist shapes the frequency content and magnitude of the exogenous

disturbances affecting the plant. For example, consider an electron microscope
as the plant. The dominant performance objective is to mechanically isolate
the microscope from outside mechanical disturbances, such as ground
excitations, sound (pressure) waves, and air currents. You can capture the
spectrum and relative magnitudes of these disturbances with the transfer
function weighting matrix Wdist.
Wperf1
Wperf1 weights the difference between the response of the closed-loop system
and the ideal model W model. Often you might want accurate matching of the
ideal model at low frequency and require less accurate matching at higher
frequency, in which case Wperf1 is flat at low frequency, rolls off at first or
5-13
second order, and flattens out at a small, nonzero value at high frequency.
The inverse of the weight is related to the allowable size of tracking errors,
when dealing with the reference commands and disturbances described by
Wcmd and Wdist.
Wperf2
Wperf2 penalizes variables internal to the process G, such as actuator states

that are internal to G or other variables that are not part of the tracking
objective.
Wact
Wact is used to shape the penalty on control signal use. Wact is a frequency
varying weighting function used to penalize limits on the deflection/position,
deflection rate/velocity, etc., response of the control signals, when dealing
with the tracking and disturbance rejection objectives defined above. Each
control signal is usually penalized independently.
Wsnois
Wsnois represents frequency domain models of sensor noise. Each sensor

measurement feedback to the controller has some noise, which is often
higher in one frequency range than another. The Wsnois weight tries to
capture this information, derived from laboratory experiments or based on
manufacturer measurements, in the control problem. For example, medium
grade accelerometers have substantial noise at low frequency and high
frequency. Therefore the corresponding Wsnois weight would be larger at low
and high frequency and have a smaller magnitude in the mid-frequency
range. Displacement or rotation measurement is often quite accurate at low
frequency and in steady state, but responds poorly as frequency increases.
The weighting function for this sensor would be small at low frequency,
gradually increase in magnitude as a first- or second-order system, and level
out at high frequency.
Hsens
Hsens represents a model of the sensor dynamics or an external antialiasing

filter. The transfer functions used to describe Hsens are based on physical
5-14
characteristics of the individual components. These models might also be

lumped into the plant model G.
This generic block diagram has tremendous flexibility and many control
performance objectives can be formulated in the H∞ framework using this
block diagram description.
Robustness in the H-Infinity Framework

Performance and robustness tradeoffs in control design were discussed in the
context of multivariable loop shaping in “Tradeoff Between Performance and
Robustness” on page 2-2. In the H∞ control design framework, you can include
robustness objectives as additional disturbance to error transfer functions —
disturbances to be kept small. Consider the following figure of a closed-loop
feedback system with additive and multiplicative uncertainty models.
The transfer function matrices are defined as:
−1
TF ( s ) z →w = TI ( s ) = KG ( I + GK )
1 1
−1
TF ( s ) z = KSO ( s ) = K ( I + GK )
2 →w2
where TI(s) denotes the input complementary sensitivity function and SO(s)
denotes the output sensitivity function. Bounds on the size of the transfer
function matrices from z1 to w1 and z2 to w2 ensure that the closed-loop
system is robust to multiplicative uncertainty, ΔM(s), at the plant input, and
additive uncertainty, ΔA(s), around the plant G(s). In the H∞ control problem
formulation, the robustness objectives enter the synthesis procedure as
additional input/output signals to be kept small. The interconnection with the
uncertainty blocks removed follows.
5-15
The H∞ control robustness objective is now in the same format as the

performance objectives, that is, to minimize the H∞ norm of the transfer
matrix from z, [z1,z2], to w, [w1,w2].
Weighting or scaling matrices are often introduced to shape the frequency and
magnitude content of the sensitivity and complementary sensitivity transfer
function matrices. Let WM correspond to the multiplicative uncertainty and
WA correspond to the additive uncertainty model. ΔM(s) and ΔA(s) are assumed
to be a norm bounded by 1, i.e., |ΔM(s)|<1 and |ΔA(s)|<1. Hence as a function
of frequency, |WM(jω)| and |WA(jω)| are the respective sizes of the largest
anticipated additive and multiplicative plant perturbations.
The multiplicative weighting or scaling WM represents a percentage error in

the model and is often small in magnitude at low frequency, between 0.05 and
0.20 (5% to 20% modeling error), and growing larger in magnitude at high
frequency, 2 to 5 ((200% to 500% modeling error). The weight will transition
by crossing a magnitude value of 1, which corresponds to 100% uncertainty
in the model, at a frequency at least twice the bandwidth of the closed-loop
system. A typical multiplicative weight is
1
s+1
WM = 0.10 5 .
1
s+1
200
By contrast, the additive weight or scaling WA represents an absolute error

that is often small at low frequency and large in magnitude at high frequency.
The magnitude of this weight depends directly on the magnitude of the plant
model, G(s).
5-16
Active Suspension Control Design

This example shows how to use robust control techniques to design an active
suspension system for a quarter car body and wheel assembly model. In this
example, you use H∞ design techniques to design a controller for a nominal
quarter-car model. Then, you use μ synthesis to design a robust controller
that accounts for uncertainty in the model.
Conventional passive suspensions use a spring and damper between the

car body and wheel assembly. The spring-damper characteristics are
adjusted to emphasize one of several conflicting objectives such as passenger
comfort, road holding, and suspension deflection. Active suspensions use
a feedback-controller hydraulic actuator between the chassis and wheel
assembly, which allows the designer to better balance these objectives.
Quarter-Car Suspension Model
This example uses the quarter-car model of the following illustration to design
active suspension control laws.
mb xb
bs
ks fs
xw
mw
kt
r
The mass, mb, represents the car chassis (body) and the mass, mw, represents
the wheel assembly. The spring, ks, and damper, bs, represent the passive
5-17
spring and shock absorber placed between the car body and the wheel
assembly. The spring, kt, models the compressibility of the pneumatic tire.
The variables xb, xw, and r are the car body travel, the wheel travel, and the
road disturbance, respectively. The force, fs, which is applied between the
body and wheel assembly, is controlled by feedback. This force represents the
active component of the suspension system.
The following state-space equations describe the quarter-car dynamics.
x1  x2 ,
1
x 2    ks  x1  x3   bs  x2  x4   fs  ,
mb 
x 3  x4 ,
1
x 4   ks  x1  x3   bs  x2  x4   kt  x3  r   fs  .
mw 
The state variables in the system are defined as x1 : xb , x2 : x b , x3 : xw ,

and x4 : x w .
Define the physical parameters of the system.
mb = 300; % kg
mw = 60; % kg
bs = 1000; % N/m/s
ks = 16000 ; % N/m
kt = 190000; % N/m
Use these equations and parameter values to construct a state-space model,

qcar, of the quarter-car suspension system.
A = [ 0 1 0 0; [-ks -bs ks bs]/mb ; ...

0 0 0 1; [ks bs -ks-kt -bs]/mw];
B = [0 0; 0 10000/mb ; 0 0; [kt -10000]/mw];
C = [1 0 0 0; 1 0 -1 0; A(2,:)];
D = [0 0; 0 0; B(2,:)];
qcar = ss(A,B,C,D);
qcar.StateName = {'body travel xb (m)';'body vel (m/s)';...
5-18
'wheel travel xw (m)';'wheel vel (m/s)'};

qcar.InputName = {'r';'fs'};
qcar.OutputName = {'xb';'sd';'ab'};
The model inputs are the road disturbance, r, and actuator force, fs. The
model outputs are the car body travel, xb , suspension deflection sd  xb  xw ,
and car body acceleration ab  
xs .
The transfer function from actuator to body travel and acceleration has an
imaginary-axis zero. Examine the zero of this transfer function.
tzero(qcar({'xb','ab'},'fs'))
ans =
-0.0000 +56.2731i
-0.0000 -56.2731i
The natural frequency of this zero, 56.27 rad/s, is called the tire-hop frequency.
The transfer function from the actuator to suspension deflection also has an
imaginary-axis zero. Examine this zero.
zero(qcar('sd','fs'))
ans =
0.0000 +22.9734i
0.0000 -22.9734i
The natural frequency of this zero, 22.97 rad/s, is called the rattlespace
frequency.
Plot the frequency response of the quarter-car model from inputs, (r,fs), to
outputs, (ab,sd). Both zeros are visible on the Bode plot.
bodemag(qcar({'ab','sd'},'r'),'b',qcar({'ab','sd'},'fs'),'r',{1 100});
legend('Road disturbance (r)','Actuator force (fs)','location','SouthWest')
title(['Gain from road dist (r) and actuator force (fs) '...
'to body accel (ab) and suspension travel (sd)'])
5-19
Road disturbances influence the motion of the car and suspension:
• Small body acceleration influences passenger comfort.

• Small suspension travel contributes to good road handling. Further, limits
on the actuator displacement constrain the allowable travel.
Because of the imaginary axis zeros, feedback control cannot improve the
response from road disturbance (r) to body acceleration (ab) at the tire-hop
frequency. Similarly, feedback control cannot improve the response from r to
suspension deflection (sd) at the rattlespace frequency. Moreover, there is an
inherent trade-off between passenger comfort and suspension deflection. Any
reduction of body travel at low frequency results in an increase of suspension
deflection. This trade-off arises because of the relationship xw = xb – sd and the
fact that xw roughly follows r at low frequency (less than 5 rad/s).
Hydraulic Actuator Model
5-20
The hydraulic actuator used for active suspension control is connected

between the body mass, mb, and the wheel assembly mass, mw. The nominal
actuator dynamics can be represented by the first-order transfer function:
1
ActNom  s   .
1
s1
60
The maximum displacement is 0.05 m.
The nominal actuator model approximates the physical actuator dynamics.

You can model variations between the actuator model and the physical
device as a family of actuator models. You can also use this approach to
model variations between the passive quarter-car model and actual vehicle
dynamics. The resulting family of models comprises a nominal model with a
frequency-dependent amount of uncertainty.
Create an uncertain model that represents this family of models.
ActNom = tf(1,[1/60 1]);

Wunc = makeweight(0.40,15,3);
unc = ultidyn('unc',[1 1],'SampleStateDim',5);
Act = ActNom*(1 + Wunc*unc);
Act.InputName = 'u';
Act.OutputName = 'fs';
At low frequency, below 3 rad/s, the model can vary up to 40% from its nominal
value. Around 3 rad/s, the percentage variation starts to increase. The
uncertainty crosses 100% at 15 rad/s, and reaches 2000% at approximately
1000 rad/sec. The weighting function, Wunc, reflects this profile and is used to
modulate the amount of uncertainty as a function of frequency. The result
Act is an uncertain state-space model of the actuator.
Examine the uncertain actuator model by plotting the frequency response of

20 randomly sampled models from Act.
bode(Act,'b',Act.NominalValue,'r+',logspace(-1,3,120))
title('Nominal and 20 random actuator models')
5-21
The plus (+) marker denotes the nominal actuator model. The blue solid lines
represent the randomly sampled models.
Design Objectives for H-Infinity Synthesis
To use H∞ synthesis algorithms, you must express your design objectives as

a single cost function to be minimized. For the quarter-car model, the main
control objectives are formulated in terms of passenger comfort and road
handling. These objectives relate to body acceleration, ab, and suspension
travel, sd. Other factors that influence the control design include:
• Characteristics of the road disturbance

• Quality of the sensor measurements for feedback
• Limits on the available control force
Use weights to model external disturbances and quantify the design

objectives, as shown in the following diagram.
5-22
xb
r xb
d1 Wroad
1/4 Car sd
Model W sd e2
fs (qcar)
u Act ab
W ab e3
fs
Wact e1 y1 W d2 d2
y2 W d3 d3
The feedback controller uses the measurements y1 and y2 of the suspension

travel, sd, and body acceleration, ab, to compute the control signal, u. This
control signal drives the hydraulic actuator. There are three external sources
of disturbance:
• The road disturbance, r, which is modeled as a normalized signal, d1, which

is shaped by a weighting function Wroad.
• Sensor noise on both measurements. This noise is modeled as normalized
signals, d2 and d3, which are shaped by weighting functions Wd2 and Wd3.
You can reinterpret the control objectives as a disturbance rejection goal.

The goal is to minimize the impact of the disturbances, d1, d2, and d3, on a
weighted combination of suspension travel (sd), body acceleration (ab), and
control effort (u). You can consider the H∞ norm (peak gain) as the measure
of the effect of the disturbances. Then, you can meet the requirements by
designing a controller that minimizes the H∞ norm from the disturbance
inputs, d1, d2, and d3, to the error signals , e1, e2, and e3.
Create the weighting functions that model the design objectives.
Wroad = ss(0.07);
Wroad.u = 'd1';
Wroad.y = 'r';
5-23
Wact = 8*tf([1 50],[1 500]);

Wact.u = 'u';
Wact.y = 'e1';
Wd2 = ss(0.01);
Wd2.u = 'd2';
Wd2.y = 'Wd2';
Wd3 = ss(0.5);
Wd3.u = 'd3';
Wd3.y = 'Wd3';
The constant weight Wroad = 0.07 models broadband road deflections of

magnitude 7 cm. Wact is a highpass filter. This filter penalizes high-frequency
content in the control signal, and thus limits control bandwidth. Wd2 and Wd3
model broadband sensor noise of intensity 0.01 and 0.5, respectively. In a
more realistic design, Wd2 and Wd3 would be frequency dependent to model
the noise spectrum of the displacement and acceleration sensors. The inputs
and outputs of all weighting functions are named to facilitate interconnection.
The notation u and y are shorthand for the InputName and OutputName
properties, respectively.
Specify target functions for the closed-loop response of the system from the
road disturbance, r, to the suspension deflection, sd, and body acceleration, ab.
HandlingTarget = 0.04 * tf([1/8 1],[1/80 1]);

ComfortTarget = 0.4 * tf([1/0.45 1],[1/150 1]);
Targets = [HandlingTarget;ComfortTarget];
Because of the actuator uncertainty and imaginary-axis zeros, the targets

attenuate disturbances only below 10 rad/s. These targets represent the
goals of passenger comfort (small car body acceleration) and adequate road
handling (small suspension deflection).
Plot the closed-loop targets and compare to the open-loop response.
bodemag(qcar({'sd','ab'},'r')*Wroad,'b',Targets,'r--',{1,1000})
grid, title('Response to road disturbance')
legend('Open-loop','Closed-loop target')
5-24
The corresponding performance weights Wsd and Wab are the reciprocals of the
comfort and handling targets. To investigate the trade-off between passenger
comfort and road handling, construct three sets of weights, (βWsd,(1 – β)Wab).
These weights use a blending parameter, β, to modulate the trade-off.
beta = reshape([0.01 0.5 0.99],[1 1 3]);
Wsd = beta/HandlingTarget;
Wsd.u = 'sd';
Wsd.y = 'e3';
Wab = (1-beta)/ComfortTarget;
Wab.u = 'ab';
Wab.y = 'e2';
Wsd and Wab are arrays of weighting functions that correspond to three
different trade-offs: emphasizing comfort (β = 0.01), balancing comfort and
handling (β = 0.5), and emphasizing handling (β = 0.99).
5-25
Connect the quarter-car plant model, actuator model, and weighting functions
to construct the block diagram of the plant model weighted by the objectives.
sdmeas = sumblk('y1 = sd+Wd2');

abmeas = sumblk('y2 = ab+Wd3');
ICinputs = {'d1';'d2';'d3';'u'};
ICoutputs = {'e1';'e2';'e3';'y1';'y2'};
qcaric = connect(qcar(2:3,:),Act,Wroad,Wact,Wab,Wsd,Wd2,Wd3,...
sdmeas,abmeas,ICinputs,ICoutputs);
qcaric is an array of three models, one for each value of the blending
parameter, β. Also, the models in qcaric are uncertain, because they contain
the uncertain actuator model Act.
Nominal H-Infinity Synthesis
Use hinfsyn to compute an H∞ controller for each value of the blending

parameter, β. hinfsyn ignores the uncertainty in the plant models and
synthesizes a controller for the nominal value of each model.
ncont = 1;
nmeas = 2;
K = ss(zeros(ncont,nmeas,3));
gamma = zeros(3,1);
for i=1:3
[K(:,:,i),~,gamma(i)] = hinfsyn(qcaric(:,:,i),nmeas,ncont);
end
The weighted plant model has one control input (ncont), the hydraulic
actuator force. The model also has two measurement outputs (nmeas), which
give the suspension deflection and body acceleration.
Examine the resulting closed-loop H∞ norms, gamma.
gamma
gamma =
0.9410
0.6724
0.8877
5-26
The three H∞ controllers achieve closed-loop H∞ norms of 0.94 (emphasizing

comfort), 0.67 (balancing comfort and handling), and 0.89 (emphasizing
handling).
Construct closed-loop models of the quarter-car plant with the synthesized

controller, corresponding to each of the three blending parameter values.
Compare the frequency response from the road disturbance to xb, sd, and ab
for the passive and active suspensions.
K.u = {'sd','ab'}; K.y = 'u';

CL = connect(qcar,Act.Nominal,K,'r',{'xb';'sd';'ab'});
clf
bodemag(qcar(:,'r'),'b', CL(:,:,1),'r-.', ...
CL(:,:,2),'m-.', CL(:,:,3),'k:',{1,140})
grid
legend('Open-loop','Comfort','Balanced','Handling','location','SouthEast')
title('Body travel, suspension deflection, and body acceleration due to road')
5-27
The solid blue line corresponds to the open-loop response. The other lines
are the closed-loop frequency responses for the different comfort and
handling blends. All three controllers reduce suspension deflection and body
acceleration below the rattlespace frequency (23 rad/s).
Time-Domain Evaluation
To further evaluate the three designs, perform time-domain simulations using

the following road disturbance signal r(t):
 a 1  cos 8 t  , 0  t  0.25
r t  
 0, otherwise.
This signal corresponds to a road bump of height 5 cm.
Create a vector that represents the road disturbance.
5-28
t = 0:0.005:1;
roaddist = zeros(size(t));
roaddist(1:51) = 0.025*(1-cos(8*pi*t(1:51)));
Build the closed-loop model using the synthesized controller.
SIMK = connect(qcar,Act.Nominal,K,'r',{'xb';'sd';'ab';'fs'});
SIMK is a model array containing three closed-loop models, one for each of
the three blending parameter values. Each model in the array represents a
closed loop that is built from the original quarter-car plant model, the nominal
actuator model, and the corresponding synthesized controller.
Simulate and plot the time-domain response of the closed-loop models to the
road disturbance signal.
p1 = lsim(qcar(:,1),roaddist,t);
y1 = lsim(SIMK(1:4,1,1),roaddist,t);
clf
subplot(221)
plot(t,p1(:,1),'b',t,y1(:,1),'r.',t,y2(:,1),'m.',t,y3(:,1),'k.',t,roaddist,'g')
title('Body travel')
ylabel('x_b (m)')
subplot(222)
title('Body acceleration')
ylabel('a_b (m/s^2)')
subplot(223)
title('Suspension deflection')
xlabel('Time (s)')
ylabel('s_d (m)')
subplot(224)
plot(t,zeros(size(t)),'b',t,y1(:,4),'r.',t,y2(:,4),'m.',t,y3(:,4),'k.',t,roaddist,'g')
title('Control force')
5-29
xlabel('Time (s)')
ylabel('f_s (N)')
legend('Open-loop','Comfort','Balanced','Suspension','Road Disturbance','location','SouthEast')
5-30
The simulations show that the body acceleration is smallest for the controller
emphasizing passenger comfort. Body acceleration is largest for the controller
emphasizing suspension deflection. The balanced design achieves a good
tradeoff between body acceleration and suspension deflection.
Robust µ Design
So far you designed H∞ controllers that meet the performance objectives for
the nominal actuator model. However, this model is only an approximation of
the true actuator. To make sure that controller performance is maintained
even with model error and uncertainty, you must design the model to have
robust performance. In this part of the example, you use μ-synthesis to design
a controller that achieves robust performance for the entire family of actuator
models that takes uncertainty into account.
Use D-K iteration to synthesize a controller for the quarter-car model with
actuator uncertainty.
[Krob,~,RPmuval] = dksyn(qcaric(:,:,2),nmeas,ncont);
The model qcaric(:,:,2) is the weighted quarter-car model for the uncertain
model that corresponds to the blending variable β = 0.5.
Examine the resulting μ-synthesis controller.
size(Krob)
State-space model with 1 outputs, 2 inputs, and 11 states.
Build the closed-loop model using the robust controller, Krob.
Krob.u = {'sd','ab'};
Krob.y = 'u';
SIMKrob = connect(qcar,Act.Nominal,Krob,'r',{'xb';'sd';'ab';'fs'});
Simulate and plot the nominal time-domain response to a road bump with
the robust controller.
p1 = lsim(qcar(:,1),roaddist,t);
y1 = lsim(SIMKrob(1:4,1),roaddist,t);
clf
5-31
subplot(221)
plot(t,p1(:,1),'b',t,y1(:,1),'r',t,roaddist,'g')
title('Body travel'), ylabel('x_b (m)')
subplot(222)
plot(t,p1(:,3),'b',t,y1(:,3),'r')
title('Body acceleration'), ylabel('a_b (m/s^2)')
subplot(223)
plot(t,p1(:,2),'b',t,y1(:,2),'r')
title('Suspension deflection'), xlabel('Time (s)'), ylabel('s_d (m)')
subplot(224)
plot(t,zeros(size(t)),'b',t,y1(:,4),'r')
title('Control force'), xlabel('Time (s)'), ylabel('f_s (N)')
legend('Open-loop','Robust design','location','SouthEast')
5-32
These responses are similar to those obtained with the balanced H∞ controller.
Examine the effect of the robust controller on variability caused by model

uncertainty. To do so, simulate the response to a road bump for 120 actuator
models randomly sampled from the uncertain model, Act. Perform this
simulation for both the H∞ and the robust controllers, to compare the results.
Compute an uncertain closed-loop model with the balanced H∞ controller, K.

Sample this model, simulate the sampled models, and plot the results.
5-33
CLU = connect(qcar,Act,K(:,:,2),'r',{'xb','sd','ab'});
nsamp = 120;
rng('default');
figure(1)
clf
lsim(usample(CLU,nsamp),'b',CLU.Nominal,'r+',roaddist,t)
title('Nominal "balanced" design')
Compute an uncertain closed-loop model with the balanced robust controller,

Krob. Sample this model, simulate the sampled models, and plot the results.
CLU = connect(qcar,Act,Krob,'r',{'xb','sd','ab'});
figure(2)
clf
lsim(usample(CLU,nsamp),'b',CLU.Nominal,'r+',roaddist,t)
title('Robust "balanced" design')
5-34
The robust controller reduces variability caused by model uncertainty, and

delivers more consistent performance.
Controller Simplification
The robust controller Krob has eleven states. It is often the case that
controllers synthesized with dksyn have high order. You can use the model
reduction functions to find a lower-order controller that achieves the same
level of robust performance. Use reduce to generate approximations of
various orders. Then, use robustperf to compute the robust performance
margin for each reduced-order approximation.
Create an array of reduced-order controllers.
5-35
NS = order(Krob);
StateOrders = 1:NS;
Kred = reduce(Krob,StateOrders);
Krob is a model array containing a reduced-order controller of every order

from 1 up to the original 11 states.
Compute the robust performance margin for each reduced controller.
CLP = lft(qcaric(:,:,2),Kred);
ropt = robustperfOptions('Sensitivity','off','Display','off','Mussv','a');
PM = robustperf(CLP,ropt);
Compare the robust performance of the reduced- and full-order controllers.
plot(StateOrders,[PM.LowerBound],'b-o',...
StateOrders,repmat(1/RPmuval,[1 NS]),'r');
title('Robust performance as a function of controller order')
legend('Reduced order','Full order')
5-36
There is no significant difference in robust performance between the 8th- and

11th-order controllers. Therefore, you can safely replace Krob by its 8th-order
approximation.
Krob8 = Kred(:,:,8);
You now have a simplified controller, Krob8, that provides robust control with
a balance between passenger comfort and handling.
5-37
Bibliography
[1] Balas, G.J., and A.K. Packard, “The structured singular value
µ-framework,” CRC Controls Handbook, Section 2.3.6, January, 1996, pp.
671-688.
[2] Ball, J.A., and N. Cohen, “Sensitivity minimization in an H∞ norm:

Parametrization of all suboptimal solutions,” International Journal of Control,
Vol. 46, 1987, pp. 785-816.
[3] Bamieh, B.A., and Pearson, J.B., “A general framework for linear periodic
systems with applications to H∞ sampled-data control,” IEEE Transactions on
Automatic Control, Vol. AC-37, 1992, pp. 418-435.
[4] Doyle, J.C., Glover, K., Khargonekar, P., and Francis, B., “State-space
solutions to standard H2 and H∞ control problems,” IEEE Transactions on
Automatic Control, Vol. AC-34, No. 8, August 1989, pp. 831-847.
[5] Fialho, I., and Balas, G.J., “Design of nonlinear controllers for active
vehicle suspensions using parameter-varying control synthesis,” Vehicle
Systems Dynamics, Vol. 33, No. 5, May 2000, pp. 351-370.
[6] Francis, B.A., A course in H∞ control theory, Lecture Notes in Control and
Information Sciences, Vol. 88, Springer-Verlag, Berlin, 1987.
[7] Glover, K., and Doyle, J.C., “State-space formulae for all stabilizing
controllers that satisfy an H∞ norm bound and relations to risk sensitivity,”
Systems and Control Letters, Vol. 11, pp. 167-172, August 1989. International
Journal of Control, Vol. 39, 1984, pp. 1115-1193.
[8] Hedrick, J.K., and Batsuen, T., “Invariant Properties of Automotive

Suspensions,” Proceedings of The Institution of Mechanical Engineers, 204
(1990), pp. 21-27.
[9] Lin, J., and Kanellakopoulos, I., “Road Adaptive Nonlinear Design of
Active Suspensions,” Proceedings of the American Control Conference, (1997),
pp. 714-718.
5-38
Bibliography
[10] Packard, A.K., Doyle, J.C., and Balas, G.J., “Linear, multivariable robust
control with a µ perspective,” ASME Journal of Dynamics, Measurements and
Control: Special Edition on Control, Vol. 115, No. 2b, June, 1993, pp. 426-438.
[11] Skogestad, S., and Postlethwaite, I., Multivariable Feedback Control:

Analysis & Design, John Wiley & Sons, 1996.
[12] Stein, G., and Doyle, J., “Beyond singular values and loopshapes,” AIAA
Journal of Guidance and Control, Vol. 14, Num. 1, January, 1991, pp. 5-16.
[13] Zames, G., “Feedback and optimal sensitivity: model reference

transformations, multiplicative seminorms, and approximate inverses,” IEEE
Transactions on Automatic Control, Vol. AC-26, 1981, pp. 301-320.
5-39

Robust Matlab

Uploaded by

Copyright:

Available Formats

Robust Matlab

Uploaded by

Document Information

Copyright

Available Formats

Share this document

Share or Embed Document

Sharing Options

Did you find this document useful?

Is this content inappropriate?

Copyright:

Available Formats

Robust Matlab

Uploaded by

Copyright:

Available Formats

Robust Control Toolbox™

Getting Started Guide

[email protected] Product enhancement suggestions

The MathWorks, Inc.

Product Requirements . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1-3

Modeling Uncertainty . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1-4

System with Uncertain Parameters ................. 1-5

Worst-Case Performance . . . . . . . . . . . . . . . . . . . . . . . . . . . 1-9

Worst-Case Performance of Uncertain System . . . . . . . . 1-10

Loop-Shaping Controller Design . . . . . . . . . . . . . . . . . . . . 1-12

Model Reduction and Approximation . . . . . . . . . . . . . . . . 1-16

Reduce Order of Synthesized Controller . . . . . . . . . . . . . 1-17

LMI Solvers . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1-21

Extends Control System Toolbox Capabilities . . . . . . . . 1-22

Norms and Singular Values . . . . . . . . . . . . . . . . . . . . . . . . . 2-4

Typical Loop Shapes, S and T Design . . . . . . . . . . . . . . . . 2-6

Using LOOPSYN to Do H-Infinity Loop Shaping . . . . . . 2-15

Loop-Shaping Control Design of Aircraft Model . . . . . . 2-16

Fine-Tuning the LOOPSYN Target Loop Shape Gd to

Mixed-Sensitivity Loop Shaping . . . . . . . . . . . . . . . . . . . . 2-22

Mixed-Sensitivity Loop-Shaping Controller Design . . . 2-24

Loop-Shaping Commands . . . . . . . . . . . . . . . . . . . . . . . . . . 2-26

Model Reduction for Robust Control

Hankel Singular Values . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3-3

Model Reduction Techniques . . . . . . . . . . . . . . . . . . . . . . . 3-5

Approximate Plant Model by Multiplicative Error

Using Modal Algorithms . . . . . . . . . . . . . . . . . . . . . . . . . . . . 3-12

Reducing Large-Scale Models . . . . . . . . . . . . . . . . . . . . . . . 3-16

Normalized Coprime Factor Reduction . . . . . . . . . . . . . . 3-17

Robust Controller Design .......................... 4-9

MIMO Robustness Analysis . . . . . . . . . . . . . . . . . . . . . . . . . 4-14

Summary of Robustness Analysis Tools . . . . . . . . . . . . . . 4-27

H-Infinity Performance . . . . . . . . . . . . . . . . . . . . . . . . . . . . 5-9

Active Suspension Control Design . . . . . . . . . . . . . . . . . . . 5-17

Tuning Fixed Control Architectures

Choosing an Approach to Fixed Control Structure

Difference Between Fixed-Structure Tuning and

Structure of Control System for Tuning With

Set Up Your Control System for Tuning with

Tune MIMO Control System for Specified

What Is hinfstruct? . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-19

Formulating Design Requirements as H-Infinity

Structured H-Infinity Synthesis Workflow . . . . . . . . . . . 6-21

Build Tunable Closed-Loop Model for Tuning with

Tune the Controller Parameters . . . . . . . . . . . . . . . . . . . . 6-29

Interpret the Outputs of hinfstruct . . . . . . . . . . . . . . . . . . 6-30

Validate the Controller Design . . . . . . . . . . . . . . . . . . . . . . 6-31

Set Up Your Control System for Tuning with systune . . 6-35

Specifying Design Requirements for systune . . . . . . . . . 6-37

Tune Controller Against Set of Plant Models . . . . . . . . . 6-39

Speed Up Tuning with Parallel Computing Toolbox

Tuning Control Systems with SYSTUNE . . . . . . . . . . . . . 6-44

Tuning Feedback Loops with LOOPTUNE . . . . . . . . . . . 6-54

Tuning Control Systems in Simulink . . . . . . . . . . . . . . . . 6-62

Building Tunable Models . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-72

Using Design Requirement Objects . . . . . . . . . . . . . . . . . . 6-81

Validating Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 6-102

Tuning Multi-Loop Control Systems . . . . . . . . . . . . . . . . . 6-112

Using Parallel Computing to Accelerate Tuning . . . . . . 6-124

PID Tuning for Setpoint Tracking vs. Disturbance

C=100ss((s+1)/(.001s+1))^3 % LTI controller