Contents

1 SYSTEM............................................................................................................................1
1.1 MODELLING............................................................................................................1
1.2 ANALYSIS...............................................................................................................11
2 CONTROL.......................................................................................................................18
2.1 PID...........................................................................................................................18
2.2 Root locus.................................................................................................................26
2.3 FREQUENCY..........................................................................................................36
2.4 STATE SPACE.........................................................................................................56
2.5 DIGITAL..................................................................................................................71
3 SIMULINK......................................................................................................................86
3.1 MODELLING..........................................................................................................86
3.2 CONTROL.............................................................................................................106
3.3 SIMSCAPE............................................................................................................110

1 SYSTEM
1.1 MODELLING
Inverted Pendulum: System Modelling
Problem setup and design requirements
The system in this example consists of an inverted pendulum mounted to a motorized cart.
The inverted pendulum system is an example commonly found in control system textbooks
and research literature. Its popularity derives in part from the fact that it is unstable without
control, that is, the pendulum will simply fall over if the cart isn't moved to balance it.
Additionally, the dynamics of the system are nonlinear. The objective of the control system is
to balance the inverted pendulum by applying a force to the cart that the pendulum is attached
to. A real-world example that relates directly to this inverted pendulum system is the attitude
control of a booster rocket at take-off.
In this case we will consider a two-dimensional problem where the pendulum is constrained
to move in the vertical plane shown in the figure below. For this system, the control input is
the force that moves the cart horizontally and the outputs are the angular position of the
pendulum and the horizontal position of the cart .
For this example, let's assume the following quantities:
(M) mass of the cart 0.5 kg
(m) mass of the pendulum 0.2 kg
(b) coefficient of friction for cart 0.1 N/m/sec
(l) length to pendulum center of mass 0.3 m
(I) mass moment of inertia of the pendulum 0.006 kg.m^2
(F) force applied to the cart
(x) cart position coordinate
(theta) pendulum angle from vertical (down)
For the PID, root locus, and frequency response sections of this problem, we will be
interested only in the control of the pendulum's position. This is because the techniques used
in these sections are best-suited for single-input, single-output (SISO) systems. Therefore,
none of the design criteria deal with the cart's position. We will, however, investigate the
controller's effect on the cart's position after the controller has been designed. For these
sections, we will design a controller to restore the pendulum to a vertically upward position
after it has experienced an impulsive "bump" to the cart. Specifically, the design criteria are
that the pendulum return to its upright position within 5 seconds and that the pendulum never
move more than 0.05 radians away from vertical after being disturbed by an impulse of
magnitude 1 Nsec. The pendulum will initially begin in the vertically upward equilibrium,
= .
In summary, the design requirements for this system are:
 Settling time for of less than 5 seconds
  Pendulum angle never more than 0.05 radians from the vertical
Employing state-space design techniques, we are more readily able to address a multi-output
system. In our case, the inverted pendulum system is single-input, multi-output (SIMO).
Therefore, for the state-space section of the Inverted Pendulum example, we will attempt to
control both the pendulum's angle and the cart's position. To make the design more
challenging in this section, we will command a 0.2-meter step in the cart's desired position.
Under these conditions, it is desired that the cart achieve its commanded position within 5
seconds and have a rise time under 0.5 seconds. It is also desired that the pendulum settle to
its vertical position in under 5 seconds, and further, that the pendulum angle not travel more
than 20 degrees (0.35 radians) away from the vertically upward position.
In summary, the design requirements for the inverted pendulum state-space example are:
 Settling time for and of less than 5 seconds
 Rise time for of less than 0.5 seconds
 Pendulum angle never more than 20 degrees (0.35 radians) from the vertical
 Steady-state error of less than 2% for and
Force analysis and system equations
Below are the free-body diagrams of the two elements of the inverted pendulum system.

Summing the forces in the free-body diagram of the cart in the horizontal direction, you get
the following equation of motion.
(1)
Note that you can also sum the forces in the vertical direction for the cart, but no useful
information would be gained.
Summing the forces in the free-body diagram of the pendulum in the horizontal direction,
you get the following expression for the reaction force .
(2)
If you substitute this equation into the first equation, you get one of the two governing
equations for this system.

(3)
To get the second equation of motion for this system, sum the forces perpendicular to the
pendulum. Solving the system along this axis greatly simplifies the mathematics. You should
get the following equation.

(4)
To get rid of the and terms in the equation above, sum the moments about the centroid of
the pendulum to get the following equation.
(5)
Combining these last two expressions, you get the second governing equation.

(6)
Since the analysis and control design techniques we will be employing in this example apply
only to linear systems, this set of equations needs to be linearized. Specifically, we will
linearize the equations about the vertically upward equillibrium position, = , and will
assume that the system stays within a small neighborhood of this equillbrium. This
assumption should be reasonably valid since under control we desire that the pendulum not
deviate more than 20 degrees from the vertically upward position. Let represent the
deviation of the pedulum's position from equilibrium, that is, = + . Again presuming a
small deviation ( ) from equilibrium, we can use the following small angle approximations of
the nonlinear functions in our system equations:
(7)
(8)

(9)
After substiting the above approximations into our nonlinear governing equations, we arrive
at the two linearized equations of motion. Note has been substituted for the input .

(10)

(11)
1. Transfer Function
To obtain the transfer functions of the linearized system equations, we must first take the
Laplace transform of the system equations assuming zero initial conditions. The resulting
Laplace transforms are shown below.

(12)

(13)
Recall that a transfer function represents the relationship between a single input and a single
output at a time. To find our first transfer function for the output and an input of we
need to eliminate from the above equations. Solve the first equation for .

(14)
Then substitute the above into the second equation.

(15)
Rearranging, the transfer function is then the following

(16)
where,

(17)
From the transfer function above it can be seen that there is both a pole and a zero at the
origin. These can be canceled and the transfer function becomes the following.

(18)
Second, the transfer function with the cart position as the output can be derived in a
similar manner to arrive at the following.

(19)
2. State-Space
The linearized equations of motion from above can also be represented in state-space form if
they are rearranged into a series of first order differential equations. Since the equations are
linear, they can then be put into the standard matrix form shown below.
(20)

(21)
The matrix has 2 rows because both the cart's position and the pendulum's position are part
of the output. Specifically, the cart's position is the first element of the output and the
pendulum's deviation from its equilibrium position is the second element of .
MATLAB representation
1. Transfer Function
We can represent the transfer functions derived above for the inverted pendulum system
within MATLAB employing the following commands. Note that you can give names to the
outputs (and inputs) to differentiate between the cart's position and the pendulum's position.
Running this code in the command window produces the output shown below.
M = 0.5;
m = 0.2;
b = 0.1;
I = 0.006;
g = 9.8;
l = 0.3;
q = (M+m)*(I+m*l^2)-(m*l)^2;
s = tf('s');
P_cart = (((I+m*l^2)/q)*s^2 - (m*g*l/q))/(s^4 + (b*(I + m*l^2))*s^3/q - ((M +
m)*m*g*l)*s^2/q - b*m*g*l*s/q);
P_pend = (m*l*s/q)/(s^3 + (b*(I + m*l^2))*s^2/q - ((M + m)*m*g*l)*s/q - b*m*g*l/q);
sys_tf = [P_cart ; P_pend];
inputs = {'u'};
outputs = {'x'; 'phi'};
set(sys_tf,'InputName',inputs)
set(sys_tf,'OutputName',outputs)
sys_tf
sys_tf =

From input "u" to output...

4.182e-06 s^2 - 0.0001025
x: ---------------------------------------------------------
2.3e-06 s^4 + 4.182e-07 s^3 - 7.172e-05 s^2 - 1.025e-05 s

1.045e-05 s
phi: -----------------------------------------------------
2.3e-06 s^3 + 4.182e-07 s^2 - 7.172e-05 s - 1.025e-05

Continuous-time transfer function.

2. State-Space
We can also represent the system using the state-space equations. The following additional
MATLAB commands create a state-space model of the inverted pendulum and produce the
output shown below when run in the MATLAB command window. Again note that the names
of the inputs, outputs, and states can be specified to make the model easier to understand.
M = .5;
m = 0.2;
b = 0.1;
I = 0.006;
g = 9.8;
l = 0.3;

p = I*(M+m)+M*m*l^2; %denominator for the A and B matrices

A = [0 1 0 0;
0 -(I+m*l^2)*b/p (m^2*g*l^2)/p 0;
 0 0 0 1;
0 -(m*l*b)/p m*g*l*(M+m)/p 0];
B=[ 0;
(I+m*l^2)/p;
0;
m*l/p];
C = [1 0 0 0;
0 0 1 0];
D = [0;
0];

states = {'x' 'x_dot' 'phi' 'phi_dot'};

inputs = {'u'};
outputs = {'x'; 'phi'};

sys_ss = ss(A,B,C,D,'statename',states,'inputname',inputs,'outputname',outputs)
sys_ss =

A=
x x_dot phi phi_dot
x 0 1 0 0
x_dot 0 -0.1818 2.673 0
phi 0 0 0 1
phi_dot 0 -0.4545 31.18 0

B=
u
x 0
x_dot 1.818
 phi 0
phi_dot 4.545

C=
x x_dot phi phi_dot
x 1 0 0 0
phi 0 0 1 0

D=
u
x 0
phi 0

Continuous-time state-space model.

The above state-space model can also be converted into transfer function form employing
the tf command as shown below. Conversely, the transfer function model can be converted
into state-space form using the ss command.
sys_tf = tf(sys_ss)
sys_tf =

From input "u" to output...

1.818 s^2 + 1.615e-15 s - 44.55
x: --------------------------------------
s^4 + 0.1818 s^3 - 31.18 s^2 - 4.455 s

4.545 s - 1.277e-16
phi: ----------------------------------
s^3 + 0.1818 s^2 - 31.18 s - 4.455
Continuous-time transfer function.

Examining the above, note the existance of some terms with very small coefficients. These
terms should actually be zero and show up due to numerical round-off errors that accumulate
in the conversion algorithms that MATLAB employs. If you set these coefficients to zero,
then the above transfer function models will match those generated earlier in the Transfer
Function section of the example.
1.2 ANALYSIS
Inverted pendulum: System Analysis
2 Contents
 Open-loop impulse response
 Open-loop step response
From the main problem, we derived the open-loop transfer functions of the
inverted pendulum system as the following.

(1)

(2)
where
(3)
Recall that the above two transfer functions are valid only for small values of
the angle , which is the angular displacement of the pendulum from the
vertically upward position. Also, the absolute pendulum angle is equal to +
.
For the original problem setup and the derivation of the above transfer
functions, please refer to the Inverted Pendulum: System Modeling page.
Considering the response of the pendulum to a 1-Nsec impulse applied to the
cart, the design requirements for the pendulum are:
 Settling time for of less than 5 seconds
 Pendulum angle never more than 0.05 radians from the vertical
Additionally, the requirements for the response of the system to a 0.2-meter step
command in cart position are:
 Settling time for and of less than 5 seconds
 Rise time for of less than 0.5 seconds
 Pendulum angle never more than 20 degrees (0.35 radians) from the
vertical
3 Open-loop impulse response
We will begin by looking at the open-loop response of the inverted pendulum
system. Create a new m-file and type in the following commands to create the
system model (refer to the main problem for the details of getting these
commands).
M = 0.5;
m = 0.2;
b = 0.1;
I = 0.006;
g = 9.8;
l = 0.3;
q = (M+m)*(I+m*l^2)-(m*l)^2;
s = tf('s');

P_cart = (((I+m*l^2)/q)*s^2 - (m*g*l/q))/(s^4 + (b*(I + m*l^2))*s^3/q - ((M +

m)*m*g*l)*s^2/q - b*m*g*l*s/q);

P_pend = (m*l*s/q)/(s^3 + (b*(I + m*l^2))*s^2/q - ((M + m)*m*g*l)*s/q -

b*m*g*l/q);

sys_tf = [P_cart ; P_pend];

inputs = {'u'};
outputs = {'x'; 'phi'};

set(sys_tf,'InputName',inputs)
set(sys_tf,'OutputName',outputs)
We can now examine the open-loop impulse response of the system.
Specifically, we will examine how the system responds to an impulsive force
applied to the cart employing the MATLAB command impulse. Add the
following commands onto the end of the m-file and run it in the MATLAB
command window to get the associated plot shown below.
t=0:0.01:1;
impulse(sys_tf,t);
title('Open-Loop Impulse Response')

As you can see from the plot, the system response is entirely unsatisfactory. In
fact, it is not stable in open loop. Although the pendulum's position is shown to
increase past 100 radians (15 revolutions), the model is only valid for small .
You can also see that the cart's position moves infinitely far to the right, though
there is no requirement on cart position for an impulsive force input.
The poles of a system can also tell us about its time response. Since our system
has two outputs and one input, it is described by two transfer functions. In
general, all transfer functions from each input to each output of a multi-input,
multi-output (MIMO) system will have the same poles (but different zeros)
unless there are pole-zero cancellations. We will specifically examine the poles
and zeros of the system using the MATLAB function zpkdata. The
parameter 'v' shown below returns the poles and zeros as column vectors rather
than as cell arrays.
The zeros and poles of the system where the pendulum position is the output are
found as shown below:
[zeros poles] = zpkdata(P_pend,'v')
zeros =
0
poles =
5.5651
-5.6041
-0.1428
Likewise, the zeros and poles of the system where the cart position is the output
are found as follows:
[zeros poles] = zpkdata(P_cart,'v')
zeros =
4.9497
-4.9497
poles =
0
5.5651
-5.6041
-0.1428
As predicted, the poles for both transfer functions are identical. The pole at
5.5651 indicates that the system is unstable since the pole has positive real part.
In other words, the pole is in the right half of the complex s-plane. This agrees
with what we observed above.
4 Open-loop step response
Since the system has a pole with positive real part its response to a step input
will also grow unbounded. We will verify this using the lsim command which
can be employed to simulate the response of LTI models to arbitrary inputs. In
this case, a 1-Newton step input will be used. Adding the following code to your
m-file and running it in the MATLAB command window will generate the plot
shown below.
t = 0:0.05:10;
u = ones(size(t));
[y,t] = lsim(sys_tf,u,t);
plot(t,y)
title('Open-Loop Step Response')
axis([0 3 0 50])
legend('x','phi')

You can also identify some important characteristics of the response using
the lsiminfo command as shown.
step_info = lsiminfo(y,t);
cart_info = step_info(1)
pend_info = step_info(2)
cart_info =
struct with fields:

SettlingTime: 9.9959
Min: 0
MinTime: 0
Max: 8.7918e+21
MaxTime: 10
pend_info =
struct with fields:

SettlingTime: 9.9959
Min: 0
MinTime: 0
Max: 1.0520e+23
MaxTime: 10
The above results confirm our expectation that the system's response to a step
input is unstable.
It is apparent from the analysis above that some sort of control will need to be
designed to improve the response of the system. Four example controllers are
included with these tutorials: PID, root locus, frequency response, and state
space. You may select a choice from the menu to the left for further details.
Note: The solutions shown in the PID, root locus, and frequency response
examples may not yield a workable controller for the inverted pendulum
problem. As stated previously, when we treat the inverted pendulum as a single-
input, single-output system, we ignore the position of the cart, . Where possible
in these examples, we will show what happens to the cart's position when a
controller is implemented on the system.
 2 CONTROL
2.1 PID
INVERTED PENDULUM: PID CONTROLLER DESIGN

Contents
 System structure
 PID control
 What happens to the cart's position?
In this page we will design a PID controller for the inverted pendulum system. In the design
process we will assume a single-input, single-output plant as described by the following
transfer function. Otherwise stated, we will attempt to control the pendulum's angle without
regard for the cart's position.

(1)
where,

(2)
More specifically, the controller will attempt to maintain the pendulum vertically upward
when the cart is subjected to a 1-Nsec impulse. Under these conditions, the design criteria
are:
 Settling time of less than 5 seconds
 Pendulum should not move more than 0.05 radians away from the vertical
For the original problem setup and the derivation of the above transfer function, please
consult the Inverted Pendulum: System Modeling page.
System structure
The structure of the controller for this problem is a little different than the standard control
problems you may be used to. Since we are attempting to control the pendulum's position,
which should return to the vertical after the initial disturbance, the reference signal we are
tracking should be zero. This type of situation is often referred to as a Regulator problem.
The external force applied to the cart can be considered as an impulsive disturbance. The
schematic for this problem is depicted below.
You may find it easier to analyze and design for this system if we first rearrange the
schematic as follows.

The resulting transfer function for the closed-loop system from an input of force to an
output of pendulum angle is then determined to be the following.

(3)
Before we begin designing our PID controller, we first need to define our plant within
MATLAB. Create a new m-file and type in the following commands to create the plant model
(refer to the main problem for the details of getting these commands).
M = 0.5;
m = 0.2;
b = 0.1;
I = 0.006;
g = 9.8;
l = 0.3;
q = (M+m)*(I+m*l^2)-(m*l)^2;
s = tf('s');
P_pend = (m*l*s/q)/(s^3 + (b*(I + m*l^2))*s^2/q - ((M + m)*m*g*l)*s/q - b*m*g*l/q);
Next we will define a PID controller.
PID control
This closed-loop transfer function can be modeled in MATLAB by copying the following
code to the end of your m-file (whether you're using the transfer function form or the state-
space representation of the plant). Specifically, we define our controller using the pid object
within MATLAB. We then use the feedback command to generate the closed-loop transfer
function as depicted in the figure above where the disturbance force is the input and
the deviation of the pendulum angle from the vertical is the output.
Kp = 1;
Ki = 1;
Kd = 1;
C = pid(Kp,Ki,Kd);
T = feedback(P_pend,C);
Now we can begin to tune our controller. First let's examine the response of the closed-loop
system to an impulse disturbance for this initial set of control gains. Enter the following code
to the end of your m-file and run in the MATLAB command window. You should generate
the response plot shown below.
t=0:0.01:10;
impulse(T,t)
title({'Response of Pendulum Position to an Impulse Disturbance';'under PID Control: Kp =
1, Ki = 1, Kd = 1'});
This response is still not stable. Let's begin to modify the response by increasing the
proportional gain. Increase the variable to see what effect it has on the response. If you
modify your m-file to the following where = 100 and run in the command window, you
should get the response plot shown below.
Kp = 100;
Ki = 1;
Kd = 1;
C = pid(Kp,Ki,Kd);
T = feedback(P_pend,C);
t=0:0.01:10;
impulse(T,t)
axis([0, 2.5, -0.2, 0.2]);
title({'Response of Pendulum Position to an Impulse Disturbance';'under PID Control: Kp =
100, Ki = 1, Kd = 1'});
Right-clicking on the resulting plot and choosing Characteristics from the resulting menu
allows you to identify important characteristics of the response. Specifically, the settling time
of the response is determined to be 1.64 seconds, which is less than the requirement of 5
seconds. Since the steady-state error approaches zero in a sufficiently fast manner, no
additional integral action is needed. You can set the integral gain constant to zero to see
for yourself that some integral control is needed. The peak response, however, is larger than
the requirement of 0.05 radians. Recall that overshoot often can be reduced by increasing the
amount of derivative control. After some trial and error it is found that a derivative gain of
= 20 provides a satisfactory response. Modifying your m-file as follows and re-running
should produce the response plot shown below
Kp = 100;
Ki = 1;
Kd = 20;
C = pid(Kp,Ki,Kd);
T = feedback(P_pend,C);
t=0:0.01:10;
impulse(T,t)
axis([0, 2.5, -0.2, 0.2]);
title({'Response of Pendulum Position to an Impulse Disturbance';'under PID Control: Kp =
100, Ki = 1, Kd = 20'});

As you can see, the overshoot has been reduced so that the pendulum does not move more
than 0.05 radians away from the vertical. Since all of the given design requirements have
been met, no further iteration is needed.
What happens to the cart's position?
At the beginning of this page, a block diagram for the inverted pendulum system was given.
The diagram was not entirely complete. The block representing the response of the cart's
position was not included because that variable is not being controlled. It is interesting
though, to see what is happening to the cart's position when the controller for the pendulum's
angle is in place. To see this we need to consider the full system block diagram as shown in
the following figure.
Rearranging, we get the following block diagram.

In the above, the block is the controller designed for maintaining the pendulum vertical.
The closed-loop transfer function from an input force applied to the cart to an output of
cart position is, therefore, given by the following.

(4)
Referring to the Inverted Pendulum: System Modeling page, the transfer function for
is defined as follows.

(5)
where,

(6)
Adding the following commands to your m-file (presuming and are still defined)
will generate the response of the cart's position to the same impulsive disturbance we have
been considering.
P_cart = (((I+m*l^2)/q)*s^2 - (m*g*l/q))/(s^4 + (b*(I + m*l^2))*s^3/q - ((M +
m)*m*g*l)*s^2/q - b*m*g*l*s/q);
T2 = feedback(1,P_pend*C)*P_cart;
t = 0:0.01:5;
impulse(T2, t);
title({'Response of Cart Position to an Impulse Disturbance';'under PID Control: Kp = 100,
Ki = 1, Kd = 20'});

As you can see, the cart moves in the negative direction with approximately constant velocity.
Therefore, although the PID controller stabilizes the angle of the pendulum, this design
would not be feasible to implement on an actual physical system.
2.2 Root locus
INVERTED PENDULUM: ROOT LOCUS DESIGN
Contents
 System structure
 Root locus design
 PID control
 What happens to the cart's position?
In this page we will design a controller for the inverted pendulum system using the root locus
design method. In the design process we will assume a single-input, single-output plant as
described by the following transfer function. Otherwise stated, we will attempt to control the
pendulum's angle without regard for the cart's position.

(1)
where,

(2)
More specifically, the controller will attempt to maintain the pendulum vertically upward
when the cart is subjected to a 1-Nsec impulse. Under these conditions, the design criteria
are:
 Settling time of less than 5 seconds
 Pendulum should not move more than 0.05 radians away from the vertical
For the original problem setup and the derivation of the above transfer function, please
consult the Inverted Pendulum: System Modeling page.
System structure
The structure of the controller for this problem is a little different than the standard control
problems you may be used to. Since we are attempting to control the pendulum's position,
which should return to the vertical after the initial disturbance, the reference signal we are
tracking should be zero. This type of situation is often referred to as a Regulator problem.
The external force applied to the cart can be considered as an impulsive disturbance. The
schematic for this problem is depicted below.
You may find it easier to analyze and design for this system if we first rearrange the
schematic as follows.

The resulting transfer function for the closed-loop system from an input of force to an
output of pendulum angle is then determined to be the following.

(3)
Before we begin designing our controller, we first need to define our plant within MATLAB.
Create a new m-file and type in the following commands to create the plant model (refer to
the main problem for the details of getting these commands).
M = 0.5;
m = 0.2;
b = 0.1;
I = 0.006;
g = 9.8;
l = 0.3;
q = (M+m)*(I+m*l^2)-(m*l)^2;
s = tf('s');
P_pend = (m*l*s/q)/(s^3 + (b*(I + m*l^2))*s^2/q - ((M + m)*m*g*l)*s/q - b*m*g*l/q);
Root locus design
We will now begin to design a controller for our system employing a root locus design
method. We can use the MATLAB command rlocus for generating the root locus plots.
Adding the following commands to your m-file and running it in the MATLAB command
window will create the root locus plot shown below. This plot displays all possible closed-
loop pole locations as a simple proportional control gain is varied from 0 to infinity. The
root locus is the same whether the multiplicative gain is in the forward or feedback path of
the closed-loop system.
rlocus(P_pend)
title('Root Locus of Plant (under Proportional Control)')

As you can see, one of the branches of the root locus is entirely in the right-half of the
complex s-plane. This means that no matter the choice of gain , there will always be a
closed-loop pole in the right-half plane making the system's impulse response unstable.
To solve this problem, we need to add a pole at the origin (an integrator) via the controller to
cancel the plant zero at the origin. This addition will produce two closed-loop poles in the
right-half plane. In our subsequent design we can then modify our controller to draw these
poles into the left-half plane, thereby stabilizing the closed-loop system. Modifying your m-
file with the following commands and re-running in the MATLAB command window will
produce the root locus plot shown below.
C = 1/s;
rlocus(C*P_pend)
title('Root Locus with Integral Control')

Let's also examine the locations of the system's open-loop poles and zeros so that we may
begin to think about how to draw the root locus branches into the left-half plane. Entering the
following commands into the MATLAB command window will generate the following
output.
zeros = zero(C*P_pend)
poles = pole(C*P_pend)
zeros =
0
poles =
0
5.5651
-5.6041
 -0.1428
As you can see, there are four poles and only one zero. This means that the root locus will
have three asymptotes: one along the real axis in the negative direction, and the other two at
120 degree angles to this one.
This configuration is also unsatisfactory because we still have branches of the root locus that
are entirely in the right-half complex plane. In general, we can pull the branches of our root
locus to the left in the complex plane by adding zeros to our system. Adding a zero to our
controller will reduce the number of asymptotes from three to two. These two asymptotes
will be parallel to the imaginary axis and will intersect the real axis at the
location s calculated from the following expression.

(4)
Therefore, for our system as described so far, we have the following assuming a minimum-
phase zero (negative).

(5)
Based on the above, the farthest we can pull the asymptotes to the left in the complex plane is
approximately -0.1 for a negligibly small zero. Recall that 2% settling time can be estimated
from the following equation.

(6)
Therefore, dominant closed-loop poles with real parts that approach -0.1 will not be sufficient
to meet the 5 second settling time that we require.
PID control
In the above discussion we demonstrated that adding a zero to our integral controller could
pull the branches of the root locus to the left in the complex plane, but we were not able to
the pull the dominant branches far enough to the left. A possible solution is to add yet another
zero. If we place both zeros on the negative real axis between the two plant poles, then the
two branches in the right-half plane will be pulled into the left-half plane and will terminate
at these two zeros. Let's specifically evaluate the root locus for a controller with an integrator
and zeros at -3 and -4. Note that this controller is actually a PID controller. We can create this
controller within MATLAB using the zpk command which creates a model by specifying the
zeros, poles, and gain of the system. Modifying your m-file with the following commands
and re-running will produce the root locus plot shown below.
z = [-3 -4];
p = 0;
k = 1;
C = zpk(z,p,k);
rlocus(C*P_pend)
title('Root Locus with PID Controller')

Examining the above root locus helps us to determine whether or not our given requirements
can be met. Specifically, since it is desired that the settling time of the system be less than 5
seconds, the real parts of our dominant closed-loop poles should be less than approximately -
4/5 = -0.8. In other words, our dominatnt closed-loop poles should be located in the
complex s-plane to the left of a vertical line at s = -0.8. Inspection of the above shows that
this is possible. Since it is also desired that the pendulum not move more than 0.05 radians
away from vertical, we also want to ensure that the closed-loop system has sufficient
damping. Placing the dominant closed-loop poles near the real axis will increase the system's
damping (small ).
To find the gain corresponding to a specific point on the root locus, we can use
the rlocfind command. Specifically, enter the command [k,poles] = rlocfind(C*P_pend) in the
MATLAB command window.
Then go to the plot and select a point on the root locus on left side of the loop, close to the
real axis as shown below with the small + marks. Selecting these poles will ensure that the
system settles sufficiently fast and, hopefully, that it has sufficient damping.
After doing this, you should see an output like the following in the MATLAB command
window.
Select a point in the graphics window

selected_point =

-3.5367 + 0.7081i

k=

20.2396

poles =
 0
-85.1333
-3.5232 + 0.7086i
-3.5232 - 0.7086i

Note that the values returned in your MATLAB command window may not be exactly the
same, but they should at least have the same order of magnitude.
Then we can check the impulse response of our closed-loop system to see if our requirements
are actually met for a gain of approximately 20. Add the following commands to your m-
file and re-run to generate a closed-loop impulse response like the one shown below.
K = 20;
T = feedback(P_pend,K*C);
impulse(T)
title('Impulse Disturbance Response of Pendulum Angle under PID Control');

Examination of the above demonstrates that all of the given requirements are met.
What happens to the cart's position?
At the beginning of this page, a block diagram for the inverted pendulum system was given.
The diagram was not entirely complete. The block representing the response of the cart's
position was not included because that variable is not being controlled. It is interesting
though, to see what is happening to the cart's position when the controller for the pendulum's
angle is in place. To see this we need to consider the full system block diagram as shown in
the following figure.

Rearranging, we get the following block diagram.

In the above, the block is the controller designed for maintaining the pendulum vertical.
The closed-loop transfer function from an input force applied to the cart to an output of
cart position is, therefore, given by the following.

(7)
Referring to the Inverted Pendulum: System Modeling page, the transfer function for
is defined as follows.
(8)
where,

(9)

Adding the following commands to your m-file (presuming and are still defined)
will generate the response of the cart's position to the same impulsive disturbance we have
been considering.
P_cart = (((I+m*l^2)/q)*s^2 - (m*g*l/q))/(s^4 + (b*(I + m*l^2))*s^3/q - ((M +
m)*m*g*l)*s^2/q - b*m*g*l*s/q);
T2 = feedback(1,P_pend*C)*P_cart;
t = 0:0.01:8.5;
impulse(T2, t);
title('Impulse Disturbance Response of Cart Position under PID Control');

As you can see, the cart's position goes unstable for this impulse disturbance. Therefore,
although the PID controller stabilizes the angle of the pendulum, this design would not be
feasible to implement on an actual physical system.
2.3 FREQUENCY
Inverted Pendulum: Frequency Domain Methods for Controller Design
Contents
 System structure
 Closed-loop response without compensation
 Closed-loop response with compensation
 What happens to the cart's position?
In this page we will design a controller for the inverted pendulum system using a frequency
response design method. In the design process we will assume a single-input, single-output
plant as described by the following transfer function. Otherwise stated, we will attempt to
control the pendulum's angle without regard for the cart's position.

(1)
where,

(2)
The controller we are designing will specifically attempt to maintain the pendulum vertically
upward when the cart is subjected to a 1-Nsec impulse. Under these conditions, the design
criteria are:
 Settling time of less than 5 seconds
 Pendulum should not move more than 0.05 radians away from the vertical
For the original problem setup and the derivation of the above transfer function, please
consult the Inverted Pendulum: System Modeling page.
Note: Applying a frequency response design approach is relatively challenging in the case of
this example because the open-loop system is unstable. That is, the open-loop transfer
function has a pole in the right-half complex plane. For this reason, attempting this example
is not recommended if you are just attempting to learn the basics of applying frequency
response techniques. This problem is better suited for more advanced students who wish to
learn about some nuances of the frequency response design approach.
System structure
The structure of the controller for this problem is a little different than the standard control
problems you may be used to. Since we are attempting to control the pendulum's position,
which should return to the vertical after the initial disturbance, the reference signal we are
tracking should be zero. This type of situation is often referred to as a regulator problem. The
external force applied to the cart can be considered as an impulsive disturbance. The
schematic for this problem is depicted below.

You may find it easier to analyze and design for this system if we first rearrange the
schematic as follows.

The resulting transfer function for the closed-loop system from an input of force to an
output of pendulum angle is then determined to be the following.

(3)
Before we begin designing our controller, we first need to define our plant within MATLAB.
Create a new m-file and type in the following commands to create the plant model (refer to
the main problem for the details of getting these commands).
M = 0.5;
m = 0.2;
b = 0.1;
I = 0.006;
g = 9.8;
l = 0.3;
q = (M+m)*(I+m*l^2)-(m*l)^2;
s = tf('s');
P_pend = (m*l*s/q)/(s^3 + (b*(I + m*l^2))*s^2/q - ((M + m)*m*g*l)*s/q - b*m*g*l/q);
As mentioned above, this system is unstable without control. We can prove this to ourselves
by employing the MATLAB command zpkdata. In this case, zpkdata returns the zeros and
poles for the transfer function. The added parameter 'v' returns the outputs in the form of
vectors instead of cell arrays and can only be employed with single-input, single-output
models. Entering the following code in the MATLAB command window generates the output
shown below.
[zeros poles] = zpkdata(P_pend,'v')
zeros =
0
poles =
5.5651
-5.6041
-0.1428

Closed-loop response without compensation

We will now examine the response of the closed-loop system without compensation before
we begin to design our controller. In this example, we will employ the Control System
Designer for examining the various analysis plots rather than employing individual
commands such as bode, nyquist, and impulse. The Control System Designer is an
interactive tool with graphical user interface (GUI) which can be launched by the MATLAB
command controlSystemDesigner as shown below, or by going to the APPS tab of the
MATLAB toolstrip and clicking on the app icon under Control System Design and
Analysis.
controlSystemDesigner('bode',P_pend)

The additional parameter 'bode' opens the Control System Designer window with the Bode
plot and closed-loop step response of the system (which was passed to the function)
as shown below.
We can then modify the system architecture being employed to reflect the fact that our
controller is in the feedback path of our system as discussed above. This is accomplished
from within the Control System Designer window by clicking on the icon labeled Edit
Architecture. Then in the new window, modify the default configuration by clicking on the
second architecture under Select Control Architecture to match the form shown below.
Next we will begin to examine some of the analysis plots for our system. Recall that we can
assess the closed-loop stability of our system based on the open-loop frequency response. In
the case of this example, the open-loop transfer function of our system is given by the
following.

(4)
Note, this is true even though the controller is in the feedback path of the system.
In other examples we have specifically employed a Bode plot representation of the open-loop
frequency response. For our system without compensation, we could employ the MATLAB
code bode(P_Pend) to generate this Bode plot. Instead, we will use the Control System
Designer that we are employing currently. The open-loop Bode plot of our system is already
open, but if it weren't, or if we wished to change the type of plot we are employing for
design, we could open a new plot from under the New Plot tab of the Control System
Designer window as shown below.
Examination of the above Bode plot shows that the magnitude is less than 0 dB and the phase
is greater than -180 degrees for all frequencies. For a minimum-phase system this would
indicate that the closed-loop system is stable with infinite gain margin and infinite phase
margin. However, since our system has a pole in the right-half complex plane, our system is
nonminimum phase and the closed-loop system is actually unstable. We will prove this to
ourselves by examining a couple of other analysis plots.
In general, when dealing with nonminimum-phase systems it is preferrable to analyze relative
stability using the Nyquist plot of the open-loop transfer function. The Nyquist plot is also
preferred when analyzing higher-order systems. This is because the Bode plot shows
frequencies that are 360 degrees apart as being different when in fact they are the same.
Since the Nyquist plot is a polar-type plot, this ambiguity is removed.
In order to generate additional plots to better understand the closed-loop performance of the
system, click on the New Plot menu in the toolstrip of the Control System
Designer window. The Control System Designer allows the user to view up to eight plots at
the same time for analysis. These plots can be viewed for a number of options, such as the
open-loop and closed-loop system. We will view the Nyquist plot for the open-loop system
and the impulse response of the closed-loop system by following the steps given below.
1. Under the New Plot menu, select New Nyquist under Create New Plot. A new window
named New Nyquist to plot should appear. In the dropdown menu Select Response to Plot,
select New Open-Loop Response. In the boxes Specify the open-loop response at the
following locations and Specify the open-loop response with the following loops open,
select . Then click on Plot. You can also select LoopTransfer_C in the dropdown Select
Response to Plot to obtain the same Nyquist plot.

Under the New Plot menu, select New Impulse under Create New Plot. A new window
named New Impulse to plot should appear. In the dropdown menu Select Response to Plot,
select IOTransfer_r2y. Then click on Plot.
The resulting figure should then appear as follows.
Examination of the above impulse response plot shows that the closed-loop system is
unstable. This also can be verified from the open-loop Nyquist plot by applying the Nyquist
stability criterion which is stated below.
(5)
Where is the number of closed-loop poles in the right-half plane, is the number of open-
loop poles in the right-half plane, and is the number of clockwise encirclements of the
point -1 by the open-loop Nyquist plot.
From our previous discussion we know that our system has one open-loop pole in the right-
half plane ( = 1) and from examination of the open-loop Nyquist plot above we can see that
there are no encirclements of the point -1 ( = 0). Therefore, = 0 + 1 = 1 and the closed-
loop system has 1 pole in the right-half plane indicating that it is indeed unstable.
Closed-loop response with compensation
Since the closed-loop system is unstable without compensation, we need to use our controller
to stabilize the system and meet the given requirements. Our first step will be to add an
integrator to cancel the zero at the origin. To add an integrator, you may right-click on the
Bode plot that is already open and choose Add Pole/Zero > Integrator from the resulting
menu. The resulting Bode plot that is generated is shown below and reflects the low
frequency behavior we would expect.
If you look at the transfer function closely, you will notice that there is a pole-zero
cancellation at the origin. However, the pole-zero cancellation is not automatically
recognized by the Control System Designer which can cause numerical errors in the
subsequent plots. Therefore, you will need exit out of the Control System Designer tool and
re-open the tool with the integrator already added to the plant as shown below.
controlSystemDesigner('bode',minreal(P_pend*(1/s)))

Since the integrator is bundled with the plant which is in the forward path, while our
controller is actually in the feedback path, we will not analyze the closed-loop response of
this system from within the Control System Designer. However, the open-loop transfer
function is unchanged by whether the controller is in the forward or feedback path, therefore,
we can still use the plots of the open-loop system for analysis and design. The resulting Bode
plot that is generated is shown below.
Even with the addition of this integrator, the closed-loop system is still unstable. We can
attempt to better understand the instability (and how to resolve it) by looking more closely at
the Nyquist plot. We can generate this analysis plot in the same manner we did previously.
Following these steps generates a figure like the one given below.
Notice that the open-loop Nyquist plot now encircles the -1 point in the clockwise direction.
This means that the closed-loop system now has two poles in the right-half plane (
). Hence, the closed-loop system is still unstable. We need to add
phase in order to get a counterclockwise encirclement. We will do this by adding a zero to our
controller. For starters, we will place this zero at -1 and view the resulting plots. This action
can be achieved graphically by right-clicking on the plot as we did previously. Instead, we
will add the zero from the Compensator Editor window of the Control System Designer,
which can be opened by right-clicking on the Nyquist plot and then selecting Edit
Compensator. Then, right-click in the Dynamics section of the Compensator
Editor window and select Add Pole/Zero > Real Zero from the resulting menu. By default,
the location of the resulting zero is -1. The resulting window should appear as shown in the
figure below.
This additional zero will automatically change the Bode and Nyquist plots that are already
open. The resulting Nyquist plot should appear as shown below.

As you can see, this change did not provide enough phase. The encirclement around -1 is still
clockwise. We will try adding a second zero at -1 in the same manner as was described above.
The resulting Nyquist diagram is shown below.
We still have one clockwise encirclement of the -1 point. However, if we add some gain we
can increase the magnitude of each point of the Nyquist plot, thereby increasing the radius of
the counterclockwise circle such that it encircles the -1 point. This results in = -1 where
is negative because the encirclement is counterclockwise. To achieve this, you can manually
enter a new gain value in the Compensator Editor window of the Control System
Designer window. Alternatively, you can modify the gain graphically from the Bode plot.
Specifically, go to the Bode plot window, click on the magnitude plot, and drag the curve up
until you have shifted the Nyquist plot far enough to the left that the counterclockwise circle
encompasses the point -1. Note that this is achieved for a gain value of approximately 3.8.
We will continue to increase the gain to a magnitude of approximately 10 as indicated in
the Preview portion of the Control System Designer window. The resulting Bode plot
should appear as shown in the figure below.
The corresponding Nyquist plot should then appear as follows.
From our previous discussion we know that = 1, now that = -1, we have = -1 + 1 = 0
closed-loop poles in the right-half plane indicating that our closed-loop system is stable. We
can verify the stability of our system and determine whether or not the other requirements are
met by examining the system's response to a unit impulse force disturbance. Since the
integrator of our controller is currently bundled with the plant, we will exit from the Control
System Designer and generate the closed-loop impulse response from the MATLAB
command line. So far, the controller we have designed has the form given below.

(6)
Adding the following code to your m-file will construct the closed-loop transfer function
from an input of to an output of . Running your m-file at the command line will then
generate an impulse response plot as shown below.
K = 10;
C = K*(s+1)^2/s;
T = feedback(P_pend,C);

t = 0:0.01:10;
impulse(T,t), grid
title({'Response of Pendulum Position to an Impulse Disturbance';'under Closed-loop
Control'});

From examination of the figure above, it is apparent that the response of the system is now
stable. However, the pendulum position overshoots past the required limit of 0.05 radians and
the settle time is on the verge of being greater than the requirement of 5 seconds. Therefore,
we will now concentrate on improving the response. We can use the Control System
Designer to see how changing the controller gain and the location of the zeros affects the
system's frequency response plots. In this case, increasing the controller gain increases the
system's phase margin which should help reduce the overshoot of the response. Furthermore,
trial and error shows that moving one of the zeros farther to the left in the complex plane
(more negative) makes the response faster. This change also increases the overshoot, but this
is offset by the increase in gain. Experimentation demonstrates that the following controller
satisfies the given requirements.

(7)
Modify your m-file as shown and re-run at the command line to generate the impulse
response plot given below.
K = 35;
C = K*(s+1)*(s+2)/s;
T = feedback(P_pend,C);
t = 0:0.01:10;
impulse(T, t), grid
title({'Response of Pendulum Position to an Impulse Disturbance';'under Closed-loop
Control'});

Our response has met our design goals. Feel free to vary the parameters further to observe
what happens. Note that it is also possible to generate the system's response to an impulse
disturbance from within the Control System Designer. Including the integrator with the
controller, rather than with the plant, as is proper will result in the numerical abberation
observed above, but will not significantly affect the impulse
response generated by MATLAB.
What happens to the cart's position?
At the beginning of this page, a block diagram for the inverted pendulum system was given.
The diagram was not entirely complete. The block representing the response of the cart's
position was not included because that variable is not being controlled. It is interesting
though, to see what is happening to the cart's position when the controller for the pendulum's
angle is in place. To see this we need to consider the full system block diagram as shown in
the following figure.
Rearranging, we get the following block diagram.

In the above, the block is the controller designed for maintaining the pendulum vertical.
The closed-loop transfer function from an input force applied to the cart to an output of
cart position is, therefore, given by the following.

(8)
Referring to the Inverted Pendulum: System Modeling page, the transfer function for
is defined as follows.

(9)
where,

(10)
Adding the following commands to your m-file (presuming and are still defined)
will generate the response of the cart's position to the same impulsive disturbance we have
been considering.
P_cart = (((I+m*l^2)/q)*s^2 - (m*g*l/q))/(s^4 + (b*(I + m*l^2))*s^3/q - ((M +
m)*m*g*l)*s^2/q - b*m*g*l*s/q);
T2 = feedback(1,P_pend*C)*P_cart;
T2 = minreal(T2);
t = 0:0.01:10;
impulse(T2, t), grid
title({'Response of Cart Position to an Impulse Disturbance';'under Closed-loop Control'});

The command minreal effectively cancels out all common poles and zeros in the closed-loop
transfer function. This gives the impulse function better numerical properties. As you can see,
the cart moves in the negative direction and stabilizes at about -0.14 meters. This design
might work pretty well for the actual controller, assuming that the cart had that much room to
move. Keep in mind that this was pure luck. We did not design our controller to stabilize the
cart's position, the fact that we have is a fortunate side effect.
2.4 STATE SPACE
INVERTED PENDULUM: STATE SPACE METHOD FOR CONTROLLER DESIGN
Contents
 Open-loop poles
 Linear Quadratic Regulation (LQR)
 Adding precompensation
 Observer-based control
From the main problem, the dynamic equations of the inverted pendulum system in state-
space form are the following:

(1)

(2)
To see how this problem was originally set up and the system equations were derived, consult
the Inverted Pendulum: System Modeling page. For this problem the outputs are the cart's
displacement ( in meters) and the pendulum angle ( in radians) where represents the
deviation of the pedulum's position from equilibrium, that is, .
The design criteria for this system for a 0.2-m step in desired cart position are as follows:
 Settling time for and of less than 5 seconds
 Rise time for of less than 0.5 seconds
 Pendulum angle never more than 20 degrees (0.35 radians) from the vertical
 Steady-state error of less than 2% for and
As you may have noticed if you went through some of the other inverted pendulum examples,
the design criteria for this example are different. In the other examples we were attemping to
keep the pendulum vertical in response to an impulsive disturbance force applied to the cart.
We did not attempt to control the cart's position. In this example, we are attempting to keep
the pendulum vertical while controlling the cart's position to move 0.2 meters to the right. A
state-space design approach is well suited to the control of multiple outputs as we have here.
This problem can be solved using full-state feedback. The schematic of this type of control
system is shown below where is a matrix of control gains. Note that here we feedback all
of the system's states, rather than using the system's outputs for feedback.

Open-loop poles
In this problem, represents the step command of the cart's position. The 4 states represent
the position and velocity of the cart and the angle and angular velocity of the pendulum. The
output contains both the position of the cart and the angle of the pendulum. We want to
design a controller so that when a step reference is given to the system, the pendulum should
be displaced, but eventually return to zero (i.e. vertical) and the cart should move to its new
commanded position. To view the system's open-loop response please refer to the Inverted
Pendulum: System Analysis page.
The first step in designing a full-state feedback controller is to determine the open-loop poles
of the system. Enter the following lines of code into an m-file. After execution in the
MATLAB command window, the output will list the open-loop poles (eigenvalues of ) as
shown below.
M = 0.5;
m = 0.2;
b = 0.1;
I = 0.006;
g = 9.8;
l = 0.3;

p = I*(M+m)+M*m*l^2; %denominator for the A and B matrices

A = [0 1 0 0;
0 -(I+m*l^2)*b/p (m^2*g*l^2)/p 0;
0 0 0 1;
0 -(m*l*b)/p m*g*l*(M+m)/p 0];
B=[ 0;
(I+m*l^2)/p;
0;
m*l/p];
C = [1 0 0 0;
0 0 1 0];
D = [0;
0];

states = {'x' 'x_dot' 'phi' 'phi_dot'};

inputs = {'u'};
outputs = {'x'; 'phi'};

sys_ss = ss(A,B,C,D,'statename',states,'inputname',inputs,'outputname',outputs);

poles = eig(A)
poles =
0
-0.1428
-5.6041
5.5651
As you can see, there is one right-half plane pole at 5.5651. This should confirm your
intuition that the system is unstable in open loop.
Linear Quadratic Regulation (LQR)
The next step in the design process is to find the vector of state-feedback control gains
assuming that we have access (i.e. can measure) all four of the state variables. This can be
accomplished in a number of ways. If you know the desired closed-loop pole locations, you
can use the MATLAB commands place or acker. Another option is to use the lqr command
which returns the optimal controller gain assuming a linear plant, quadratic cost function, and
reference equal to zero (consult your textbook for more details).
Before we design our controller, we will first verify that the system is controllable.
Satisfaction of this property means that we can drive the state of the system anywhere we like
in finite time (under the physical constraints of the system). For the system to be completely
state controllable, the controllability matrix must have rank where the rank of a matrix is
the number of linearly independent rows (or columns). The controllability matrix of the
system takes the form shown below. The number corresponds to the number of state
variables of the system. Adding additional terms to the controllability matrix with higher
powers of the matrix will not increase the rank of the controllability matrix since these
additional terms will just be linear combinations of the earlier terms.

(3)
Since our controllability matrix is 4x4, the rank of the matrix must be 4. We will use the
MATLAB command ctrb to generate the controllability matrix and the MATLAB
command rank to test the rank of the matrix. Adding the following additional commands to
your m-file and running in the MATLAB command window will produce the following
output.
co = ctrb(sys_ss);
controllability = rank(co)
controllability =
4
Therefore, we have verified that our system is controllable and thus we should be able to
design a controller that achieves the given requirements. Specifically, we will use the linear
quadratic regulation method for determining our state-feedback control gain matrix . The
MATLAB function lqr allows you to choose two parameters, and , which will balance the
relative importance of the control effort ( ) and error (deviation from 0), respectively, in the
cost function that you are trying to optimize. The simplest case is to assume ,
and . The cost function corresponding to this and places equal importance on
the control and the state variables which are outputs (the pendulum's angle and the cart's
position). Essentially, the lqr method allows for the control of both outputs. In this case, it is
pretty easy to do. The controller can be tuned by changing the nonzero elements in the
matrix to achieve a desirable response. To observe the structure of , enter the following into
the MATLAB command window to see the output given below.
Q = C'*C
Q=
1 0 0 0
0 0 0 0
0 0 1 0
0 0 0 0
The element in the (1,1) position of represents the weight on the cart's position and the
element in the (3,3) position represents the weight on the pendulum's angle. The input
weighting will remain at 1. Ultimately what matters is the relative value of and , not
their absolute values. Now that we know how to interpret the matrix, we can experiment to
find the matrix that will give us a "good" controller. We will go ahead and find the
matrix and plot the response all in one step so that changes can be made in the control and
seen automatically in the response. Add the following commands to the end of your m-file
and run in the MATLAB command window to get the following value for and the response
plot shown below.
Q = C'*C;
R = 1;
K = lqr(A,B,Q,R)

Ac = [(A-B*K)];
Bc = [B];
Cc = [C];
Dc = [D];

states = {'x' 'x_dot' 'phi' 'phi_dot'};

inputs = {'r'};
outputs = {'x'; 'phi'};

sys_cl = ss(Ac,Bc,Cc,Dc,'statename',states,'inputname',inputs,'outputname',outputs);

t = 0:0.01:5;
r =0.2*ones(size(t));
[y,t,x]=lsim(sys_cl,r,t);
[AX,H1,H2] = plotyy(t,y(:,1),t,y(:,2),'plot');
set(get(AX(1),'Ylabel'),'String','cart position (m)')
set(get(AX(2),'Ylabel'),'String','pendulum angle (radians)')
title('Step Response with LQR Control')
K=
-1.0000 -1.6567 18.6854 3.4594

The curve in red represents the pendulum's angle in radians, and the curve in blue represents
the cart's position in meters. As you can see, this plot is not satisfactory. The pendulum and
cart's overshoot appear fine, but their settling times need improvement and the cart's rise time
needs to be reduced. As I'm sure you have noticed, the cart's final position is also not near the
desired location but has in fact moved in the opposite direction. This error will be dealt with
in the next section and right now we will focus on the settling and rise times. Go back to your
m-file and change the matrix to see if you can get a better response. You will find that
increasing the (1,1) and (3,3) elements makes the settling and rise times go down, and lowers
the angle the pendulum moves. In other words, you are putting more weight on the errors at
the cost of increased control effort . Modifying your m-file so that the (1,1) element of is
5000 and the (3,3) element is 100, will produce the following value of and the step
response shown below.
Q = C'*C;
Q(1,1) = 5000;
Q(3,3) = 100
R = 1;
K = lqr(A,B,Q,R)

Ac = [(A-B*K)];
Bc = [B];
Cc = [C];
Dc = [D];

states = {'x' 'x_dot' 'phi' 'phi_dot'};

inputs = {'r'};
outputs = {'x'; 'phi'};

sys_cl = ss(Ac,Bc,Cc,Dc,'statename',states,'inputname',inputs,'outputname',outputs);

t = 0:0.01:5;
r =0.2*ones(size(t));
[y,t,x]=lsim(sys_cl,r,t);
[AX,H1,H2] = plotyy(t,y(:,1),t,y(:,2),'plot');
set(get(AX(1),'Ylabel'),'String','cart position (m)')
set(get(AX(2),'Ylabel'),'String','pendulum angle (radians)')
title('Step Response with LQR Control')
Q=
5000 0 0 0
0 0 0 0
0 0 100 0
0 0 0 0
K=
 -70.7107 -37.8345 105.5298 20.9238

You may have noted that if you increased the values of the elements of even higher, you
could improve the response even more. The reason this weighting was chosen was because it
just satisfies the transient design requirements. Increasing the magnitude of more would
make the tracking error smaller, but would require greater control force . More control effort
generally corresponds to greater cost (more energy, larger actuator, etc.).
Adding precompensation
The controller we have designed so far meets our transient requirements, but now we must
address the steady-state error. In contrast to the other design methods, where we feedback the
output and compare it to the reference input to compute an error, with a full-state feedback
controller we are feeding back all of the states. We need to compute what the steady-state
value of the states should be, multiply that by the chosen gain , and use a new value as our
"reference" for computing the input. This can be done by adding a constant gain after the
reference. The schematic below shows this relationship:
We can find this factor by employing the used-defined function rscale.m as shown below.
The matrix is modified to reflect the fact that the reference is a command only on cart
position.
Cn = [1 0 0 0];
sys_ss = ss(A,B,Cn,0);
Nbar = rscale(sys_ss,K)
Nbar =
-70.7107
Note that the function rscale.m is not a standard function in MATLAB. You will have to
download it here and place it in your current directory. More information can be found
here, Extras: rscale.m. Now you can plot the step response by adding the above and following
lines of code to your m-file and re-running at the command line.
sys_cl = ss(Ac,Bc*Nbar,Cc,Dc,'statename',states,'inputname',inputs,'outputname',outputs);

t = 0:0.01:5;
r =0.2*ones(size(t));
[y,t,x]=lsim(sys_cl,r,t);
[AX,H1,H2] = plotyy(t,y(:,1),t,y(:,2),'plot');
set(get(AX(1),'Ylabel'),'String','cart position (m)')
set(get(AX(2),'Ylabel'),'String','pendulum angle (radians)')
title('Step Response with Precompensation and LQR Control')
Now, the steady-state error is within our limits, the rise and settle times are met, and the
pendulum's overshoot is within range of the design criteria.
Note that the precompensator employed above is calculated based on the model of the
plant and further that the precompensator is located outside of the feedback loop. Therefore,
if there are errors in the model (or unknown disturbances) the precompensator will not correct
for them and there will be steady-state error. You may recall that the addition of integral
control may also be used to eliminate steady-state error, even in the presence of model
uncertainty and step disturbances. For an example of how to implement integral control in the
state space setting, see the Motor Position: State-Space Methods example. The tradeoff with
using integral control is that the error must first develop before it can be corrected for,
therefore, the system may be slow to respond. The precompensator on the other hand is able
to anticipitate the steady-state offset using knowledge of the plant model. A useful technique
is to combine the precompensator with integral control to leverage the advantages of each
approach.
Observer-based control
The response achieved above is good, but was based on the assumption of full-state feedback,
which is not necessarily valid. To address the situation where not all state variables are
measured, a state estimator must be designed. A schematic of state-feedback control with a
full-state estimator is shown below, without the precompensator .
Before we design our estimator, we will first verify that our system is observable. The
property of observability determines whether or not based on the measured outputs of the
system we can estimate the state of the system. Similar to the process for verifying
controllability, a system is observable if its observability matrix is full rank. The observability
matrix is defined as follows.

(4)
We can employ the MATLAB command obsv to contruct the observability matrix and
the rank command to check its rank as shown below.
ob = obsv(sys_ss);
observability = rank(ob)
observability =
4
Since the observability matrix is 8x4 and has rank 4, it is full rank and our system is
observable. The observability matrix in this case is not square since our system has two
outputs. Note that if we could only measure the pendulum angle output, we would not be able
to estimate the full state of the system. This can be verified by the fact
that obsv(A,C(2,:)) produces an observability matrix that is not full rank.
Since we know that we can estimate our system state, we will now describe the process for
designing a state estimator. Based on the above diagram, the dynamics of the state estimate
are described by the following equation.

(5)
The spirit of this equation is similar to that of closed-loop control in that last term is a
correction based on feedback. Specifically, the last term corrects the state estimate based on
the difference between the actual output and the estimated output . Now let's look at the
dynamics of the error in the state estimate.

(6)
Therefore, the state estimate error dynamics are described by
(7)
and the error will approach zero ( will approach ) if the matrix is stable (has
negative eigenvalues). As is with the case for control, the speed of convergence depends on
the poles of the estimator (eigenvalues of ). Since we plan to use the state estimate as
the input to our controller, we would like the state estimate to converge faster than is desired
from our overall closed-loop system. That is, we would like the observer poles to be faster
than the controller poles. A common guideline is to make the estimator poles 4-10 times
faster than the slowest controller pole. Making the estimator poles too fast can be problematic
if the measurement is corrupted by noise or there are errors in the sensor measurement in
general.
Based on this logic, we must first find the controller poles. To do this, copy the following
code to the end of your m-file. If you employed the updated matrix, you should see the
following poles in the MATLAB command window.
poles = eig(Ac)
poles =
-8.4910 + 7.9283i
-8.4910 - 7.9283i
-4.7592 + 0.8309i
-4.7592 - 0.8309i
The slowest poles have real part equal to -4.7592, therefore, we will place our estimator poles
at -40. Since the closed-loop estimator dynamics are determined by a matrix ( ) that
has a similar form to the matrix that determines the dynamics of the state-feedback system (
), we can use the same commands for finding the estimator gain as we can for
finding the state-feedback gain . Specifically, since taking the transpose of leaves
the eigenvalues unchanged and produces a result that exactly matches the form
of , we can use the acker or place commands. Recalling that the place command
cannot place poles of multiplicity greater than one, we will place the observer poles as
follows. Add the following commands to your m-file to calculate the matrix and generate
the output shown below.
P = [-40 -41 -42 -43];
L = place(A',C',P)'
L=
1.0e+03 *
0.0826 -0.0010
1.6992 -0.0402
-0.0014 0.0832
-0.0762 1.7604
We are using both outputs (the angle of the pendulum and the position of the cart) to design
the observer.
Now we will combine our state-feedback controller from before with our state estimator to
get the full compensator. The resulting closed-loop system is described by the following
matrix equations.

(8)

(9)
The closed-loop system described above can be implemented in MATLAB by adding the
following commands to the end of your m-file. After running the m-file the step response
shown will be generated.
Ace = [(A-B*K) (B*K);
zeros(size(A)) (A-L*C)];
Bce = [B*Nbar;
zeros(size(B))];
Cce = [Cc zeros(size(Cc))];
Dce = [0;0];
states = {'x' 'x_dot' 'phi' 'phi_dot' 'e1' 'e2' 'e3' 'e4'};
inputs = {'r'};
outputs = {'x'; 'phi'};

sys_est_cl = ss(Ace,Bce,Cce,Dce,'statename',states,'inputname',inputs,'outputname',outputs);

t = 0:0.01:5;
r = 0.2*ones(size(t));
[y,t,x]=lsim(sys_est_cl,r,t);
[AX,H1,H2] = plotyy(t,y(:,1),t,y(:,2),'plot');
set(get(AX(1),'Ylabel'),'String','cart position (m)')
set(get(AX(2),'Ylabel'),'String','pendulum angle (radians)')
title('Step Response with Observer-Based State-Feedback Control')

This response is almost identical to the response achieved when it was assumed that we had
full access to the state variables. This is because the observer poles are fast, and because the
model we assumed for the observer is identical to the model of the actual plant (including the
same initial conditions). Therefore, all of the design requirements have been met with the
minimal control effort expended. No further iteration is needed.
This example demonstrates that it is much easier to control multi-input, multi-output systems
with the state-space method than with the other methods we have presented.

2.5 DIGITAL
Inverted Pendulum: Digital Controller Design

Contents
 Discrete state-space
 Controllability and observability
 Control design via pole placement
 Precompensator design
 Observer design
In this digital control version of the inverted pendulum problem, we will use the state-space
method to design the digital controller. If you refer to the Inverted Pendulum: System
Modeling page, the linearized state-space equations were derived as:

(1)
(2)
where:
(M) mass of the cart 0.5 kg
(m) mass of the pendulum 0.2 kg
(b) coefficient of friction for cart 0.1 N/m/sec
(l) length to pendulum center of mass 0.3 m
(I) mass moment of inertia of the pendulum 0.006 kg.m^2
(F) force applied to the cart
(x) cart position coordinate
(theta) pendulum angle from vertical (down)
For this problem the outputs are the cart's displacement ( in meters) and the pendulum angle
( in radians) where represents the deviation of the pedulum's position from equilibrium,
that is, = + .
The design criteria for this system for a 0.2-m step in desired cart position are as follows:
 Settling time for and theta of less than 5 seconds
 Rise time for of less than 0.5 seconds
 Pendulum angle never more than 20 degrees (0.35 radians) from the vertical
 Steady-state error of less than 2% for and
Discrete state-space
Our first step in designing a digital controller is to convert the above continuous state-space
equations to a discrete form. We will accomplish this employing the MATLAB function c2d.
This function requires that we specify three arguments: a continuous system model, the
sampling time (Ts in sec/sample), and the 'method'. You should already be familiar with how
to construct a state-space system from , , , and matrices.
In choosing a sample time, note that it is desired that the sampling frequency be fast
compared to the dynamics of the system. One measure of a system's "speed" is its closed-loop
bandwidth. A good rule of thumb is that the sampling frequency be at least 30 times larger
than the closed-loop bandwidth frequency which can be determined from the closed-loop
Bode plot.
Assuming that the closed-loop bandwidth frequencies are around 1 rad/sec for both the cart
and the pendulum, let the sampling time be 1/100 sec/sample. The discretization method we
will use is the zero-order hold ('zoh'). For further details, refer to the Introduction: Digital
Controller Design page. Now we are ready to use c2d function. Enter the following
commands into an m-file. Running this m-file in the MATLAB command window gives you
the following four matrices representing the discrete time state-space model.
M = 0.5;
m = 0.2;
b = 0.1;
I = 0.006;
g = 9.8;
l = 0.3;

p = I*(M+m)+M*m*l^2; %denominator for the A and B matrices

A = [0 1 0 0;
0 -(I+m*l^2)*b/p (m^2*g*l^2)/p 0;
0 0 0 1;
0 -(m*l*b)/p m*g*l*(M+m)/p 0];
B=[ 0;
(I+m*l^2)/p;
0;
m*l/p];
C = [1 0 0 0;
0 0 1 0];
D = [0;
0];

states = {'x' 'x_dot' 'phi' 'phi_dot'};

inputs = {'u'};
outputs = {'x'; 'phi'};
sys_ss = ss(A,B,C,D,'statename',states,'inputname',inputs,'outputname',outputs);

Ts = 1/100;

sys_d = c2d(sys_ss,Ts,'zoh')
sys_d =

A=
x x_dot phi phi_dot
x 1 0.009991 0.0001336 4.453e-07
x_dot 0 0.9982 0.02672 0.0001336
phi 0 -2.272e-05 1.002 0.01001
phi_dot 0 -0.004544 0.3119 1.002

B=
u
x 9.086e-05
x_dot 0.01817
phi 0.0002272
phi_dot 0.04544

C=
x x_dot phi phi_dot
x 1 0 0 0
phi 0 0 1 0

D=
u
x 0
 phi 0

Sample time: 0.01 seconds

Discrete-time state-space model.

Now we have obtained the discrete state-space model of the form:

(3)

(4)
Controllability and observability
The next step is to check the controllability and the observability of the system. For the
system to be completely state controllable, the controllability matrix

(5)
must have the rank of n. The rank of the matrix is the number of linearly independent rows
(or columns). In the same token, for the system to be completely state observable, the
observability matrix

(6)
must also have the rank of n. These tests for controllability and observability are identical to
the situation of continuous control except that now the state space model is discrete.
Since the number of state variables in our system is 4, the rank of both matrices must be 4.
The function rank can give you the rank of each matrix. Adding the following commands to
your m-file and running in the MATLAB command window will generate the results shown
below.
co = ctrb(sys_d);
ob = obsv(sys_d);

controllability = rank(co)
observability = rank(ob)
controllability =
4
observability =
4
This proves that our discrete system is both completely state controllable and completely
state observable.
Control design via pole placement
The schematic of a full-state feedback control system is shown below.

The next step is to assume that all four states are measurable and design the control gain
matrix . If you refer to the continuous Inverted Pendulum: State-Space Methods for
Controller Design page the Linear Quadratic Regulator (LQR) method was used to find
the control gain matrix . In this digital version, we will use the same LQR method. This
method allows you to find the control gain that results in the optimal balance between system
errors and control effort. Please consult your control textbook for details. To use this LQR
method, we need to specify two parameters, the performance index matrix and the state-
cost matrix . For simplicity, we will initially choose the performance index matrix equal
to 1, and the state-cost matrix equal to . The relative weightings of these two matrices
will then be tuned by trial and error. The state-cost matrix has the following structure.

(7)
The element in the (1,1) position of represents the weight on the cart's position and the
element in the (3,3) position represents the weight on the pendulum's angle.
Now we are ready to find the control gain matrix and observe the resulting closed-loop
response of the system. Since we are designing a digital controller, we will specificially
employ the MATLAB function dlqr. Add the following commands to a your m-file and run it
in the MATLAB command window. Note that in the following we are overwriting the values
of the state-space matrices , , , and with their discrete-time equivalents using the
model derived with the c2d command above.
A = sys_d.a;
B = sys_d.b;
C = sys_d.c;
D = sys_d.d;
Q = C'*C
R = 1;
[K] = dlqr(A,B,Q,R)

Ac = [(A-B*K)];
Bc = [B];
Cc = [C];
Dc = [D];

states = {'x' 'x_dot' 'phi' 'phi_dot'};

inputs = {'r'};
outputs = {'x'; 'phi'};

sys_cl = ss(Ac,Bc,Cc,Dc,Ts,'statename',states,'inputname',inputs,'outputname',outputs);

t = 0:0.01:5;
r =0.2*ones(size(t));
[y,t,x]=lsim(sys_cl,r,t);
[AX,H1,H2] = plotyy(t,y(:,1),t,y(:,2),'plot');
set(get(AX(1),'Ylabel'),'String','cart position (m)')
set(get(AX(2),'Ylabel'),'String','pendulum angle (radians)')
title('Step Response with Digital LQR Control')
Q=
1 0 0 0
0 0 0 0
0 0 1 0
0 0 0 0
K=
-0.9384 -1.5656 18.0351 3.3368

The curve in red represents the pendulum's angle in radians, and the curve in blue represents
the cart's position in meters. As you can see, this plot is not satisfactory. The pendulum and
cart's overshoot appear fine, but their settling times need improvement and the cart's rise time
needs to be reduced. As I'm sure you have noticed, the cart's final position is also not near the
desired location but has in fact moved in the opposite direction. This error will be dealt with
in the next section and right now we will focus on the settling and rise times. Go back to your
m-file and change the matrix to see if you can get a better response. You will find that
increasing the (1,1) and (3,3) elements makes the settling and rise times go down, and lowers
the angle the pendulum moves. In other words, you are putting more weight on the errors at
the cost of increased control effort . Modifying your m-file so that the (1,1) element of is
5000 and the (3,3) element is 100, will produce the following value of and the step
response shown below.
A = sys_d.a;
B = sys_d.b;
C = sys_d.c;
D = sys_d.d;
Q = C'*C;
Q(1,1) = 5000;
Q(3,3) = 100
R = 1;
[K] = dlqr(A,B,Q,R)

Ac = [(A-B*K)];
Bc = [B];
Cc = [C];
Dc = [D];

states = {'x' 'x_dot' 'phi' 'phi_dot'};

inputs = {'r'};
outputs = {'x'; 'phi'};

sys_cl = ss(Ac,Bc,Cc,Dc,Ts,'statename',states,'inputname',inputs,'outputname',outputs);

t = 0:0.01:5;
r =0.2*ones(size(t));
[y,t,x]=lsim(sys_cl,r,t);
[AX,H1,H2] = plotyy(t,y(:,1),t,y(:,2),'plot');
set(get(AX(1),'Ylabel'),'String','cart position (m)')
set(get(AX(2),'Ylabel'),'String','pendulum angle (radians)')
title('Step Response with Digital LQR Control')
Q=
5000 0 0 0
0 0 0 0
0 0 100 0
0 0 0 0
K=
-61.9933 -33.5040 95.0597 18.8300

From this plot, we see that all design requirements are satisfied except the steady-state error
of the cart position . We can easily correct this by introducing a feedforward scaling
factor .
Precompensator design
Unlike other design methods, the full-state feedback system does not compare the output
directly to the reference, rather, it compares the state vector multiplied by the control matrix (
) to the reference (see the schematic shown above). Thus, we should not expect the output
to converge to the commanded reference. To obtain the desired output, we need to scale the
reference input so that the output equals the reference. This can be easily done by introducing
a feedforward scaling factor . The basic full state-feedback schematic with scaling factor is
shown below.
Unfortunately, we cannot use our user-defined function rscale to find because this function
was defined for continuous-time single-output systems. We can, however, find the scaling
factor by trial and error. After several trials, equal to -61.55 provided a satisfactory
response. Adding the following commands to your m-file and running in the command
window will generate the response shown below.
Nbar = -61.55;
sys_cl = ss(Ac,Bc*Nbar,Cc,Dc,Ts,'statename',states,'inputname',inputs,'outputname',outputs);

t = 0:0.01:5;
r =0.2*ones(size(t));
[y,t,x]=lsim(sys_cl,r,t);
[AX,H1,H2] = plotyy(t,y(:,1),t,y(:,2),'plot');
set(get(AX(1),'Ylabel'),'String','cart position (m)')
set(get(AX(2),'Ylabel'),'String','pendulum angle (radians)')
title('Step Response with Digital LQR Control and Precompensation')
Notice that the steady-state error of the cart's position has been eliminated. Now we have
designed a system that satisfies all of the design requirements. Note, however, that the scaling
factor was designed based on a model of the system. If our model is in error or there are
unknown disturbances, then the steady-state error will no longer be driven to zero.
Observer design
The above response satisfies all of the design requirements, however, it was found assuming
all state variables of the system are measurable. This assumption may not be valid for all
systems. In this section, we develop a technique for estimating the state of the system based
on the measured outputs and a model of the plant. The object that estimates the state of
system is called an observer. Thus, in this section we will design a full-order state observer
to estimate all of the system's state variables, including those that are measured. For further
explanation on how an observer works, please consult your control textbook.
A basic schematic of the observer-based state-feedback system is shown below.
Designing the observer equates to finding the observer gain matrix . To accomplish this, we
need to first determine the closed-loop poles of the system without the observer (the
eigenvalues of ). This can be achieved using the MATLAB command eig as shown
below.
poles = eig(A-B*K)
poles =
0.9157 + 0.0728i
0.9157 - 0.0728i
0.9535 + 0.0079i
0.9535 - 0.0079i
Since the observer is attempting to estimate the values of state variables which are themselves
changing, it is desired that the dynamics of the observer be significantly faster than the
dynamics of the closed-loop system without the observer. A common guideline is to make the
estimator poles (eigenvalues of ) 4-10 times faster than the slowest controller pole
(eigenvalue of ). Making the estimator poles too fast can be problematic if the
measurement is corrupted by noise or there are errors in the sensor measurement in general.
Based on the poles found above, we will place the observer poles at [-0.2 -0.21 -0.22 -0.23].
These poles can be modified later, if necessary. We will use the MATLAB function place to
find the matrix. Add the following code to your m-file and re-run in the command window
to generate the observer gain matrix shown below.
P = [-0.2 -0.21 -0.22 -0.23];
L = place(A',C',P)'
L=
2.4308 -0.0104
147.6324 -1.2418
-0.0131 2.4305
-1.8079 147.9057
Now we will obtain the overall system response including the observer. Add the following
commands to your m-file and run to generate the response shown below.
Ace = [(A-B*K) (B*K);
zeros(size(A)) (A-L*C)];
Bce = [B*Nbar;
zeros(size(B))];
Cce = [Cc zeros(size(Cc))];
Dce = [0;0];

states = {'x' 'x_dot' 'phi' 'phi_dot' 'e1' 'e2' 'e3' 'e4'};

inputs = {'r'};
outputs = {'x'; 'phi'};

sys_est_cl =
ss(Ace,Bce,Cce,Dce,Ts,'statename',states,'inputname',inputs,'outputname',outputs);

t = 0:0.01:5;
r = 0.2*ones(size(t));
[y,t,x]=lsim(sys_est_cl,r,t);
[AX,H1,H2] = plotyy(t,y(:,1),t,y(:,2),'plot');
set(get(AX(1),'Ylabel'),'String','cart position (m)')
set(get(AX(2),'Ylabel'),'String','pendulum angle (radians)')
title('Step Response with Digital Observer-Based State-Feedback Control')

This response is almost identical to the response achieved when it was assumed that we had
full access to the state variables. This is because the observer poles are fast, and because the
model we assumed for the observer is identical to the model of the actual plant (including the
same initial conditions). Therefore, all of the design requirements have been met with the
minimal control effort expended. No further iteration is needed.
 3 SIMULINK
3.1 MODELLING
Contents
 Physical setup and system equations
 Building the nonlinear model with Simulink
 Building the nonlinear model with Simscape
 Generating the open-loop response
 Extracting a linear model from the simulation
In this page we outline how to build a model of our inverted pendulum system for the
purposes of simulation using Simulink and its add-ons. A great advantage of simulation, as
will be demonstrated in this example, is that it can generate numerical solutions to nonlinear
equations for which closed-form solutions cannot be generated. The nonlinear simulation can
then be employed to test the validity of a linearized version of the model. The simulation
model can also be used to evaluate the performance of the control scheme designed based on
the linearized model.
Physical setup and system equations
In this example we will consider a two-dimensional version of the inverted pendulum system
with cart where the pendulum is constrained to move in the vertical plane shown in the figure
below. For this system, the control input is the force that moves the cart horizontally and
the outputs are the angular position of the pendulum and the horizontal position of the cart
.
For this example, let's assume the following quantities:
(M) mass of the cart 0.5 kg
(m) mass of the pendulum 0.2 kg
(b) coefficient of friction for cart 0.1 N/m/sec
(l) length to pendulum center of mass 0.3 m
(I) mass moment of inertia of the pendulum 0.006 kg.m^2
(F) force applied to the cart
(x) cart position coordinate
(theta) pendulum angle from vertical (down)
Below are the two free-body diagrams of the system.

This system is challenging to model in Simulink because of the physical constraint (the pin
joint) between the cart and pendulum which reduces the degrees of freedom in the system.
Both the cart and the pendulum have one degree of freedom ( and , respectively). We will
generate the differential equations for these degrees of freedom from first principles
employing Newton's second law ( ) as shown below.

(1)

(2)
It is necessary, however, to include the interaction forces and between the cart and the
pendulum in order to fully model the system's dynamics. The inclusion of these forces
requires modeling the - and -components of the translation of the pendulum's center of
mass in addition to its rotational dynamics. In the Inverted Pendulum: System
Modeling tutorial, the interaction forces and were solved for algebraically.
In general, we would like to exploit the modeling power of Simulink to take care of the
algebra for us. Therefore, we will model the additional - and -component equations for the
pendulum as shown below.

(3)
(4)

(5)

(6)
However, the position coordinates and are exact functions of . Therefore, we can
represent their derivatives in terms of the derivatives of . First addressing the -component
equations we arrive at the following.
(7)

(8)

(9)
Then addressing the -component equations gives us the following.
(10)

(11)

(12)
These expressions can then be substituted into the expressions for and from above as
follows.

(13)

(14)
We can now represent these equations within Simulink. Simulink can work directly with
nonlinear equations, so it is unnecessary to linearize these equations as was done in
the Inverted Pendulum: System Modeling page.
Building the nonlinear model with Simulink
We can build the inverted pendulum model in Simulink employing the equations derived
above by following the steps given below.
 Begin by typing simulink into the MATLAB command window to open the Simulink
environment. Then open a new model window in Simulink by choosing New >
 Simulink > Blank Model of the open Simulink Start Page window or by
pressing Ctrl-N.
 Insert four Fcn Blocks from the Simulink/User-Defined Functions library. We will
build the equations for , , , and employing these blocks.
 Change the label of each Fcn block to match its associated function.
 Insert four Integrator blocks from the Simulink/Continuous library. The output of each
Integrator block is going to be a state variable of the system: , , , and .
 Double-click on each Integrator block to add the State Name: of the associated state
variable. See the following figure for an example. Also change the Initial
condition: for (pendulum angle) to "pi" to represent that the pendulum begins
pointing straight up.

 Insert four Multiplexer (Mux) blocks from the Simulink/Signal Routing library, one
for each Fcn block.
 Insert two Out1 blocks and one In1 block from the Simulink/Sinks and
Simulink/Sources libraries, respectively. Then double-click on the labels for the
blocks to change their names. The two outputs are for the "Position" of the cart and
the "Angle" of the pendulum, while the one input is for the "Force" applied to the cart.
  Connect each output of the Mux blocks to the input of the corresponding Fcn block.
 Connect the output of the and Fcn blocks to two consecutive integrators to
generate the cart's position and the pendulum's angle. Your current model should now
appear as follows.

Now we will enter each of the four equations (1), (2), (13), and (14) into a Fcn block. Let's
start with equation (1) which is repeated below.

(15)
 This equation requires three inputs: , , and . Double-click
on the corresponding Mux block and change the Number of inputs: to "3".
 Connect these three inputs to this Mux block in the order prescribed in the previous
step.
 Double-click on the first Fcn block and enter the equation for xddot as shown below.
Now, let's enter equation (2) which is repeated below.

(16)
 This equation also requires three inputs: , , and .
 Enter the above equation into the Fcn block, change the number of inputs of the Mux
block, and connect the correct signals to the Mux block in the correct order.
 Repeat this process for equations (13) and (14) repeated below.

(17)

(18)
When all of these steps are completed, the resulting model should appear as follows.
In order to save all of these components as a single subsystem block, first select all of the
blocks, then select Create Subsystem from Selection from the menu after right-clicking on
the selected portion. Your model should appear as follows. You can also download the file for
this system by right-clicking here and then selecting Save link as ....
 Building the nonlinear model with Simscape
In this section, we alternatively show how to build the inverted pendulum model using the
physical modeling blocks of the Simscape extension to Simulink. The blocks in the Simscape
library represent actual physical components; therefore, complex multi-body dynamic models
can be built without the need to build mathematical equations from physical principles as was
done above by applying Newton's laws.
Open a new Simulink model and follow the steps below to create the inverted pendulum
model in Simscape. In order to orient oneself, we will assume a coordinate system where the
cart moves in the -direction (positive to the right) and the positive -direction is directed up.
Following standard convention, the positive -direction is then pointed out of the plane of
motion.
 Insert a Body block from the Simscape/Multibody/First Generation(1G)/Bodies
library to represent the cart. Following the system parameters given at the top of this
page, double-click on the block and set the Mass: to "0.5" with units of kg. The Body
block by default includes two ports. Since we need ports to define where the
pendulum is connected to the cart and where the external force and the frictional force
are applied, a third port must be added. This can be accomplished from the button on
the right-side of the Position tab. Since the cart will only move in one dimension, the
two forces must be co-located along that dimension (the -direction). Since we are in
essence modeling the cart as a point mass that can only translate, you do not have to
 change any of the other default parameters. However, we plan to use Simscape to
animate the motion of the system and hence will create additional ports to define the
four corners of the cart (2-dimensional only) with respect to its center of gravity (CG).
The following figure shows a possible definition of the cart body.

 Insert a second Body block to represent the pendulum. Double-click on the block and
set the Mass: to "0.2" with units of kg. Since the pendulum can only rotate about the
-axis, the inertia associated with that principle direction is the only one that needs to
be defined. For simplicity, define the Inertia: equal to "0.006*eye(3)" with units
of kg*m^2. Since we are modeling the pendulum as a rigid body that has size as well
as mass, the body can rotate and it is important to define the location of the
pendulum's attachment to the cart and its CG correctly. Specifically, define the point
of attachment CS1 to have a position "[0 0 0]" and an origin that is Adjoining and
define the CG to be 0.3 meters away from the attachment CS1 (as defined above).
Also define the four corners of the pendulum. Make sure to show the port defining the
attachment point. Under the Visualization tab, you can also change the pendulum's
color to make it stand out from the cart.
 Next add a Revolute block from the Simscape/Multibody/First Generation(1G)/Jjoints
library to define the joint connecting the pendulum to the cart. By default, the joint
will be defined to rotate about the -axis which matches the situation we are
modeling. Connect the Body block corresponding to the cart to the base port (B) of
the joint and the Body block corresponding to the pendulum to the follower port (F)
of the joint. Double-click on the Revolute block and set the Number of sensor /
actuator ports: to "2".
 Then add a Joint Initial Condition block and a Joint Sensor block from the
Simscape/Multibody/First Generation(1G)/Sensors & Actuators library and connect
these blocks to the Revolute block. Double-click on the Joint Initial Condition block
and check the Enable box. We can use the default values for initial position and
velocity of the joint. Employing an initial position of 0 degrees corresponds to the
pendulum being pointed vertically upward based on the definition of the pendulum
body above. This isn't consistent with the original definition of , but it will make the
response results consistent with those generated from the linearized model in the other
pages of this example. Next double-click on the Joint Sensor block and change the
units on the Angle measurement to rad. Angular position is the only measurement that
is needed for this joint, the other boxes may remain unchecked.
 Add two Prismatic blocks from the Simscape/Multibody/First Generation(1G)/Joints
library to define the translational degree of freedom of the cart and the application of
the forces to the cart. Since the cart is technically a point mass we need only one
Prismatic block, but by employing two we can apply the forces at different locations.
 Double-click on each Prismatic block and change the Axis of Action to "[1 0 0]" to
reflect the fact that the two forces act in the -direction. Then connect the follower
port (F) of each block to the ports for the applied force (CS1) and frictional force
(CS2) on the Body block representing the cart.
 Next add two Ground blocks from the Simscape/Multibody/First
Generation(1G)/Bodies library to define the base for the motion of the cart.
Specifically, connect the output of each ground block to the base port (B) of each
Prismatic block.
 For one of the Ground blocks you just created, double-click on the block and check
the Show Machine Environment port box. Then add a Machine Environment block
from the Simscape/Multibody/First Generation(1G)/Bodies library and connect it to
the Ground block for which you just added the port. The Machine Environment block
allows us to define the gravitational force in the simulation. In this case the default
direction (negative -direction) and magnitude ("9.81") for units of m/s^2 are correct.
This block also allows us to define the parameters for visualization and the numerical
solver. The default parameters are fine for this example.
 Next add two Joint Actuator blocks and one Joint Sensor block from
Simscape/Multibody/First Generation(1G)/Sensors & Actuators library. The Joint
Actuator blocks will be employed for generating the external applied force and the
frictional force, while the Joint Sensor block will sense the motion of the cart. Note,
there is also a Translational Friction block that is available, but we will calculate the
frictional force ourselves since we are employing only a simple viscous model.
Double-click on one of the Prismatic blocks and set the Number of sensor / actuator
ports: to "1" (for the force actuator). For the other Prismatic block, set the Number of
sensor / actuator ports: to "2" (one for the force actuator and one for the cart sensor).
Then connect the Joint Actuator and Joint Sensor blocks as described. The default
values for the Joint Actuator blocks are sufficient for this case, but we must change
the Joint Sensor block to output position and velocity since the velocity is needed for
calculating the frictional force. Double-click on the Joint Sensor block and check the
box for Velocity while leaving the box for Position checked. The default metric units
do not need to be changed. Also uncheck the box for Output selected parameters as
one signal.
 Add a Gain block from the Simulink/Math Operations library to represent the viscous
friction coefficient . Set the Gain to "0.1" as defined at the top of this page and
connect the input to the velocity output of the Joint Sensor block for the cart and
connect the output of the gain to the Joint Actuator for the frictional force.
 Next, add two Out1 blocks and one In1 block from the Simulink/Ports & Subsystems
library. Connect the Out1 blocks to the remaining Joint Sensor block outputs and the
In1 block to the remaining Joint Actuator input.
  Finally, connect and label the components as shown in the following figure. You can
rotate a block in a similar manner to the way you flipped blocks, that is, by right-
clicking on the block then selecting Rotate Block from the Rotate & Flip menu.

You can also save this model as a single subsystem block as described in the previous section.
You can change the color of the subsystem by right-clicking on the block and
choosing Background Color from the resulting menu. You can download the complete
model file by right-clicking here, but note that you will need the Simscape addition to
Simulink in order to run the file. We use this model in the Inverted Pendulum: Simulink
Controller Design page.
Generating the open-loop response
We will now simulate the response of the inverted pendulum system to an impulsive force
applied to the cart. This simulation requires an impulse input. Since there is no such block in
the Simulink library, we will use the Pulse Generator block to approximate a unit impulse
input. We could use either of the models we generated above, however, we will use the
Simscape model in this case because it will allow us to visualize the motion of the inverted
pendulum system. Follow the steps given below.
 Open the inverted pendulum simscape model generated above.
 Add a Pulse Generator block from the Simulink/Sources library. Double-click on the
block and change the parameters as shown below. In particular, change the Period: to
 "10". Since we will run our simulation for 10 seconds, this will ensure that only a
single "pulse" is generated. Also change the Amplitude to "1000" and the Pulse
Width (% of period): to "0.01". Together, these settings generate a pulse that
approximates a unit impulse in that the magnitude of the input is very large for a very
short period of time and the area of the pulse equals 1.

 Add a Scope block from the Simulink/Sinks library.

 In order display two inputs on the scope, right-click on the Scope block, choose
the Signals & Ports and change the Number of Input Ports to "2".
Connect the blocks and label the signals connected to the Scope block as shown.
Save this system as Pend_Openloop.slx, or download by right-clicking here and
selecting Save link as ....
Before we start the simulation, we would like to enable the visualization of the inverted
pendulum system. From the menus at the top of the model window choose Simulation >
Model Configuration Parameters. Then from the directory on the left-side of the window
choose Simscape Multibody 1G. Then check the box for Show animation during
simulation as shown in the figure below.
Now, start the simulation (select Run from the Simulation menu or enter Ctrl-T). As the
simulation runs, an animation of the inverted pendulum like the one shown below will
visualize the system's resulting motion.
Then open the Scope. You will see the following output for the pendulum angle and the cart
position.
Notice that the pendulum repeatedly swings through full revolutions where the angle rolls
over at radians. Furthermore, the cart's position grows unbounded, but oscillates under the
influence of the swinging pendulum. These results differ quite a bit from the results of the
open-loop simulation shown in the Inverted Pendulum: System Analysis page. This is due of
course to the fact that this simulation employed a fully nonlinear model, while the previous
analysis had relied on a linear approximation of the inverted pendulum model. In order to
compare the results of the simulation model more directly to the prior results, we will extract
a linear model from our simulation model.
Extracting a linear model from the simulation
Aside from comparing our simulation model to our prior results, it may also be desirable to
extract a linear model for the purposes of analysis and design. Much of the analytical
techniques that are commonly applied to the analysis of dynamic systems and the design of
their associated control can only be applied to linear models. Therefore, it may be desirable to
extract an approximate linear model from the nonlinear simulation model. We will
accomplish this from within Simulink.
 To begin, open either of the Simulink models generated above, Pend_Model.slx or
Pend_Model_Simscape.slx.
 If you generated your simulation model using variables, it is necessary to define the
physical constants in the MATLAB workspace before performing the linearization.
This can be accomplished by entering the following commands in the MATLAB
command window.
M = 0.5;
m = 0.2;
b = 0.1;
I = 0.006;
g = 9.8;
l = 0.3;
 Next choose from the menus at the top of the model window Analysis > Control
Design > Linear Analysis. This will cause the Linear Analysis Tool window to open.
 In order to perform our linearization, we need to first identify the inputs and outputs
for the model and the operating point that we wish to perform the linearization about.
First right-click on the signal representing the Force input in the Simulink/Simscape
model. Then choose Linear Analysis Points > Open-loop Input from the resulting
menu. Similarly, right-click on each of the two output signals of the model (pendulum
angle and cart position) and select Linear Analysis Points > Open-loop Output from
the resulting menu in each case. The resulting inputs and outputs should now be
identified on your model by arrow symbols as shown in the figure below.
 Next we need to identify the operating point to be linearized about. From
the Operating Point: menu choose Trim Model as shown in the figure below. This
will open the Trim the model window. Within this window, select the Start
trimming button indicated by the green triangle. This will create the operating
point op_trim1.
 Since we wish to examine the impulse response of this system, return to the LINEAR
ANALYSIS tab and choose Impulse as shown in the figure below.

 Finally, choose op_trim1 from the Operating Point: drop-down menu and click
the Impulse button indicated by the small green triangle. This automatically generates
an impulse response plot and the linearized model linsys1.
 In order to compare the results to those plots generated in the Inverted Pendulum:
System Analysis page, it is necessary to change the -axis scaling. This can be
achieved from by choosing Properties from the right-click menu. The resulting
window should then appear as follows, where the top plot is response of the pendulum
angle and the bottom plot is the response of the cart position.
These plots are very similar, though not exactly the same, as those generated in the Inverted
Pendulum: System Analysis page.
We can also export the resulting linearized model into the MATLAB workspace for further
analysis and design. This can be accomplished by simply right-clicking on the linsys1 object
in the Linear Analysis Workspace to copy the object. Then right-click within the MATLAB
Workspace to paste the object.
3.2 CONTROL
INVERTED PENDULUM: SIMULINK CONTROLLER DESIGN

Contents
 Problem setup and design requirements
 Implementing PID control for the nonlinear model
 Nonlinear closed-loop response
Problem setup and design requirements
In this problem, the cart with an inverted pendulum, shown below, is "bumped" with an
impulse force, .

For this example, let's assume that

(M) mass of the cart 0.5 kg
(m) mass of the pendulum 0.2 kg
(b) friction of the cart 0.1 N/m/sec
(l) length to pendulum center of mass 0.3 m
(I) inertia of the pendulum 0.006 kg*m^2
(F) force applied to the cart N
(x) cart position coordinate m
(theta) pendulum angle from vertical radians
In the design process we will develop a PID controller and apply it to a single-input, single-
output plant. More specifically, the controller will attempt to maintain the pendulum
vertically upward when the cart is subjected to a 1-Nsec impulse. The cart's position will be
ignored. Under these conditions, the design criteria are:
 Settling time of less than 5 seconds
 Pendulum should not move more than 0.05 radians away from the vertical
Implementing PID control for the nonlinear model
In the Inverted Pendulum: PID Controller Design page a PID controller was designed with
proportional, integral, and derivative gains equal to 100, 1, and 20, respectively. To
implement this closed-loop system, we will start with one of our plant models from
the Inverted Pendulum: Simulink Modeling page. Following the steps below, we will build a
closed-loop model with reference input of pendulum position and a disturbance force applied
to the cart.
 To begin, open either of the Simulink models generated previously. You can download
them by right-clicking on Pend_Model.slx or Pend_Model_Simscape.slx and then
selecting Save link as .... We will use the Simscape model to leverage the animation
capabilities it offers.
 Insert two Add blocks from Simulink/Math Operations library.
 Change the List of signs: of one of the Add blocks to "+-".
 Insert a Constant block from Simulink/Sources library. Change its value to 0. This is
the reference input that corresponds to the pendulum vertically upward. Note, the non-
Simscape model (and the rest of the pages in this example) define the pendulum angle
to equal pi when pointing straight up.
 Insert a PID Controller block from the Simulink/Continuous library.
 Edit the PID block by doubleclicking on it. Change the Proportional (P): gain to
"100", leave the Integral (I): gain as "1", and change the Derivative (D): gain to
"20".
 Now connect the blocks as they appear in the following figure:
You can download our version by right-clicking here and then selecting Save link as ....
Nonlinear closed-loop response
We can now simulate the closed-loop system. Be sure that the physical parameters are set as
follows.
M = 0.5;
m = 0.2;
b = 0.1;
I = 0.006;
g = 9.8;
l = 0.3;
Now, start the simulation (select Run from the Simulation menu or enter Ctrl-T). As the
simulation runs, an animation of the inverted pendulum will visualize the system's resulting
motion. Recall that the Show animation during simulation option must be checked under
the Simulation > Model Configuration Parameters menu. After the simulation has run you
should see the following response.
This response is almost identical to the closed-loop response obtained in the MATLAB
tutorials (for example, in the Inverted Pendulum: PID Controller Design page). Note that the
PID controller handles the nonlinear system very well because the deviation of the angle from
the operating point is very small (approximately .05 radians).
3.3 SIMSCAPE
Inverted Pendulum: Simscape Modeling

Contents
 Physical Setup
 Create world frame and basic configuration
 Assemble base plane and cart
 Pendulum subsystem and connecting the cart to the pendulum
 Selecting outputs for controller and angle conversion
 Create subsystem for pendulum and cart
 Closed-loop setup
 Controller implementation
Physical Setup
In this section we show how to build the inverted pendulum model using the physical
modeling blocks of Simscape Multibody. The blocks in the Simscape library represent actual
physical components; therefore, complex multibody dynamic models can be built without the
need to build mathematical equations from physical principles as was done by applying
Newton's laws to generate the model implemented in Inverted Pendulum: Simulink
Modeling.

The system parameters are defined as follows:

(M) mass of the cart 0.5 kg
(m) mass of the pendulum 0.2 kg
(b) coefficient of friction for cart 0.1 N/m/sec
(l) length to pendulum center of mass 0.3 m
(I) mass moment of inertia of the pendulum 0.006 kg.m^2
Create world frame and basic configuration
Open a new Simscape Multibody model by typing smnew in the MATLAB command
window. A new model opens, as shown below, with a few commonly used blocks already in
the model. The PS-Simulink and Simulink-PS blocks define the boundary between Simulink
input/output models where the blocks are evaluated sequentially and Simscape models where
the equations are evaluated simultaneously.

To configure basic settings in the model, carry out the following:

 Double-click on the Mechanism Configuration block and set Gravity to "[0 0 -9.81]",
this represents an acceleration due to gravity of acting along the global
-Z direction
 Open the Solver Configuration block and ensure that the Use local solver checkbox is
not selected.
 Type CTRL-E to open the Configration Parameters dialog
 On the Solver pane, ensure that Type is set to "Variable-Step" and that Solver is set to
"auto", also set the Stop time to "10"
 Assemble base plane and cart
We will model the cart as a point mass moving along an axis. We use the World Frame to
define the axis along which the cart will travel.
 Connect the B port of the Rigid Transform block to the W port of the World Frame
 Double-click on the Rigid Transform block
 In group Rotation, set Method to "Standard Axis", Axis to "+Y" and Angle to "90
deg"
 Rename the Rigid Transform block to "Transform Vehicle Axis"

Tips for adding blocks:

1. Use Quick Insert to add the blocks. Click in the diagram and type the name of the
block (use the letters in bold below). A list of blocks will appear and you can select
the block you want from the list.
2. After the block is entered, a prompt will appear for you to enter the parameter. Enter
the variable names as shown below.
3. To rotate a block or flip blocks, right-click on the block and select an option from
the Rotate & Flip menu.
4. To show the parameter below the block name, see Set Block Annotation Properties in
the Simulink documentation.
Add the following blocks:
* Prismatic Reference
* Scope
* Pulse Generator
For the Pulse Generator, double-click on the block and set Period to "10", Amplitude to
"1000", Pulse Width to "0.01", and Phase delay to "1".
The Prismatic Joint allows only one translational degree of freedom. The Prismatic Joint will
be actuated by a force input.
 Double-Click on the Prismatic Joint to open the dialog box
 In group Internal Mechanics, set Damping Coefficient to "0.1 N/(m/s)"
  Under group Actuation, set the Force to "Provided by Input"
 In group Sensing, select Position and Velocity
 Rename the Prismatic Joint to "Prismatic Cart"

The Prismatic Joint next needs to be connected to the rest of the model.
 Connect the B port of Prismatic Cart to the F port of block Transform Vehicle Axis
 Connect the F port of Prismatic Cart to the R port of the solid block
 Rename the Pulse Generator block to "Disturbance" and connect the output of the
"Disturbance" to the Simulink-PS Converter block that is already in the diagram
 Connect the output of the Simulink-PS Converter block to the force input of Prismatic
Cart
 Double-click on this signal and name it "Force"
 Double-click on the Simulink-PS Converter block and set Input signal unit to "N"
for newtons
 Make a copy of the PS-Simulink block
 Double-click on one PS-Simulink block and set Output signal unit to "m", and
connect that block to the p port of the Prismatic Cart block
 Double-click on the other PS-Simulink block, set the Output signal unit to "m/s" and
connect it to the v port of the Prismatic Cart block
  Connect both PS-Simulink blocks to the Scope
 Name the signal from the p port "x Cart"
 Name the signal from the v port "v Cart"
Since we are modeling the sliding cart as a point mass, only the mass will affect the
simulation results.
 Double-click on the solid block
 Set parameter Dimensions to "[0.2,0.04,0.6]" so that the longest dimension points
along the direction of travel.
 In group Inertia, set Type to "Point Mass" and Mass to "0.5 kg"
 In group Graphic, under Visual Properties set Color to "[0.8 0.45 0]"
 Rename the block to "Cart"

The model at this point should now appear as follows.

Running a simulation (type CTRL-T or press the green arrow run button), the resulting plot
shows both the distance traveled by the cart as well as its velocity.

4 Pendulum subsystem and connecting the cart to the pendulum

Add the following blocks to the model:
1. A Rigid Transform block
2. A Brick Solid block (Solid block prior to R2019b)
 3. A Revolute Joint
To define the axis of rotation for the pendulum:
 Connect the B port of the new Rigid Transform block to the F port of Prismatic Cart
 Connect the F port of the new Rigid Transform block to the B port of the new
Revolute Joint block
 Double-click on the new Rigid Transform block
 In group Rotation, set Method to "Standard Axis", Axis to "+X", and Angle to "90
deg"
 Rename the block "Transform Pendulum Pivot" the revolute

To define the degree of freedom of rotation for the pendulum:

 Rename Revolute Joint to "Revolute Pendulum"
 Connect F port of Revolute Pendulum to R port of Brick Solid block (Solid block
prior to R2019b)
To model the pendulum:
 Double-click on the Brick Solid block
 In group Geometry, set Dimensions to "[0.6 0.03 0.05] m"
 In group Inertia, set Type to "Point Mass" and Mass to "0.2 kg"
 In group Graphic, under Visual Properties set Color to "[0.25 0.4 0.7]"
 Press OK, and rename the block to "Pendulum"
To model the connection point to the cart:
 Double-click on the Pendulum block
 In group Frames, click on the + next to New Frame, this opens the frame definition
interface
 Click on the small face of the brick facing you (along positive x direction) to select it
 In section Frame Origin, select radio button Based on Geometric Feature
 Click the Save button at the bottom of the dialog box
 Change the name of "Frame1" to "B"
 Unselect the checkbox for Show Port R
 Click OK
 Connect the B port of Pendulum to the F port of "Transform Pendulum Pivot"
The resulting model should appear as follows:

Running a simulation (type CTRL-T or press the green arrow run button), the following plot
is generated, where one can see that the addition of the pendulum alters the cart behavior both
in its distance traveled as well as its velocity.
5 Selecting outputs for controller and angle conversion
We now need to measure the angle and angular velocity of the pendulum:
 Double-click on "Revolute Pendulum"
 In group Z Revolute Primitive (Rz) under Sensing, select checkboxes
for Position and Velocity
 Make two copies of the PS-Simulink converter block
 Double-click on one PS-Simulink block and set Output signal units to "rad"
 Connect that PS-Simulink block to the q port on Revolute Pendulum
 Double-click on the other PS-Simulink block and set Output signal units to "rad/s"
 Connect that PS-Simulink block to the w port on Revolute Pendulum
We need to limit the measured angle to stay between -pi and pi radians.
 Add a Subsystem block to your model
 Rename the Subsytem block "Wrap Angle"
 Double-click to enter the Wrap Angle subsystem
In this subsystem we will add pi radians to the measurement, find the remainder when the
signal is divided by 2*pi, and then subtract pi radians.
 Rename the input block "q"
 Delete the signal connection between the inport and the outport
 Add a Bias block to the model
 Set parameter Bias to "pi"
 Connect the q block to the Bias block
 Add a Math Function block to the model
 Double-click on the Math Function block and set Function to "rem"
 Connect the output of Bias to the first input of the Math Function block
 Add a Constant block to the model
 Set parameter Constant to "2*pi"
 Connect the Constant block to the second input of the Math Function
 Make a copy of the Bias block
 Set parameter Bias in the new block to "-pi"
 Connect Math Function output to the input of the new Bias block
 Connect the output of the new Bias block to the outport
 Rename the outport "qwrap"
 Go up one level in the diagram and rename the subsytem "Wrap Angle"
Here you can see the resulting subsystem for wrapping the angle

We will plot the new signals on a Scope

 Connect PS-Simulink output for the q measurement of Revolute Pendulum to the
input of Wrap Angle
 Add a new Scope
  Connect the qwrap output of Wrap Angle to the new Scope and change the name of
this signal to "q pendulum"
 Connect the PS-Simulink output for the w measurement of Revolute Pendulum to the
new Scope and change the name of this signal to "w pendulum"
The resulting model should appear as follows:

Running a simulation (type CTRL-T or press the green arrow run button), the followings
plot are generated. The motion of the cart is the same as before, but now we can see the
motion of the pendulum. The associated animation shown below is also generated.
 Create subsystem for pendulum and cart
We have now successfully created all the elements of the inverted pendulum system. We will
convert this into a subsystem.
 Click once in the diagram (but not on a block) and press CTRL-A to select all blocks
 Hold the Shift key and click on the Disturbance block and each Scope to unselect
those blocks
 Press CTRL-G to create a subsystem
 Rename the subsystem "Inverted Pendulum on Cart"
Your model should now appear as shown:
 Closed-loop setup
We will now add blocks for open and closed-loop testing.
Add the following blocks:
 Subtract
 PID Controller
 Constant
 Manual Switch
 Sum
The Disturbance also will be added to the control signal.
 Delete the signal connecting Disturbance to the cart subsystem
 Connect the output of the sum block to the Force input of the cart subsystem
 Connect Disturbance to the bottom + port of the Sum block
We wish to manually select open or closed-loop behavior.
 Connect the Manual Switch output to the + input of the Sum block
 Connect the Constant block to the lower input of the Manual Switch, then set the
parameter Constant to "0" and rename block "No Force"
 Connect the output of the PID Controller to the upper input of Manual Switch
 Connect the Subtract block output to the input of PID Controller
We will close the loop to control the pendulum angle.
 Connect q pendulum output of the cart subsystem to the - port of the Subtract block
 Make a copy of the Constant block and connect it to the + port of the Subtract block
and rename the block "Desired Pendulum Angle"
Your model should now appear as follows. The output of the simulation is unchanged from
prior results when in open-loop mode.

Controller implementation
We will now implement the PID control gains developed in the Inverted Pendulum: PID
Controller Design page.
 Double-click on the PID Controller block
 Set parameter Proportional (P) to "100"
 Set parameter Integral (I) to "1"
 Set parameter Derivative (D) to "20"
Double-click on Manual Switch until the input from the PID Controller is selected.

Running a simulation produces the following two plots show the controlled response of the
system. After the initial impact the controller was able to quickly bring down the pendulum
angle to zero and the pendulum velocity is also zero. The cart moves slowly and with a
constant velocity in the negative X direction to keep the pendulum balanced.