Digital Control Lecture Notes
Digital Control Lecture Notes
Digital Control Lecture Notes
CONTROL
CONTENTS:
2 LYAPUNOV STABILITY 34
2.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 34
2.2 Equilibrium points . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35
i
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau ii
3 DIGITAL CONTROL 47
3.1 Introduction to digital control systems . . . . . . . . . . . . . . . . . . . . . . . . 47
3.1.1 System configuration and terminology . . . . . . . . . . . . . . . . . . . . 47
3.1.2 Advantages and disadvantages of digital control . . . . . . . . . . . . . . . 48
3.2 The z-transform . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 49
3.2.1 Definition of the z-transform . . . . . . . . . . . . . . . . . . . . . . . . . 49
3.2.2 z-transforms of common functions . . . . . . . . . . . . . . . . . . . . . . 52
3.2.3 Theorems of z-transform . . . . . . . . . . . . . . . . . . . . . . . . . . . . 53
3.2.4 Inverse z-transform . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 57
3.3 Difference equations and pulse transfer functions . . . . . . . . . . . . . . . . . . 59
3.3.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 59
3.3.2 Solving difference equations by recursion . . . . . . . . . . . . . . . . . . 60
3.3.3 Pulse transfer functions and weighting sequence . . . . . . . . . . . . . . 61
3.4 Mapping between the s-plane and the z-plane . . . . . . . . . . . . . . . . . . . . 62
3.4.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 62
3.4.2 Mapping of the jω axis to the z-plane . . . . . . . . . . . . . . . . . . . . 63
3.4.3 Lines of constant attenuation . . . . . . . . . . . . . . . . . . . . . . . . . 63
3.4.4 Lines of constant frequency . . . . . . . . . . . . . . . . . . . . . . . . . . 64
3.4.5 Important note . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66
3.5 Data holds . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66
3.5.1 Introduction . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 66
3.5.2 The Zero-Order-Hold (ZOH) . . . . . . . . . . . . . . . . . . . . . . . . . . 67
1 − e−T s
3.5.3 z-transform of functions involving the term . . . . . . . . . . 69
s
3.6 Block diagram analysis . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
3.6.1 Open-loop systems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70
3.6.2 Closed-loop systems . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 74
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau iii
In the previous course ELC 412, we learnt that if we can measure the n internal state
variables x1 , x2 · · · xn of a system, then it is possible to feed back each state variable
xi through a gain ki to the input of the system. By choosing the feedback gains K =
[k1 k2 · · · kn ], we can arbitrarily locate the poles of the closed-loop system. However,
the system must be completely state controllable, and all the states must be accessible for
feedback.
It is assumed without loss of generality that the system is a Single Input Single Output
(SISO) system with a state space representation:
The state matrix of the closed-loop system is therefore (A − BK). The block diagram of
the closed-loop system is shown in Figure 1.1.
When designing by pole placement, it is important to keep in mind that the control effort
required is proportional to how far the system poles are moved by the state feedback
- the further the poles move, the greater the control effort required (the K matrix will
have high gains).
1
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 2
rag replacements
ẋ(t) R x(t) y(t)
B C
+
+
A
SYSTEM
u(t) = −Kx(t)
−K STATE
FEEDBACK
There are three methods of determining the feedback gain matrix K and they are out-
lined below without proofs – a reader interested in the proofs can refer to the handout
for ELC 412.
Suppose the desired closed-loop poles are µ1 , µ2 , . . . µn . Then the desired char-
acteristic equation is
Then
h i −1
K = α0 − a 0 α1 − a 1 α2 − a 2 ··· αn−1 − an−1 MW ,
(1.6)
where M is the controllability matrix:
M = B AB A2 B An−1 B
··· (1.7)
and
a1 a2 a3 · · · an−1 1
a2 a3 a4 ··· 1 0
W= .. .. .. .. .. ..
. (1.8)
. . . . . .
an−1 1 0 ··· 0 0
1 0 0 ··· 0 0
where
φ (A) = An + αn−1 An−1 + αn−2 An−2 + · · · + α1 A + α0 I
(compare this with the desired characteristic equation (1.5)). This method is
known as Ackermann’s formula.
Method 3: Direct substitution method From equation (1.3), the characteristic equa-
tion of the closed loop system is
sI − (A − BK) = sI − A + BK = 0. (1.10)
For the direct substitution method, the above equation is equated to the desired
characteristic equation (1.5), i.e.
This is known as the Direct Substitution Method. The elements of the matrix K
are obtained by equating coefficients of equal powers of s in (1.11) above. Note
that this method becomes tedious if the number of states n is large - it is mostly
used for cases where n ≤ 3.
Example
Given is a continuous-time system:
" # " #" # " #
ẋ1 (t) 0 1 x1 (t) 0
= 2
+ u(t).
ẋ2 (t) −ω0 0 x2 (t) 1
This system has open-loop poles at s = ±jω0 . We would like to use state feedback law
u = −Kx to locate both of the closed-loop poles at s = −2ω0 (i.e. the closed-loop
system should have double poles at s = −2ω0 ). Determine K.
Solution:
The system is first tested for controllability. The controllability matrix is
" #
0 1
M = [B AB] = .
1 0
|M| = −1 6= 0,
hence the system is completely state controllable and design by pole placement is pos-
sible.
|sI − A| = s2 + ω02 = s2 + a1 s + a0 ,
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 4
α1 − a1 ] (MW)−1 = 3ω02
K = [α0 − a0 4ω0 .
s −1
= s2 + sk2 + ω02 + k1 = s2 + 4ω0 s + 4ω02 .
sI − A + BK =
ω02 + k1 s + k2
Equating coefficients of equal powers of s gives k1 = 3ω02 and k2 = 4ω0 , hence
K = 3ω02
4ω0 .
1.2.1 Introduction
For design by state variable feedback (design by pole placement), it is assumed that all
the states are available for feedback. However, there are situations where not all states
are available for feedback for the following reasons:
When some or all the states are not available, we need to estimate the unavailable state
variables. This is done using a state observer (also known as a state estimator). A state
observer is an electronic circuit or a computer program which estimates the states of
the system from the measured input and output signals.
An observer which estimates all the state variables of a system regardless of whether
some of the state variables are available for direct measurement is known as a full-order
state observer.
An observer which estimates fewer than n state variables, where n is the order of the sys-
tem is called a reduced order state observer. If the reduced order observer only estimates
those states which cannot be measured, it is called a minimum-order state observer.
A state observer can only be designed if the system is observable. The design also
requires an accurate model of the system.
A block diagram of a system connected to a full-order state observer is shown in Fig-
ure 1.2.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 6
SYSTEM
ẋ(t) x(t) y(t)
C
R
B
+
+
replacements
A
u(t)
x̂(t)
−+
˙
x̂(t) x̂(t) ŷ(t)
C
R
B
++ +
+
−K
Ke y(t) − ŷ(t) y(t) − ŷ(t)
Ke
FULL ORDER OBSERVER
Note that the full-order observer has two inputs, u(t) and y(t). The feedback term
Ke (y(t) − ŷ(t)) is needed to ensure that
lim x̂(t) −→ x(t), (1.13a)
t→∞
lim ŷ(t) −→ y(t), (1.13b)
t→∞
i.e.
lim x̂(t) − x(t) = 0, (1.14a)
t→∞
lim ŷ(t) − y(t) = 0. (1.14b)
t→∞
From Figure 1.2, we can write the state equation of the full-order observer (the lower
part of the block diagram) as:
˙
x̂(t) = Ax̂(t) + Bu(t) + Ke y(t) − ŷ(t) . (1.15)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 7
˙
x̂(t) = Ax̂(t) + Bu(t) + Ke Cx(t) − Cx̂(t) . (1.16)
For the system (the upper part of the block diagram), the state equation is given by
Subtracting equation (1.16) from equation (1.17) gives the observer error equation, i.e.
the equation which gives the error between the actual state x(t) and the estimated state
x̂(t):
˙
ẋ(t) − x̂(t) = Ax(t) + Bu(t) − Ax̂(t) − Bu(t) − Ke Cx(t) − Cx̂(t)
= A x(t) − x̂(t) − Ke C x(t) − x̂(t) = A − Ke C x(t) − x̂(t) . (1.18)
Defining the error vector e(t) = (x(t) − x̂(t)), then the above equation can be written as
ė(t) = A − Ke C e(t). (1.19)
From (1.19) above, it can be seen that the behaviour of the error vector depends on
the eigenvalues of the matrix (A − Ke C). If the matrix (A − Ke C) is a stable matrix,
then the error vector e(t) will converge to zero for any initial error vector e(0), i.e. x̂(t)
will converge to x(t) regardless of the x̂(0) and x(0). The eigenvalues of the matrix
(A − Ke C) also determine how fast the estimated state x̂(t) converge to the state x(t).
The characteristic equation for the error vector system is
Now recall dual systems discussed in ELC 412. For a system given by the state space
representation
where z(t) is the state vector of the dual system, v(t) the input and w(t) the output.
Consider system (1.21). For the pole placement problem, the control law is given by
equation (1.2), reproduced here for convenience:
Now consider the dual system (1.22). For the pole-placement problem, the correspond-
ing control law is:
v(t) = −Kz(t). (1.26)
The closed-loop system is therefore given by
ż(t) = AT z(t) + CT − Kz(t) = A − CT K z(t). (1.27)
The state matrix of the closed-loop system is (AT −CT K) and this matrix determines the
eigenvalues of the closed-loop dual system. The characteristic equation of the closed-
loop dual system is
sI − AT + CT K = 0. (1.28)
From matrix theory, the eigenvalues of a matrix F are the same as the eigenvalues of the
transpose of that matrix FT . It follows that the eigenvalues of the matrix (AT − CT K)
are the same as those of (AT − CT K)T .
But T
AT − C T K = A − KT C.
From the above equation, the characteristic equation of the closed-loop dual system
given in (1.28) can be written as
sI − A + KT C = 0. (1.29)
Comparing equations (1.20) and (1.29), it can be seen that the matrix Ke can
be determined by pole placement approach for the dual of the original system
then using Ke = KT .
For pole placement of the dual system, the dual system (1.22) must be completely state
controllable. The controllability matrix of the dual system (1.22) is
h i
T T T T 2 T n−1
T T
N= C A C A C ··· A C . (1.30)
We recognize (1.30) as the observability matrix of the original system (1.21), which
leads us to the necessary condition for observer design:
A system (1.21) must be observable for it to be possible to design an observer.
Conversely, if the system (1.21) is observable, then it is possible to choose a
matrix Ke such that the matrix (A − Ke C) has the desired eigenvalues.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 9
We now look at the methods of determining the matrix Ke . There are three methods
which are derived from the methods given in Section 1.1.
Determination of Ke : Method 1
Since the characteristic equations of the original system and its dual are identical, then
a1 a2 a3 · · · an−1 1
a2 a3 a4 · · · 1 0
. .. .. . . .. ..
W= .. . . . . . (1.35)
an−1 1 0 · · · 0 0
1 0 0 ··· 0 0
But Ke = KT hence
h i T
−1
Ke = α0 − a 0 α1 − a 1 α2 − a 2 ··· αn−1 − an−1 (NW) , (1.38)
or
α0 − a 0
α1 − a 1
−1
Ke = WNT α2 − a 2 . (1.39)
..
.
αn−1 − an−1
Determination of Ke : Method 2
But Ke = KT hence
h h 2 n−1 T i−1 iT
Ke = [0 0 · · · 0 1] CT AT CT A T CT · · · A T C φ AT
T h T i−1 T
T 2 T n−1
T T T
[0 0 · · · 0 1]T .
T T
= φ A C A C A C ··· A C
(1.41)
T
From matrix theory, φ AT = φ (A), hence
−1
C 0
CA 0
CA2
Ke = φ (A) 0 , (1.42)
.. ..
. .
CAn−1 1
where
φ (A) = An + αn−1 An−1 + αn−2 An−2 + · · · + α1 A + α0 I
(compare this with the desired characteristic equation (1.36)). This method is known
as Ackermann’s formula for observer design.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 11
Determination of Ke : Method 3
The location of the observer poles is a trade-off between the time it takes for the esti-
mated state x̂(t) to converge with the actual state x(t), the control effort required and
the effects of noise and disturbances.
As a rule of thumb, the observer poles (eigenvalues of the matrix (A−Ke C)) are chosen
to be a factor of 2 to 6 times faster than the controller poles (eigenvalues of the matrix
(A − BK) from pole placement).
If the observer poles are chosen faster than 6 times the controller poles, the estimated
state x̂(t) converges with the actual state x(t) in a short time but a large control effort
is required and the effects of noise and disturbances are highly amplified.
If the observer poles are chosen slower than 2 times the controller poles, the estimated
state x̂(t) takes a long time to converge with the actual state x(t), but the control effort
required is low and the effect of noise and disturbances is minimal.
1.2.5 Example
x3 (t)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 12
suppose that the states of the system are not available for feedback. Design a full order
state observer with the eigenvalues of the matrix (A − Ke C) are located at s = −4,
s = −3 ± j.
Solution:
Observability matrix
C 0 0 1
N = CA = 0 2 0 .
CA2 6 −2 2
N = 12 6= 0
hence the system is observable and observer design is possible.
Method 1:
Since the given system is a third order system,
−1 α0 − a0
Ke = WNT α1 − a 1 .
α2 − a 2
The desired characteristic equation is
(s + 4)(s + 3 + j)(s + 3 − j) = s3 + 10s2 + 34s + 40 = s3 + α2 s2 + α1 s + α0 = 0,
hence α2 = 10, α1 = 34 and α0 = 40.
Also
sI − A = s3 − 9s + 2 = s3 + a2 s2 + a1 s + a0 ,
hence a2 = 0, a1 = −9 and a0 = 2.
h i 0 0 6
2
N = CT A T CT AT CT = 0 2 −2 .
1 0 2
For a third order system,
a1 a2 1 −9 0 1
W = a2 1 0 = 0 1 0 .
1 0 0 1 0 0
From the above working,
−1
−1 α 0 − a 0 −9 0 1 0 0 1 38
T
Ke = WN α1 − a 1 = 0 1 0 0 2 0 43
α2 − a 2 1 0 0 6 −2 2 10
1 1 7
6 6 6
38 25 61
= 0 12 0 43 = 21 12 .
0 0 1 10 10
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 13
Method 2:
Ackermann’s formula for observer design for a third order system is
−1
C 0
Ke = φ (A) CA 0
CA2 1
where
φ (A) = A3 + α2 A2 + α1 A + α0 I.
From the working above, α2 = 10, α1 = 34 and α0 = 40, hence
151 86 20
φ (A) = A3 + 10A2 + 34A + 40I = 129 85 33 .
60 66 58
Also,
C 0 0 1
CA = 0 2 0 .
CA2 6 −2 2
It follows that
−1
151 86 20 0 0 1 0
Ke = 129 85 33 0 2 0 0
60 66 58 6 −2 2 1
151 86 20 − 31 16 1
6
0 25 1
6
= 129 85 33 0 12 0 0 = 21 12 .
60 66 58 1 0 0 1 10
Method 3:
For the direct substitution method for observer design,
Let
k e1
Ke = k e 2 .
k e3
Then
s 0 0 1 2 0 0 0 k e1
sI − A + Ke C = 0 s 0 − 3 −1 1 + 0 0 k e2
0 0 s 0 2 0 0 0 k e3
s − 1 −2 k e1
= −3 s + 1 −1 + ke2 .
0 −2 s + k e3
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 14
10
1.3.1 Introduction
A combined state feedback controller and observer is shown in Figure 1.3. The dif-
ference between Figures 1.2 and 1.3 is that in Figure 1.2, the estimated state is not
connected to anything but in Figure 1.3, the estimated state is used for state feedback
and there is no external input.
From Figure 1.3, we can write for the system (the upper part of the block diagram)
Defining the error vector e(t) = (x(t) − x̂(t)), then the above equation can be written as
ẋ(t) = A − BK x(t) + BKe(t). (1.47)
For the differential equation for the observer (the lower part of Figure 1.3), recall the
observer error equation (1.19), reproduced here for convenience:
ė(t) = A − Ke C e(t). (1.48)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 15
SYSTEM
ẋ(t) x(t) y(t)
C
R
B
+
+
replacements A
u(t) = −Kx̂(t)
−K
˙
x̂(t) x̂(t) ŷ(t) +
−
C
R
B
+ +
+ +
y(t) − ŷ(t)
Ke
Ke y(t) − ŷ(t)
FULL ORDER OBSERVER
sI − A + BK −BK
= 0. (1.50)
0 sI − A + Ke C
sI − A + BK sI − A + Ke C = 0. (1.52)
The equation |sI − A + BK| is the characteristic equation of the closed-loop system
with state feedback (compare with equation (1.11)). On the other hand, |sI − A + Ke C|
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 16
is the characteristic equation of the observer (compare with equation (1.20)). This leads
us to the following conclusion:
Pole placement and observer design are independent of each
other, i.e. K and Ke are determined independently.
Note also that the poles of the closed-loop system shown in Figure 1.3 are the poles of
the controller (eigenvalues of (A − BK)) and the poles of the observer (eigenvalues of
(A − Ke C)). These poles can also be obtained by solving equation (1.52).
From Figure 1.3, the differential equation for the observer (the lower part of the block
diagram) is
˙
x̂(t) = Ax̂(t) + Bu(t) + Ke y(t) − ŷ(t) . (1.53)
But
ŷ(t) = Cx̂(t) (1.54)
and
u(t) = −Kx̂(t). (1.55)
Substituting (1.54) and (1.55) into (1.53) gives
˙
x̂(t) = A − BK − Ke C x̂(t) + Ke y(t). (1.56)
To obtain the transfer function of the controller-observer, we take the Laplace trans-
forms of (1.56) and (1.55) and eliminate the state variable.
The Laplace transform of (1.56) is
sX̂(s) = A − BK − Ke C X̂(s) + Ke Y (s), (1.57)
We make X̂(s) the subject of the formula by pre-multiplying both sides of (1.58) by
(sI − A + BK + Ke C)−1 , resulting in
−1
X̂(s) = sI − A + BK + Ke C Ke Y (s). (1.59)
U (s) −1
= K sI − A + BK + Ke C Ke . (1.62)
−Y (s)
Exercise:
Given a continuous-time system
" # " #" # " #
ẋ1 (t) 0 1 x1 (t) 0
= + u(t)
ẋ2 (t) 0 0 x2 (t) 1
" #
h i x (t)
1
y(t) = 1 0 + [0]u(t),
x2 (t)
(i) Determine K such that the poles of the state feedback controller (roots of the
equation |sI − A + BK|) are located at s = −2 ± j2.
(ii) Determine Ke such that the observer poles (roots of the equation
|sI − A + Ke C|) are located at s = −4 ± j4.
Solutions:
(i) K = [8 4]
" #
8
(ii) Ke =
32
−1 64(3s + 4)
(iii) K sI − A + BK + Ke C Ke = .
s2+ 12s + 72
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 18
1.4.1 Introduction
So far, we have considered systems which do not have external inputs. Such systems
are known as regulator systems. A regulator control system has no external input and
the purpose of the system is to return all state variables to zero when the states have
been perturbed.
However, many systems require that the output track an input signal. Such a control
signal, in which the output must follow command signals, is known as a servo system.
In a servo system it is generally required that the system have one or more integrators
within the closed-loop. If the plant to be controlled does not have an integrator, it is
usually necessary to add one or more integrators within the loop to eliminate/minimize
steady state errors.
A
SYSTEM
STATE
K FEEDBACK
In this case, we assume that the system has at least one integrator. Since we are feeding
back the states and not the output, the input gain required for the output to track
the input is not necessarily equal to 1 so the input gain is set to kr . To simplify the
derivations, we assume that the external input is a unit step.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 19
The system shown in Figure 1.5 is an n-th order system, but equation (1.64) has (n + 1)
unknowns k1 , k2 , . . . kn and kr , so we cannot uniquely determine all the unknowns in
equation (1.64).
To be able to solve for all the gains, we reconfigure the closed-loop system so that it
behaves like a regulator. From the design point of view, we desire
hence
lim y(t) − r(t) = 0. (1.66)
t→∞
Equations (1.65) and (1.66) are used to help us convert the system shown in Figure 1.5
to a system which behaves like a regulator.
To convert the system into a regulator system, we let the output y(t) be equal to one of
the states, usually state x1 (t), and we also let kr be equal to one of the feedback gains
(the gain corresponding to the state equal to the output).
Let
hence the system behaves like a regulator system, and the methods used to determine
the feedback gains for regulator systems can now be applied to this system.
Equation (1.68) can be rewritten as
PLANT
replacements x1 (t)
r(t) u(t) x2 (t)
y = x1 (t)
k1 ẋ(t) = Ax(t) + Bu(t) x3 (t) C
+ +
xn (t)
k2
k3
kn
Figure 1.6: State feedback where there is an external input – alternative configuration
which is the state equation of the closed-loop system. The characteristic equation of the
closed-loop system is
sI − A + BK = 0 (1.74)
(compare this with the characteristic equation of the regulator system, equation (1.10)).
For the system configuration shown in Figure 1.6, determine the feedback gain matrix
K so that the control law u(t) = −Kx(t) + k1 r(t), where r(t) is a unit step input, results
in a closed-loop system with poles located at s = −4 ± j4.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 21
Solution
Y (s) 1
= 2
U (s) s
Hence
s2 Y (s) = U (s).
ÿ(t) = u(t).
Since this is a system with an external input, the states are selected as x1 (t) = y(t) and
x2 (t) = ẏ(t).
The state equations can then be written as:
|M| = −1 6= 0,
hence the system is completely state controllable.
Since no pole placement method is specified, we can use any of the three methods in
Section 1.1 – in this case, Ackermann’s formula is chosen.
The desired characteristic equation is
(s + 4 + j4)(s + 4 − j4) = s2 + 8s + 32 = 0 = s2 + α1 s + α0 = 0,
hence
" #" # " # " # " #
0 1 0 1 0 1 1 0 32 8
φ(A) = A2 + 8A + 32I = +8 + 32 = .
0 0 0 0 0 0 0 1 0 32
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 22
A block diagram of the closed-loop system is shown in Figure 1.7 (compare this with
replacementsthe block diagram shown in Figure 1.6).
PLANT
r(t) u(t) 1 x2 (t) 1 y(t) = x1 (t)
32
+ + s s
k1
8
k2
and the transfer function of the closed-loop system can be obtained using the equation
Y (s)
= C (sI − A + BK)−1 Bk1 .
R(s)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 23
Exercise:
A plant is given by the transfer function
Y (s) 1
= .
U (s) s(s + 1)(s + 2)
For the system configuration shown in Figure 1.6, use the direct substitution method
to determine the feedback gain matrix K so that the control law
where r(t) is a unit step input, results in a closed-loop system with poles located at
√
s = −2 ± j2 3, s = −10, and draw a block diagram of the resulting system.
Solution:
h i
K= 160 54 11 .
State space design does not produce integral control unless special steps are taken in
the design process. If the system does not have an integrator (i.e. is type 0), the system
will have finite steady-state errors for step inputs and infinite steady-state errors for
ramp and parabolic inputs.
replacements
To simplify the mathematical analysis, we restrict our system to step inputs. To eliminate
steady-state errors for step inputs, we need to add an integrator to a type 0 system to
convert it to a type 1 system. The integrator is added to the closed-loop system as shown
in Figure 1.8.
SYSTEM
r(t) ζ̇(t) R ζ(t) u(t) ẋ(t) x(t) y(t) = x1 (t)
kI
R
+ + B +
C
INTEGRATOR +
STATE
K FEEDBACK
The closed-loop system shown in Figure 1.8 is of the order (n + 1). From the figure, we
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 24
At t = ∞,
" # " #" # " # " #
ẋ(∞) A 0 x(∞) B 0
= + u(∞) + r(∞), (1.78)
ζ̇(∞) −C 0 ζ(∞) 0 1
and
u(∞) = −Kx(∞) + kI ζ(∞). (1.79)
Subtracting (1.78) from (1.76), keeping in mind that r(t) = r(∞) = 1, gives
" # " #
ẋ(t) − ẋ(∞) A 0 x(t) − x(∞) B
= + u(t) − u(∞) . (1.80)
ζ̇(t) − ζ̇(∞) −C 0 ζ(t) − ζ(∞) 0
Let
and " #
h i xe (t)
ue (t) = −Kxe (t) + kI ζe (t) = − K −kI (1.84)
ζe (t)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 25
respectively.
and
ue (t) = −K̂e(t) (1.86)
respectively.
Equations (1.85) and (1.86) are similar to the regulator equations (1.1a) and (1.2),
which leads us to the following conclusion:
Example
A plant is given by the transfer function
Y (s) 1
= .
U (s) s+3
Use Ackermann’s formula to design a type 1 servo system with closed-loop poles at
√
s = −2.5 ± j2.5 3.
Solution
The given system is a first order system, and after adding an integrator, the final system
should be a type 1 second order system which gives zero steady-state errors to step
inputs.
Since the system has an external input, the output should be chosen to be equal to one
of the states.
From the given transfer function above, we can write
This is a first order system so there is only one state, which we can name x(t).
Letting
x(t) = y(t),
then
ẋ(t) = ẏ(t) = −3y(t) + u(t) = −3x(t) + u(t).
The state space representation of the system is
M̂ = −1 6= 0,
hence the system is completely state controllable.
Also " #
1 −3
M̂−1 = .
0 −1
The desired characteristic equation is
√ √
s + 2.5 + j2.5 3 s + 2.5 − j2.5 3 = s2 + 5s + 25.
Hence " #
19 0
φ(Â) = Â2 + 5Â + 25I = .
−2 25
From the working above,
" #" #
h i h i 1 −3 19 0 h i
K̂ = k1 −kI = 0 1 = 2 −25 ,
0 −1 −2 25
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 27
PLANT
replacements
−3
+
r(t) 1 u(t) + 1 y(t) = x(t)
25
+ s + s
kI
2
k1
Y (s) 1
= .
U (s) (s + 1)(s + 2)
Use the direct substitution method to design a type 1 servo system with closed-loop
√
poles at s = −10, s = −2 ± j2 3, and draw a block diagram of the resulting system.
Solution:
h i h i
K̂ = k1 k2 −kI = 54 11 −160 .
1.5.1 Introduction
A system with state feedback control is superior to a PID controller since it allows the
placing of all roots of the closed-loop system at the desired locations.
The roots of the closed-loop system are selected to satisfy the performance specifica-
tions. It thus follows that the first task in state-feedback design is to translate all per-
formance specifications into desired root locations. The control law is then obtained by
feeding back all the state variables with constant gains. The values of the gains depend
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 28
When selecting poles, it is useful to remember that the control effort required is related
to how far the open-loop poles are moved by the state feedback – the further the poles
are moved, the higher the control effort required.
Two commonly-used methods of pole selection criteria are
prototype design
In this case, it is assumed that the system has a dominant second order transient which
has complex poles at a radius ωn and a damping factor ζ. The rest of the poles are
selected to have real parts corresponding to sufficiently damped modes of the system so
that the system approximates to a second order system (see Figure 1.10). This method
is commonly used when the performance specifications are given in terms of transient
response parameters e.g. rise time, overshoot, settling time.
This method is used when the system cannot be usefully approximated as a second
order system.
There are different prototypes and they have different characteristics. One category
of prototypes is based on performance indices. A performance index is a number that
indicates the “goodness” of a system performance.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 29
j ω axis
s−plane
ωn sin θ = ζ
θ
real axis
One commonly used index is the Integral of Time Multiplied by Absolute Error (ITAE)
criterion. This method minimizes the performance index
Z ∞
J= t|e(t)|dt. (1.87)
0
The ITAE criterion has been applied to a closed-loop transfer function of the form
Y (s) a0
= n (1.88)
U (s) s + an−1 sn−1 + · · · + a1 s + a0
and the coefficients which minimize the performance index J in equation (1.87) have
been determined. The ITAE characteristic equations for systems of order 1 to order 6
are given below (in each case, a0 = (ω0 )n ):
Design of a system to meet the ITAE criterion results in a system with a small rise time and
zero steady-state error to a step input.
The designer tries different values of ω0 and selects the value which gives the best
response. In general, increasing the value of ω0 makes the system faster (smaller rise
time and smaller settling time).
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 30
The Bessel prototype is based on Bessel functions. For a closed-loop transfer function
of the form of equation (1.88), the Bessel characteristic equations for systems of order
1 to order 6 are given below (in each case, a0 = (ω0 )n ):
Again, the designer tries different values of ω0 and selects the value which gives the best
response. In general, increasing the value of ω0 makes the system faster (smaller rise
time and smaller settling time).
Systems designed using the Bessel prototypes has no overshoots and zero steady-state error
for step inputs. In addition, for a given value of ω0 , the step response of a Bessel prototype
is slower than the step response of an ITAE prototype of the same order.
Figure 1.11 shows the step responses of the ITAE prototypes given by equations (1.89a)
to (1.89f) while Figure 1.12 shows the step responses of the Bessel prototypes given by
equations (1.90a) to (1.90f). In all cases, a value of ω0 = 1 is used. As an example, to
simulate the 6th order ITAE step response, the (closed-loop) transfer function used is
Y (s) ω06
= 6
U (s) s + 3.25ω0 s5 + 6.6ω02 s4 + 8.6ω03 s3 + 7.45ω04 s2 + 3.95ω05 s + ω06
As can be seen from Figure 1.12, the Bessel prototypes result in a closed-loop system
with no overshoots.
1.5.6 Example:
1.2 1.2
0.8 0.8
0.6 0.6
0.4 0.4
0.2 0.2
0 0
0 2 4 6 8 10 12 0 2 4 6 8 10 12
time −−> (sec) time −−> (sec)
Figure 1.11: Step response – ITAE prototypes Figure 1.12: Step response – Bessel prototypes
Solution:
Y (s) 10
= 2 (1.91)
U (s) s + 2s
Hence
s2 Y (s) + 2sY (s) = 10U (s) (1.92)
Since this is a system with an external input, the states are selected as x1 (t) = y(t) and
x2 (t) = ẏ(t).
The state equations can then be written as:
s2 + 2.8s + 4 (1.98)
which is the state equation of the closed-loop system. The output equation of the
closed-loop system is the same as the output equation of the open-loop system (equation
(1.95b).
Note that the transfer function of the closed-loop system is given by
Y (s) 4
= C (sI − A + BK)−1 Bk1 = 2 (1.104)
R(s) s + 2.8s + 4
A Simulink model to simulate resulting system is shown in Figure 1.13 (the transfer
4
function block s2 +2.8s+4 is included for comparison of the system response with the
expected response).
Exercise
Repeat the above example, this time designing for a closed-loop system with poles
at the second order Bessel locations for ω0 = 2. Show that the closed-loop transfer
function of the resulting system is given by
Y (s) 4
= 2 .
R(s) s + 3.46s + 4
Chapter 2
LYAPUNOV STABILITY
2.1 Introduction
Stability plays a major role in control system design. There are many methods for
determining the stability of a given system. In this course, it is assumed that the reader
is familiar with the following methods of determining the stability of a system:
The above methods are known as classical methods, and they are only applicable if the
system being considered is Linear-Time-Invariant (LTI).
However, if the system is nonlinear, or linear but time varying, the methods above
cannot be used to determine the stability of such systems.
The Lyapunov stability methods are more general methods of determining the stabil-
ity of systems, and they are applicable to both linear and non-linear, as well as time-
invariant and time-varying systems. These methods are named after a Russian scientist,
A. M. Lyapunov, who studied dynamic systems developed the methods in the 19th cen-
tury.
Lyapunov developed 2 methods of determining the stability of systems: the first method
(also known as the indirect method) and the second method (a.k.a. the direct method).
First method (indirect method): This method involves solving of the differential equa-
tions (linear or nonlinear) describing the system to determine the stability of the
system. This method is not commonly used because it is usually difficult and
time-consuming to solve nonlinear differential equations.
34
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 35
Second method (the direct method): Does not require the solution of the differential
equations describing the system hence by using this method, the stability of the
system can be determined without having to solve the state equations.
The second method of Lyapunov is the most general method for determination of
the stability of nonlinear and/or time-varying systems as well as hybrid/switched
systems. This method also applies to linear-time-invariant systems.
ẋ = f (x, t) (2.1)
i.e. since the states have settled down to constant values, ẋe = 0.
Example
Find the equilibrium state of the autonomous LTI system described by the state equation
" # " #" #
ẋ1 (t) −1 2 x1 (t)
= (2.3)
ẋ2 (t) 1 −1 x2 (t)
Solution
At the equilibrium state xe , equation (2.2) applies hence
which gives a solution x1 = x2 = 0, hence the equilibrium state for this system is at the
origin, i.e. (x1e , x2e ) = (0, 0).
In general, a linear time-invariant autonomous system
ẋ = Ax (2.5)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 36
Solution
Let the equilibrium point be (x1 , x2 ) = (x1e , x2e ). At the equilibrium point, ẋ1e = ẋ2e = 0
hence we can write
ẋ1 = x2 (2.10a)
ẋ2 = −x1 − x21 − x2 (2.10b)
Solution
At equilibrium,
0 = x2e (2.11a)
0 = −x1e − x21e − x2e (2.11b)
hence x2e = 0 and x1e =0 or -1, hence the equilibrium states for this system are (x1e , x2e ) =
(0, 0) and (x1e , x2e ) = (−1, 0).
In general, a non-linear system may have more than one equilibrium points.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 37
Fixed
surface
R
θ
Example
Consider the pendulum shown in Figure 2.1. It is assumed that the weight at the bottom
of the pendulum is joined to the fixed surface using a rigid body (not a string!) of
negligible mass.
The dynamics of the pendulum are given by the nonlinear differential equation
where R is the pendulum’s length, M its mass, b the frictional coefficient at the hinge,
and g the gravitational acceleration. The angle θ is in radians.
Defining two states x1 and x2 such that x1 = θ and x2 = θ̇, the state equations of the
system are given by
ẋ1 = x2 (2.13a)
g b
ẋ2 = − sin x1 − x2 (2.13b)
R M R2
At equilibrium, x2e = 0, and sin x1e = 0. This gives two equilibrium states (x1e , x2e ) =
(0, 0) and (x1e , x2e ) = (π, 0). Physically, these points correspond to the pendulum resting
in the bottom position and the upright position respectively.
x2 2−dimensional
state space
x0
0 x1
PSfrag replacements
system
trajectory
Suppose that the state vector x lies within a hyperspherical region of radius R, i.e.
kxk < R. (In a 2-dimensional system, the hyperspherical region would be a circle of
radius R, for a 3-dimensional system, the hyperspherical region would be a sphere of
radius R, etc.)
Define two regions with radius δ and ε within R such that δ < ε < R (see Figure 2.3 for
a 2-dimensional state-space).
x2
2−dimensional
state space
R
ε
0 δ x1
PSfrag replacements
radius ε, there exists a radius δ < ε so that every trajectory starting within the region of
radius δ will always remain within radius ε with continually increasing time.
,→ If an equilibrium state is Lyapunov stable, the state trajectory will remain within a
given neighbourhood if the initial state is chosen close to the equilibrium point.
Definition: Asymptotic stability
If in addition to satisfying the Lyapunov stability condition the trajectory starting within a
radius δ converges to the origin with continuing time, the system is said to be asymptotically
stable.
Definition: Asymptotic stability in the large
If a system is Lyapunov stable and all trajectories starting anywhere in the state space
converge to the origin with continuing time, the system is said to be asymptotically stable
in the large. Sometimes this referred to as global stability.
Note that for linear-time-invariant systems, asymptotic stability and asymptotic stability
in the large are the same.
Figure 2.4 shows three trajectories: trajectory 1 is a Lyapunov-stable trajectory, trajec-
tory 2 is asymptotically stable, while trajectory 3 is unstable.
x2
3 2−dimensional
state space
1 δ x1
R
PSfrag replacements
Lyapunov’s linearization method is concerned with the local stability of a system. To use
this method, it is assumed that the equilibrium state of the system is at the origin. The
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 40
nonlinear system is then approximated by a linear model (i.e. the system is linearized)
and the behaviour of the linear system analyzed.
Consider a scalar function f (h). For small x, Taylor’s approximation states that
∂f (h) 1 ∂ 2 f (h) 2
f (h + x) ≈ f (h) + x+ x +... (2.16a)
∂h 2! ∂h2
∂f (h)
≈ f (h) + x + higher order terms. (2.16b)
∂h
If h = 0, equation (2.16b) can be written as
∂f (h)
f (x) ≈ f (0) + x + higher order terms. (2.17)
∂h h=0
Equation (2.20) is used to approximate the nonlinear system with a linear one. Note
that since Taylor’s theorem is only valid for small variations around a point, the model
obtained using equation (2.20) is only useful in analyzing the stability in the neighbour-
hood of the equilibrium point at the origin.
Note that equation (2.20) is a vector matrix equation, and
∂f1 ∂f1 ∂f1
∂x1 ∂x2 · · · ∂xn
∂f2 ∂f2 ∂f2
∂f (x) ···
∂x1 ∂x2 ∂xn ,
= (2.21)
∂x ... .. .. ..
. . .
∂fn ∂fn ∂fn
···
∂x1 ∂x2 ∂xn
which is known as a Jacobian.
For a second order system, the approximation can therefore be written as
∂f1 ∂f1
" # " #
ẋ1 ∂x1 ∂x2 x1
≈ ∂f ∂f2 (2.22)
ẋ2 2 x2
∂x1 ∂x2 x1 =x2 =0
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 41
If the eigenvalues of the linearized system lie strictly on the left half plane, the
equilibrium point of the nonlinear system is asymptotically stable.
If at least one of the eigenvalues of the linearized system lies on the right half
plane, the equilibrium point of the nonlinear system is unstable.
If the eigenvalues of the linearized system lie on the left half plane, but at least
one of them lies on the imaginary axis, no conclusion can be drawn about the
stability of the equilibrium point of the nonlinear system.
Example
Consider the nonlinear system
The linearized system has eigenvalues on the right half plane so the equilibrium point
of the nonlinear system is unstable.
Let us define a scalar function V (x) (this is a scalar function of a vector x).
Definition: Positive definite function
A scalar function V (x) is positive definite in a region if
V (x) ≥ 0 if x 6= 0 (2.26a)
V (x) = 0 only when x = 0. (2.26b)
The examples are fairly simple, and it is fairly straight-forward to tell which functions
are positive/negative definite/semi-definite. However, in general, there is no general
way to tell if a function is positive or negative definite. However, if V (x) is in quadratic
form, we can use Sylvester’s criterion to determine if V (x) is positive definite.
Example
Use Sylvester’s criterion to show that the function below is positive definite:
Solution
For a 3rd order system, the quadratic form is given by
p11 p12 p13 x1
V (x) = xT Px = [x1 x2 x3 ] p12 p22 p23 x2
It follows that
10 1 −2
P= 1 4 −1 ,
−2 −1 1
and this matrix satisfies Sylvester’s criterion for positive-definiteness so the function
given in equation (2.30) is positive definite.
1. V (0) = 0.
3. V (x) is continuous and has continuous derivatives with respect to all components
of x.
4. V̇ (x) ≤ 0 (i.e. the time derivative of V (x) is less or equal to zero) along the
trajectories of the system.
• Failure to find a Lyapunov function for a particular system does not mean that a
Lyapunov function for that system does not exist, and in such cases, no conclusions
about the stability of the system can be made.
ẋ = Ax. (2.31)
Let V (x) be in quadratic form, i.e. V (x) = xT Px. By use of the chain rule, we can write
that
V̇ (x) = ẋT Px + xT Pẋ. (2.32)
Substituting (2.31) in (2.32) gives
where
Q = − AT P + PA .
(2.34)
Solution
Let
" # " #
p11 p12 1 0
P= (a real-symmetric matrix), and Q = . (2.36)
p12 p22 0 1
Then
" #" # " #" # " #
T
−1 1 p11 p12 p11 p12 −1 2 −1 0
A P + PA = + = −Q =
2 −1 p12 p22 p12 p22 1 −1 0 −1
(2.37)
or " # " #
−2p11 + 2p12 2p11 − 2p12 + p22 −1 0
= (2.38)
2p11 − 2p12 + p22 4p12 − 2p22 0 −1
which gives the simultaneous equations
Solving these equations gives p11 = −0.25, p12 = −0.75 and p22 = −1, hence
" #
−0.25 −0.75
P= . (2.40)
−0.75 −1
The matrix P above does not meet Sylvester’s criteria for positive-definiteness hence we
conclude that the given system is not stable, and a Lyapunov function for this system
does not exist.
Example
Repeating the above example for the system
" # " #" #
ẋ1 −5 −4 x1
= (2.41)
ẋ2 2 1 x2
gives " #
1 7
3 12
P= 7 11
. (2.42)
12 6
In this case, p11 > 0 and |P| > 0 (satisfying Sylvester’s criterion for positive-definiteness),
hence xT Px is a positive definite function and the system is asymptotically stable. The
Lyapunov function is given by
" #" #
1 7
x 1 1 2 7 11 2
V (x) = xT Px = [x1 x2 ] 37 12 11
= x 1 + x 1 x 2 + x2 . (2.43)
12 6
x 2 3 6 6
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 46
DIGITAL CONTROL
Digital control refers to the use of digital devices or digital computers to control dynamic
systems. The theory of digital control is very wide, and cannot be covered in a single
course. This chapter covers the basic materials of digital control systems. For more
details, the reader is referred to the references at the end of the handout.
PSfrag replacements
C(s)
G(s)
47
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 48
replaced with a digital controller. Since the digital controller works in discrete-time,
replacementssuch a system is also referred to as a discrete-time control system.
plant
r(t) e(t) sample e(kT ) u(kT )
D/A
u(t) y(t)
and A/D digital smoothing G(s)
+ hold converter controller converter filter
G(s)
Sample-and-Hold Circuit: used to hold the analog input constant when the Analog-
to-Digital converter is performing the conversion
Analog-to-Digital Converter: Converts the analog input into a binary (digital) code.
Note that many Analog-to-Digital converters have built-in sample-and-hold cir-
cuits.
Digital-to-Analog converter: Converts the control signal given in the form of a binary
code to an analog signal.
• Digital controllers have high repeatability since digital devices do not suffer the
drift associated with analog devices.
• It is easier and more economical to handle large control laws (i.e. high order
controllers).
Disadvantages:
• Mathematical analysis and design of digital control systems is generally more com-
plicated than the analysis of continuous-time control systems.
The z-transform serves the same role for discrete-time control systems that the s-transform
(Laplace transform) serves for continuous-time control systems.
Consider the continuous-time signal f (t) shown in Figure 3.3
f (t)
PSfrag replacements
0 t
Suppose the signal f (t) is sampled every T seconds using an ideal sampler. The ideal
sampler is represented as a switch as shown in Figure 3.4.
ideal sampler
interval T
The output of the ideal sampler can be thought of as a train of impulses whose respective
impulse areas are equal to the magnitudes of f (t) at the corresponding sampling instants
t = kT , k = 0, 1, 2, . . . (see Figure 3.5).
The output of the ideal sampler at various time instants is tabulated below:
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 50
PSfrag replacements
f ∗ (t)
0 T 2T 3T 4T 5T 6T t
X∞ Z ∞
= f (kT ) δ(t − kT )e−st dt. (3.3b)
k=0 0
The exponential term e−skT in equation (3.5 introduces nonlinear functions in the
Laplace transform of discrete-time systems. This problem is resolved by the change
of variable
1
z = esT =⇒ s= ln z. (3.6)
T
With this change of variable, equation (3.5) can then be written as
∞
1 X
F ∗ (s = ln z) = F (z) = f (kT )z −k (3.7)
T k=0
which involves positive and negative powers of z. Equation (3.9) is known as the two-
sided or bilateral z-transform.
In general, any continuous functions which have a Laplace transform also have a z-
transform. It follows that when working with discrete-time systems, we can solve the
system equations using Laplace transform techniques but with a change of variable from
s to z so that the computations remain algebraic.
For a function f (t), the functions f (kT ) and F (z) are said to be a z-transform pair, i.e.
z
f (kT ) ←→ F (z). (3.10)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 52
An important formula
The following formula from the theory of geometric series is useful in obtaining the z
transform of many functions:
N −1
2 N −1
X A − AxN
A + Ax + Ax + · · · + Ax = Axn = . (3.11)
n=0
1−x
Step function
(
1 t≥0
f (t) = (3.13)
0 t<0
hence f (kT ) = 1. It follows that
F (z) = 1 + z −1 + z −2 + · · · + z −k + · · · (3.14a)
1 z
= = for |z −1 | < 1 or |z| > 1. (3.14b)
1−z −1 z−1
In other words,
z
z
1 ←→ for |z| > 1. (3.15)
z−1
The region of convergence is the outside of a circle of unit radius centered at the origin.
Geometric sequence ak
∞ ∞
X
k −k
X k
F (z) = a z = az −1 (3.16a)
k=0 k=0
2 −2
= 1 + az + a z + · · · + ak z −k + · · ·
−1
(3.16b)
1 z
= = for |az −1 | < 1 (or z > a) (3.16c)
1 − az −1 z−a
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 53
∞ ∞ ∞
X X X k
F (z) = f (kT )z −k
= e z
−akT −k
= e−aT z −1 (3.18a)
k=0 k=0 k=0
= 1 + e−aT z −1 + e z +···+e z +···
−2aT −2 −akT −k
(3.18b)
1 z
= = for |e−aT z −1 | < 1 (or |z| > e−aT ).(3.18c)
1−e z−aT −1 z−e −aT
z
z ejωt |z| > ejωT ,
= for (3.21a)
z − ejωT
z
z e−jωt = for |z| > e−jωT , (3.21b)
z − e−jωT
hence
1 z z z sin ωT
z (sin ωt) = − = for |z| > 1. (3.22)
2j z − ejωT z − e−jωT z 2 − 2z cos ωT + 1
Similarly
z 2 − z cos ωT
z (cos ωt) = for |z| > 1. (3.23)
z 2 − 2z cos ωT + 1
Multiplication by a constant
If
z (x(t)) = X(z), (3.25)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 54
then
z (ax(t)) = aX(z). (3.26)
Proof: ∞ ∞
X X
z (ax(t)) = ax(k)z −k
=a x(k)z −k = aX(z). (3.27)
k=0 k=0
Linearity
If
x(k) = αf (k) + βg(k), (3.28)
then
X(z) = αF (z) + βG(z). (3.29)
Multiplication by ak
If
z (x(k)) = X(z), (3.30)
then
z ak x(k) = X(a−1 z).
(3.31)
Proof:
∞
X ∞
X z
k
k
−k
z a x(k) = a x(k)z −k
= x(k) a z
−1
= X(a z) = X
−1
. (3.32)
k=0 k=0
a
If
z (x(t)) = X(z), (3.33)
then
z (x(t − nT )) = z −n X(z). (3.34)
Proof: ∞ ∞
X X
z (x(t − nT )) = x(kT − nT )z −k
= x((k − n)T )z −k . (3.35)
k=0 k=0
A one sided z-transform of x(mT ) is zero for m < 0, hence the limit can be changed
from m = −n to m = 0, giving
∞
X ∞
X
z (x(t − nT )) = z −n x(mT )z −m = z −n x(mT )z −m = z −n X(z). (3.37)
m=−n m=0
Example
Find the z-transform of the signal f (t) shown in Figure 3.6.
f(t)
0 2T t
Figure 3.6:
Solution
Figure 3.6 shows a step signal which has been shifted to the right by two sample inter-
vals.
1
A step signal starting at t = 0, denoted as 1(t), has a z-transform of 1−z −1
.
Using the real translation theorem, it follows that the z-transform of f (t) = 1(t − 2T )
shown in Figure 3.6 is given by
z −2
.
1 − z −1
Example
Find the z-transform of the signal f (t) shown in Figure 3.7.
f(t)
0 2T 4T t
Figure 3.7:
Solution
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 56
z −2 z −4
F (z) = −
1 − z −1 1 − z −1
If
z (x(t)) = X(z), (3.38)
then
z e−at x(t) = X zeaT .
(3.39)
Proof:
∞
X ∞
X ∞
X −k
x(kT ) eaT z = X(zeaT ).
z e −at
x(t) = e −akT
x(kT )z −k
= x(kT ) e z
−akT −k
=
k=0 k=0 k=0
(3.40)
If
z (x(t)) = X(z), (3.41)
then
x(0) = lim X(z). (3.42)
z→∞
Proof:
∞
X
X(z) = x(kT )z −k = x(0) + x(T )z −1 + x(2T )z −2 + · · · + x(kT )z −k + · · · . (3.43)
k=0
If
z (x(k)) = X(z) (3.44)
where x(k) = 0 for k < 0, then
lim x(k) = lim 1 − z −1 X(z) . (3.45)
k→∞ z→1
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 57
Proof:
∞
X
z (x(k)) = X(z) = x(k)z −k , (3.46a)
k=0
∞
X
z (x(k − 1)) = z −1
X(z) = x(k − 1)z −k , (3.46b)
k=0
hence ∞ ∞
X X
x(k)z −k − x(k − 1)z −k = X(z) − z −1 X(z). (3.47)
k=0 k=0
But " #
∞
X ∞
X ∞
X ∞
X
lim x(k)z −k
− x(k − 1)z −k
= x(k) − x(k − 1). (3.49)
z→1
k=0 k=0 k=0 k=0
Recalling that x(k) = 0 for k < 0, we can expand the right hand side of equation (3.49)
to get
∞
X ∞
X
x(k) − x(k − 1) = [x(0) − x(−1)] + [x(1) − x(0)] + [x(2) − x(1)] + · · · = x(∞),
k=0 k=0
(3.50)
hence
lim x(k) = lim 1 − z −1 X(z) . (3.51)
k→∞ z→1
Introduction
This method involves the direct division of the numerator polynomial of X(z) by the
denominator polynomial of X(z). Before the division is done, both the numerator and
the denominator polynomials must be written in increasing powers of z −1 . By the divi-
sion, X(z) is expanded into an infinite power series in z −1 , and then the values of x(kT )
can be determined by inspection.
Example
Given
10z + 5
X(z) = , (3.52)
(z − 1)(z − 0.2)
find x(kT ) for k = 0, 1, 2, 3, 4.
Solution
First, X(z) is rewritten as a ratio of polynomials in z −1 :
10z −1 + 5z −2
X(z) = . (3.53)
1 − 2z −1 + 0.2z −2
Dividing the numerator by the denominator using long division yields
It follows that x(0) = 0, x(T ) = 10, x(2T ) = 17, x(3T ) = 18.4 and x(4T ) = 18.68.
Note that in general, the direct division method does not yield a solution in closed form.
The linearity property of the z-transform enables us to apply the partial fraction expan-
sion method to obtain the inverse z-transform of X(z).
Note that if X(z) has one or more zeros at the origin, it is easier to obtain the inverse
z-transform of X(z) by partial fraction expansion of X(z)
z
instead of the partial fraction
of X(z).
Example
Obtain the inverse z-transform of
10z
X(z) = . (3.55)
(z − 1)(z − 0.2)
X(z) has one zero at the origin. To show the difference between expanding X(z) and
X(z)
z
into partial fractions, both functions are expanded into partial fractions, starting
with X(z).
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 59
z −1 1
z −1
= z −1
z −1
(3.57a)
1 − z −1 1 − z −1
(
1 k = 1, 2, 3, . . .
= (3.57b)
0 k<1
Similarly, (
0.2k−1 k = 1, 2, 3, . . .
z −1
z−1 = (3.58)
1 − 0.2z −1 0 k<1
From (3.57) and (3.58) above,
X(z)
Let us now expand z
into partial fractions:
X(z) 10 1 1
= = 12.5 − . (3.60)
z (z − 1)(z − 0.2) z − 1 z − 0.2
Hence
z z 1 1
X(z) = 12.5 − = 12.5 − . (3.61)
z − 1 z − 0.2 1−z −1 1 − 0.2z −1
x(k) = 12.5 1 − (0.2)k
for k ≥ 0. (3.62)
It can be seen from above that if a function has one or more finite zeros, it is easier to
obtain the inverse z-transform by expanding X(z)
z
into partial fractions than expanding
X(z) into partial fractions.
3.3.1 Introduction
u(kT ) (sometimes written simply as u(k)) is a known input sequence, while y(kT ) (or
y(k)) is the output.
Consider a system described by the following linear difference equation:
Equation (3.63) (or (3.64)) is an nth order constant coefficient difference equation or an
nth order recursion equation.
y(kT ) (or y(k))is the value of the output at the current time, y(kT − T ) is the value of
the output at time (KT − T ) (the immediate past sampling instant), y(kT − 2T ) is the
value of the output at (kT −2T ) (the second-last sampling instant), and so on. Similarly,
u(kT ) is the value of the input at the current time, u(kT − T ) is value of the input at the
immediate past sampling instant, and so on.
Equation (3.63) (or (3.64)) can easily be solved using a digital computer by a procedure
known as recursion. To illustrate this procedure, suppose u(k) = 0 for k < 0 (i.e. the
input is zero for t < 0) as is usually the case.
Substituting k = 0 in (3.64), we can write
Recall that the input sequence is known, so, u(0) in the above equation is known. But
to determine y(0), the values y(−1), y(−2) · · · y(−n) must be known. These values are
known as the initial conditions of the difference equation (3.64).
Now substituting k = 1 into equation (3.64),
(Note that the value of y(0) determined from (3.65) is used in (3.66) above).
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 61
and so on. This procedure of solving a difference equation is known as recursion. The
procedure can easily be programmed into a digital computer.
Note that in the special case where u(k) = 0 for all k (input is identically zero), the
difference equation (3.64) is known as a homogeneous difference equation, and the
solution to the equation is known as the zero input response.
Also note that if the initial conditions y(−1), y(−2) · · · y(−n) are not known, we can
only obtain a general solution to the difference equation (3.64), not a particular solu-
tion.
or
1 + a1 z −1 + · · · + an z −n Y (z) = b0 + b1 z −1 + · · · + bm z −m U (z). (3.69)
Hence
Y (z) b0 + b1 z −1 + · · · + bm z −m
= G(z) = . (3.70)
U (z) 1 + a1 z −1 + · · · + an z −n
G(z) is known as the pulse transfer function of the discrete-time system. The pulse
transfer function relates the input u(kT ) to the output y(kT ) at the sampling instants
only. It does not give any information about the nature of output between the sampling
instants.
The inverse z-transform of G(z), i.e. z−1 [G(z)] = g(k), is known as the weighting se-
quence of the system.
Example
Obtain the weighting sequence of the system described by the difference equation:
Solution:
Taking the z-transform of (3.71) above gives
3.4.1 Introduction
jω
s = σ + jω
s−plane
0 σ
PSfrag replacements
Consider the s-plane shown in Figure 3.9. Every point in the s-plane can be expressed
as a sum of a real number σ and an imaginary number jω, i.e.
s = σ + jω. (3.76)
z = 1∠ωT, (3.80)
Equation (3.80) is the equation of a unit circle centered at the origin. The imaginary
axis on the s-plane is therefore mapped onto a unit circle centered at the origin in the
z-plane.
Note that the imaginary axis on the s-plane winds itself around the unit circle in the z-
plane infinitely many times - the origin in the s-plane is mapped to the point (1,0) in the
z-plane - then the positive part of the jω axis winds itself infinitely many times around
the unit circle in the anti-clockwise direction, while the negative part of the jω axis
winds itself around the unit circle in the clockwise direction. Points along the jω axis
in the s-plane whose frequencies differ by integral multiples of 2π T
are mapped to the
jπ j3π j5π jπ
same point in the z-plane, e.g. points T , T , T and − T are all mapped onto the point
1∠π = (−1, 0) in the z-plane.
Figure3.10 shows two lines of constant, attenuation in the s-plane. The line to the right
of the jω axis is mapped onto a circle of radius eσ1 T > 1 centered at the origin in the
z-plane. The line to the left of the jω axis is mapped onto a circle of radius eσ2 T < 1
centered at the origin in the z-plane as shown in Figure 3.11.
From the above, we can deduce that the left half plane in the s-plane is mapped onto the
inside of the unit circle centered at the origin in the z-plane, while the right half-plane
in the s-plane is mapped onto the outside of the unit circle in the z-plane
,→ A linear-time-invariant discrete-time system is stable if all its poles lie inside of the unit
circle centered at the origin in the z-plane.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 64
jω
s−plane
PSfrag replacements
σ2 0 σ1 σ
e σ1 T
0
PSfrag replacements
e σ2 T
in the z-plane. This is a straight line starting at the origin in the z-plane at an angle ω1 T
radians to the horizontal axis as shown in Figure 3.13. Similarly, the line at a frequency
jω2 is mapped to z = eσT ejω2 T in the z-plane, while the line at the frequency −jω1 is
mapped to z = eσT e−jω1 T in the z-plane.
Example:
Consider the region in the s-plane that is shaded in Figure 3.14.
The shaded region is mapped into the z-plane as shown in Figure 3.15.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 65
jω
s−plane
0 σ
−jω1
eσT
eσT
z−plane
PSfrag replacements
ω2 T
ω1 T
0 −ω1 T
eσT
jω
s−plane
jω
PSfrag replacements
2
0 σ
σ
σ
2 1
jω1
Figure 3.14:
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 66
eσT
z−plane
ω T
2
0 ω T
1
PSfrag replacements
eσT
unit circle
Figure 3.15:
Note that points in the s-plane with the same attenuation factors but whose frequencies
differ by integral multiples of 2π
T
are mapped onto the same locations in the z-plane.
This means that there are infinitely many values of s for each value of z, although a single
point in the s-plane corresponds to a single point in the z-plane.
3.5.1 Introduction
?
x(t) x(kT )
PSfrag replacements t t
Ideal Hold
sampler
The ideal sampler converts a continuous time signal to a train of pulses occurring
at the sampling instants.
The data hold is a device for generating a continuous-time signal from a discrete-
time sequence x(kT )
The signal h(t) during the time interval kT ≤ t ≤ (k + 1)T may be approximated by the
polynomial in τ as
If the data hold is an n-th order extrapolator, it is called an n-th order hold:
An n-th order hold uses the past (n + 1) discrete data samples x(kT − nT ), x(kT − nT +
T ), x(kT − nT + 2T ), . . . x(kT ) to generate the signal h(kT + τ ).
Because the higher data order hold uses past samples to extrapolate a continuous-time
signal between the present sampling instant and the next sampling instant, the accu-
racy of the approximation improves as the number of past samples used is increased.
However, this better accuracy is obtained at the cost of greater time delays since a large
number of past samples must be stored and the computation time increases with an
increase in the number of samples and this leads to increased time delays. In closed-
loop control systems, any added time delay in the loop will decrease the stability of the
system, and in some cases, it may lead to instability.
PSfrag replacements t t t
This implies that the data hold circuit holds the amplitude of the sample constant from
one sampling instant to the next. Such a data hold is known as a ZERO ORDER HOLD
(ZOH).
The ZOH holds the discretized input at a constant value throughout the sample period
until a new sample is obtained, at which time the output changes. This results in a
discontinuous signal with steps at each sampling instant. Sometimes, the output of a
ZOH is passed through a low-pass filter to remove the discontinuities.
For the ZOH,
It is known that
1
L [1(t)] = (3.90)
s
and from the properties of Laplace transform, we can write that
e−kT s
L [1(t − kT )] = (3.91)
s
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 69
hence
∞ ∞
e−kT s − e−(k+1)T s 1 − e−T s −kT s
X X
L(h(t)) = H(s) = x(kT ) = x(kT ) e (3.92)
k=0
s k=0
s
1 − e−T s
In other words, a ZOH block can be replaced by the transfer function .
s
1 − e−T s
3.5.3 z-transform of functions involving the term
s
Suppose a transfer function G(s) follows a ZOH as shown in Figure 3.18.
PSfrag replacements ZOH System
1 − e−T s
G(s)
s
The overall transfer function is the product of the two transfer functions, and can be
written as:
1 − e−T s G(s)
X(s) = G(s) = 1 − e−T s (3.94)
s s
The z-transform of the overall transfer function X(s), denoted by X(z), is given by
G(s)
X(z) = z [X(s)] = 1 − z −1
z . (3.95)
s
Proof:
The overall transfer function is given by
1 − e−T s G(s)
G(s) = 1 − e−T s = 1 − e−T s G1 (s)
X(s) = (3.96)
s s
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 70
G(s)
where G1 (s) = .
s
Let X1 (s) = e−T s G1 (s). Since X1 (s) is a product of two Laplace-transformed functions,
the inverse Laplace transform is the following convolution integral:
Z ∞
x1 (t) = g0 (t − τ )g1 (τ )dτ (3.97)
0
−T s
where g0 (t) = L−1 e = δ(t − T ) and g1 (t) = L−1 (G1 (s)), hence
Z ∞
x1 (t) = δ(t − T − τ )g1 (τ )dτ = g1 (t − T ). (3.98)
0
or
G(s)
X(z) = z [X(s)] = 1 − z −1
z (3.100)
s
End of Proof
PSfrag replacements
3.6 Block diagram analysis
Figure 3.19:
Proof:
Converting equation (3.101) into time domain,
From the properties of the Dirac function, the above equation reduces to
∞
X
y(t) = g(t − nT )x(nT ). (3.106)
n=0
At time t = kT ,
∞
X
y(kT ) = g(kT − nT )x(nT ). (3.107)
n=0
x(t) y(t)
G(s)
X(s) Y (s)
Figure 3.20:
hence
Y (z) = z [GX(s)] . (3.113)
Note that the z-transform for Y (s) in Figure 3.19 (equation (3.103))is different from the
z-transform of Y (s) in Figures 3.20 (equation (3.113)): for the system in Figure 3.19,
Y (z) is the product of the z-transform of X(s) and the z-transform of G(s), while for
the system in Figure 3.20, Y (z) is the z-transform of the product G(s)X(s) (abbreviated
as GX(s)).
Example:
Suppose that in both Figures 3.19 and 3.20,
1 1
X(s) = (step input), and G(s) = .
s s+a
For the system in Figure 3.19, Y (z) = G(z)X(z).
1 1 1 1
G(z) = z = and X(z) = z = ,
s+a 1 − z e−aT
−1 s 1 − z −1
hence
1
Y (z) = .
(1 − z −1 e−aT ) (1 − z −1 )
Consider the system shown in Figure 3.21. It is assumed that all the samplers are
synchronized.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 73
Figure 3.21:
But
U (s) = G(s)X ∗ (s) =⇒ U ∗ (s) = [G(s)X ∗ (s)]∗ = G∗ (s)X ∗ (s), (3.116)
hence
PSfrag replacements
Y ∗ (s) = H ∗ (s)G∗ (s)X ∗ (s) (3.117)
or
Y (z) = H(z)G(z)X(z). (3.118)
The pulse transfer function of the block diagram shown in Figure 3.21 is
Y (z)
= H(z)G(z). (3.119)
X(z)
Consider the system shown in Figure 3.22. Again, it is assumed that all the samplers
are synchronized.
Figure 3.22:
It follows that
Y (z) = GH(z)X(z), (3.122)
and the pulse transfer function of the system is
Y (z)
= GH(z). (3.123)
X(z)
H(s)
Figure 3.23:
hence
R∗ (s)
E ∗ (s) = . (3.128)
1 + (GH(s))∗
But Y (s) = G(s)E ∗ (s), so
Y ∗ (s) = G∗ (s)E ∗ (s). (3.129)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 75
H(s)
Figure 3.24:
E(s) = R(s) − H(s)G∗ (s)E ∗ (s) =⇒ E ∗ (s) = R∗ (s) − H ∗ (s)G∗ (s)E ∗ (s). (3.134)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 76
Figure 3.25:
To obtain the pulse transfer function of the closed-loop system, we first replace the ZOH
block with its transfer function as shown in Figure 3.26.
G1 (s)
Figure 3.26:
G1 (s) is the combined transfer function of the ZOH and G(s), and the z-transform of
G1 (s) can then obtained as described in Section 3.5.3, i.e.
G(s)
G1 (z) = z [G1 (s)] = 1 − z −1
z (3.138)
s
The block-diagram of the closed-loop system reduces to what is shown in Figure 3.27.
The pulse transfer function of the closed-loop system is
Y (z) G1 (z)
= . (3.139)
R(z) 1 + G1 (z)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 77
Figure 3.27:
Exercise
For the system shown in Figure 3.25, suppose
K
G(s) =
s+1
where K is a constant.
(ii) Show that the pulse transfer function of the closed-loop system is
Y (z) Kz −1 1 − e−T
=
R(z) (1 − z −1 e−T ) + Kz −1 (1 − e−T )
replacements
Figure 3.28 shows a typical digital control system.
plant
r(t) e(t) sample e(kT ) u(kT )
D/A u(t) y(t)
and A/D digital
G(s)
+ hold converter controller converter
G(s)
To obtain the pulse transfer function of the closed-loop system, the sample-and-hold-
A/D-converter block is replaced with an ideal sampler of sampling interval T , the digital
controller is represented using a transfer function C(z) (or C ∗ (s)). Unless otherwise
specified, the D/A converter is assumed to be a Zero-Order-Hold (ZOH).
The new block diagram of the system is shown in Figure 3.29.
The block diagram shown in Figure 3.29 can be further reduced as shown in Figure 3.30,
rag replacements
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 78
G1 (s)
Figure 3.29:
where
1 − eT s
G1 (s) = G(s).
s
PSfrag replacements
Figure 3.30:
Also,
E(s) = R(s)−Y (s) =⇒ E ∗ (s) = R∗ (s)−Y ∗ (s) = R∗ (s)−E ∗ (s)C ∗ (s)G∗1 (s) (3.141)
hence
R∗ (s)
E ∗ (s) = . (3.142)
1 + G∗ (s)C ∗ (s)
Substituting (3.142) into (3.140) gives
Exercise
For the system shown in Figure 3.29, suppose
K
G(s) =
s+1
where K is a constant. Show that the pulse transfer function of the closed-loop system
is
Y (z) Kz −1 1 − e−T
=
R(z) 1 − z −1 (K − Ke−T − e−T − 1) − z −2 e−T
3.7.1 Introduction
Sampling is made necessary because of the discrete nature of digital systems. Sampling
means that the continuous-time signal is replaced by a train of impulses. This handout
examines if it is possible to make an exact reconstruction the original signal from its
samples.
Consider a signal which has a Fourier transform that is zero outside the interval (−B, B)
Hz (or (−2πB, 2πB) radians per second). Such a signal is said to be bandlimited to B
(see Figure 3.31).
PSfrag replacements
−B 0 B f (Hz)
Consider the continuous-time signal g(t) shown in Figure 3.32. Assume that this signal
is bandlimited.
g(t)
PSfrag replacements
0 t
PSfrag replacements
Figure 3.32: A bandlimited signal
Also consider a train of unit impulses shown in Figure 3.33. It can be written that
∞
X
s(t) = δ(t − kT ). (3.145)
k=−∞
s(t)
Multiplication of g(t) by the train of unit impulses s(t) yields the sampled signal gs (t)
shown in Figure 3.34. It follows that
∞
X
gs (t) = g(t) · δ(t − kT ). (3.146)
k=−∞
We now obtain the Fourier transforms of the signals g(t), s(t) and gs (t) above. It was
assumed that the signal g(t) is bandlimited and we can assume without loss of gen-
erality that its Fourier transform G(ω) is bandlimited to 2πB radians per second, with
a maximum amplitude A, and its Fourier transform G(ω) can be drawn as shown in
Figure 3.35.
PSfrag replacements
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 81
gs (t)
0 T 2T 3T 4T 5T 6T t
PSfrag replacements
Figure 3.34: The sampled signal
The Fourier
PSfrag transform of the train of unit impulses (s(t)) is given by
replacements
∞ ∞
2π X 2πk 2π X
F (s(t)) = S(ω) = δ(ω − )= δ(ω − kωs ) (3.147)
T k=−∞ T T k=−∞
where ωs = 2π T
(this is derived using the theory of Fourier transforms for periodic sig-
nals). This signal is illustrated in Figure 3.36
S(ω)
Recall that multiplication in the time domain is equivalent to convolution in the fre-
quency domain. Formally, this can written as
1 1
Z ∞
F
gs (t) = s(t) · g(t) ←→ Gs (ω) = [S(ω) ∗ G(ω)] = S(ω − τ )G(τ )dτ (3.148)
2π 2π −∞
written that ∞
A X
Gs (ω) = G(ω − kωs ) (3.149)
T k=−∞
2π
where ωs = T
.
If ws > 2(2πB), there is no overlap in the shifted replicas of G(ω) and hence g(t) can
recovered exactly from gs (t) by means of a low pass filter with a gain T and a cut-
off frequency ωc greater than 2πB and less than (ωs − 2πB). This proves Shannon’s
sampling theorem.
1
The maximum allowable sampling interval Ts = 2B is known as the Nyquist interval.
The corresponding sampling rate (2B samples per second) is known as the Nyquist
sampling rate.
3.7.4 The
PSfrag ideal low pass filter
replacements
The ideal low pass filter required for Shannon reconstruction is shown in Figure 3.38.
H(ω)
T primary
component
This ideal low pass filter shown in Figure 3.38 has a transfer function
(
T |ω| ≤ ωc
H(ω) = (3.150)
0 |ω| > ωc
where ωc is the cut-off frequency of the filter. This filter will perform exact interpolation
of the signal between the sampling instants.
Performing the inverse Fourier transform on H(ω) gives
Z ωc
1 1
Z ∞
jωt
T ejωt dω
h(t) = H(ω)e dω = (3.151)
2π −∞ 2π −ωc
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 83
1
T ωc
π
0.8
0.6
0.4
0 t
−0.2
−0.4
−5 −4 −3 −2 −1 0 1 2 3 4 5
From Section 3.5.2, it was shown that the transfer function of the ZOH is given by
1 − e−sT
Gh (s) =.
s
Replacing s with jω in the above equation gives
jωT
! jωT jωT
!
1 − e−T s
−jωT −jωT
1 − e−jωT jωT e 2 −e 2 2e 2 e 2 −e 2
Gh (jω) = = = 2e 2 =
s s=jω jω 2jω ω 2j
ωT
jωT sin
2e 2 ωT jωT 2
= sin = Te 2 . (3.153)
ω 2 ωT
2
The magnitude response of the ZOH is given by
ωT ωT
sin sin 2
jωT 2
|Gh (jω)| = T e 2 =T . (3.154)
ωT ωT
2 2
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 84
The magnitude curve of the ZOH is shown in Figure 3.40. The frequency response of
the ideal low pass filter is drawn in dotted lines for comparison.
Gh (ω)
1 T
0.8
0.6
0.637T
PSfrag replacements
0.4
0.2 0.212T
0.127T
0
− 6π − 4π − 2π 0 2π 4π 6π ω
T T T T T T
−0.2
−20 −15 −10 −5 0 5 10 15 20
From the figure, it can be seen that the frequency response of the Zero-Order-Hold is
quite different from the frequency response of the ideal low pass filter. In particular,
the spectrum of the ZOH does not go to zero beyond the frequency ω = ω2s = Tπ , and
this causes distortion in the signal reconstructed using the ZOH. However, the ZOH is
physically realizable and is simple to implement and for these reasons it is commonly
used in digital control systems.
When sampling is done at frequencies higher than 2B samples per second, the frequency
spectrum of a sampled bandlimited signal is as shown in Figure 3.37. However, if the
sampling frequency is less than 2B samples per second, the frequency spectrum of the
sampled signal looks as shown in Figure 3.41.
In this case, the replicated spectra of G(ω) overlap and g(t) is no longer recoverable by a
low pass filter. This overlapping effect is known as aliasing. (One can think of aliasing as
follows: from Figure 3.41, it can be seen that the high frequency tail of the component
2 overlaps with a lower frequency part of component 1, hence we say that the higher
frequencies of component 2 take on an alias of a lower frequency of component 1).
PSfrag replacements
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 85
A
T
2 1
When aliasing takes place, the reconstructed signal will have distortion for two reasons:
• Loss of the high frequency tail of Gs (ω) beyond |ω| > ωs /2.
• The same high frequency tail appears inverted or folded onto the spectrum at
frequencies |ω| < ωs /2.
This tail inversion is known as spectral folding and hence aliasing is sometimes referred
to as frequency folding.
In practice, aliasing always occurs regardless of the sampling frequency since bandlim-
ited signals do not exist in reality. It can be shown that if a signal is bandlimited, it
cannot be time-limited and vice-versa. All physical signals exist over a finite time dura-
tion hence they are not bandlimited.
To reduce the effects of aliasing, the following measures are usually adopted:
• Filtering the signals before sampling (this reduces the high frequency tail).
In practice, the rule of thumb is to choose the sampling time T such that, either
1
= 30 × closed-loop bandwidth (3.155)
T
or
smallest time constant of system
T = . (3.156)
20
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 86
3.8.1 Introduction
Discrete-time controller design methods can be broadly divided into two categories:
Discrete-time controller design only considers the system behaviour at the sam-
pling instants. Problems of inter-sampling behaviour can only be solved if the
continuous-time controller is designed first and then discretized.
In the discretization, one may wish to preserve some of the following properties:
d.c. gain.
Not all these properties can be preserved at the same time after discretization and the
choice of the discretization method is normally based on what the designer wants to
preserve.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 87
Introduction
Cross-multiplying gives
sU (s) + aU (s) = bE(s). (3.158)
Converting (3.158) into time domain gives the differential equation:
Similarly, at time t = kT − T ,
Z kT −T
u(kT − T ) = (−au(t) + be(t)) dt. (3.162)
0
The numerical integration methods depend on how the integral on the right of equation
(3.163) is approximated, i.e. how the area under the curve −au(t) + be(t) in the region
between t = (kT − T ) and t = kT is approximated.
By the forward rectangular rule, the area between t = (kT − T ) and t = kT is approxi-
mated by
Z kT
(−au(t) + be(t)) dt ≈ T (−au(kT − T ) + be(kT − T )) . (3.164)
KT −T
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 88
−au(t)+be(t)
kT−T kT t
hence
U (z) bT z −1 b
= = . (3.167)
E(z) 1 − z + aT z
−1 −1 1−z −1
+a
T z −1
Comparing (3.167) with (3.157) shows that
1 − z −1 z−1
s= = . (3.168)
Tz −1 T
By the backward rectangular rule, the area between t = (kT − T ) and t = kT is approx-
imated by
Z kT
(−au(t) + be(t)) dt ≈ T (−au(kT ) + be(kT )) . (3.170)
KT −T
Taking z-transforms,
U (z) = z −1 U (z) − aT U (z) + bT E(z) = z −1 − aT U (z) + bT E(z) (3.172)
hence
U (z) bT b
= = . (3.173)
E(z) 1 − z −1 + aT 1−z −1
+a
T
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 90
z−plane
Figure 3.45:
The entire left half of the s-plane is mapped into a circle of radius 21 , centered at ( 12 , 0) in
the z-plane, (see Figure 3.45), hence stable s-plane controllers produce stable z-plane
controllers (it is even possible for an unstable s-plane controller to be mapped into
the stable region in the z-plane). However, since the entire left half of the s-plane is
squeezed into a circle within the unit circle in the z-plane, there is a lot of distortion in
the transient and frequency response of the controller obtained using this method. To
minimize the distortion, a small sampling time should be used.
U (z) b
= , (3.180)
E(z) 2 1 − z −1
+a
T 1 + z −1
hence
2 1 − z −1 2 z−1
s= = (3.181)
T 1 + z −1 T z+1
Equation (3.181) is known by several names: trapezoidal rule, Tustin’s method, Tustin’s
transformation or bilinear transformation.
For this method,
sT
1+
z= 2 , (3.182)
Ts
1−
2
hence for s = jω, |z| = 1. The left half s-plane is mapped into the unit circle in the
z-plane, hence stable s-plane controllers will produce stable z-plane controllers.
For this method, the impulse response of the continuous-time controller is equal to the
impulse response of the discrete-time controller at the sampling instants.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 92
Given a continuous controller G(s), the discrete-time equivalent using the impulse in-
variance method is given by
For the step-invariance method, the step responses of the continuous and the discrete-
time controllers agree at the sampling instants.
Given a continuous-time controller G(s), the discrete-equivalent using the step-invariance
method is given by
1 − e−T s
G(s)
G(z) = z G(s) = 1−z −1
z (3.184)
s s
(compare equation (3.184) with equation (3.95) – the controller discretized using this
method is sometimes referred to as the ZOH equivalent).
Using this method, the poles and finite zeros of the continuous controller are mapped
onto the z-plane using the rule z = esT . This method is easy to apply if the transfer
function of the continuous controller is given in poles and zeros.
The complete procedure of discretizing a continuous-time controller using the matched
pole-zero mapping method is as follows:
1. The controller G(s) must be factored before this method can be applied.
2. The poles of G(s) are mapped onto the z-plane according to the relation z = esT ,
e.g. a pole at s = −a is mapped to z = e−aT .
5. The gain of the discrete-time controller is adjusted to match the gain of the con-
tinuous controller.
If the controller is a low-pass filter, the gain of the discrete-time controller at z = 1
should be the same as the gain of the continuous controller at s = 0.
For a high pass filter, match the gains at z = −1 and s = ∞.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 93
Example
Use the matched pole-zero mapping method to obtain the discrete-time equivalent of
the continuous-time controller
a
G(s) = . (3.185)
s+a
Solution
The given controller is a low-pass filter with one pole at s = −a and one zero at infinity.
The pole and zero are mapped to z = e−aT and z = −1 respectively, giving a discrete
time equivalent
z+1
G(z) = K , (3.186)
z − e−aT
where K is the gain.
Since this is a low-pass filter,
G(z)|z=1 = G(s)|s=0 , (3.187)
which gives
1 − e−aT
K= , (3.188)
2
hence
1 − e−aT (z + 1)
G(z) = . (3.189)
2 (z − e−aT )
Exercise
Use the matched pole-zero mapping method to show that the discrete-time equivalent
of the continuous-time controller
s+b
G(s) = , (3.190)
s+a
where our interest is in the low-frequency range, is given by
b 1 − e−aT z − e−bT
G(z) = (3.191)
a (1 − e−bT ) (z − e−aT )
Example
Use the matched pole-zero mapping method to obtain the discrete-time equivalent of
the high-pass filter
s
G(s) = . (3.192)
s+a
Solution
The pole at s = −a and the zero at s = 0 are mapped to z = e−aT and z = 1 respectively,
giving a discrete time equivalent
z−1
G(z) = K , (3.193)
z − e−aT
where K is the gain.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 94
3.9.1 Introduction
where (3.197a) is the state equation while (3.197b) is the output equation.
If these equations can be linearized, then it can be written that
The state equation (3.200a) gives the next state x(kT + T ) (or x(k + 1)) given the
current state x(kT ) and input u(kT ), while the output equation (3.200b) gives the
current output y(kT ) given the current state and input.
In some books, equation (3.200) is written as
x(k) + y(k)
u(k) x(k+1) +
H z I
−1 C
+
+
Figure 3.47 is similar to the one for a continuous-time LTI state-space model encoun-
tered earlier, the only difference being that the integrator in the continuous-time model
is replaced by the z −1 operator in the discrete-time model. The z −1 operator is known as
the unit delay operator and it is analogous to the integrator in continuous-time models.
Introduction
Just as in continuous-time systems, discrete time models may be expressed in the canon-
ical forms, namely
Y (z) b2 + b1 z −1 + b0 z −2
= . (3.204)
U (z) 1 + a1 z −1 + a0 z −2
The pulse transfer function (3.204) can be represented as shown in Figure 3.48
Figure 3.48:
hence
X(z) = U (z) − a1 z −1 X(z) + a0 z −2 X(z). (3.206)
Also,
Y (z) = b2 X(z) + b1 z −1 X(z) + b0 z −2 X(z). (3.207)
Equations (3.206) and (3.207) can be modelled as shown in Figure 3.49.
With the states labelled as shown in Figure 3.49, the state-space representation in con-
trol canonical form is:
" # " #" # " #
x1 (k + 1) 0 1 x1 (k) 0
= + u(k) (3.208a)
x2 (k + 1) −a0 −a1 x2 (k) 1
" #
h i x (k)
1
y(k) = b0 − a 0 b0 b1 − a 1 b2 + [b2 ] u(k) (3.208b)
x2 (k)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 97
PSfrag replacements b2
b1
+
U (z) X(z) x2 x1 +
+ Y (z)
+ + z −1 z −1 b0
+
−a1
−a0
With the states labelled as shown in Figure 3.50, the state-space representation in ob-
servable canonical form is:
" # " #" # " #
x1 (k + 1) 0 −a0 x1 (k) b0 − a 0 b0
= + u(k) (3.211a)
x2 (k + 1) 1 −a1 x2 (k) b1 − a 1 b2
" #
h i x (k)
1
y(k) = 0 1 + [b2 ] u(k) (3.211b)
x2 (k)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 98
U (z)
b0 b1 b2
−a0 −a1
where n > m (a strictly proper pulse transfer function). Assuming that the characteristic
equation of the pulse transfer function (3.212) has n distinct poles p1 , p2 , . . . , pn , the
pulse transfer function can be expanded using partial fraction expansion as
Y (z) c1 c2 cn
= + +···+ . (3.213)
U (z) z − p1 z − p2 z − pn
PSfrag replacements
A block diagram of the system is shown in Figure 3.51.
1 X1 (z)
c1
z − p1
1 Xn (z)
cn
z − pn
X1 (z) 1
= . (3.214)
U (z) z − p1
Cross-multiplying gives
zX1 (z) − p1 X1 (z) = U (z). (3.215)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 99
Similarly
x2 (k + 1) = p2 x2 (k) + u(k)
x3 (k + 1) = p3 x3 (k) + u(k)
.. . (3.217)
.
xn (k + 1) = pn xn (k) + u(k)
The output equation is
In general, discrete-time equations are easier to solve than differential equations be-
cause difference equations can be solved easily by means of a recursion procedure. The
recursion procedure can easily be programmed for solution with a digital computer.
Consider the LTI discrete-time system
Equation (3.221a) can be solved for any positive integer k using a recursion procedure
by substituting k = 0, k = 1, k = 2 and so on.
Substituting k = 0 in equation (3.221a) gives
For k = 1,
x(2) = Gx(1) + Hu(1). (3.223)
Substituting (3.222) in (3.223) gives
Similarly
x(3) = G3 x(0) + G2 Hu(0) + GHu(1) + Hu(2) (3.225)
and so on. By repeating this procedure, we obtain
k−1
X
k
x(k) = G x(0) + Gk−j−1 Hu(j). (3.226)
j=0
x(k) consists of two parts, the term Gk x(0) which is the contribution by the initial state
Pk−1 k−j−1
x(0), and j=0 G Hu(j) which is the contribution by the input u(j). The term
G x(0) is also referred to as the zero input response while the term k−1
k k−j−1
P
j=0 G Hu(j)
is referred to as the zero initial state response. The continuous-time analogy to equation
(3.226) is Z t
x(t) = eAt x(0) + eA(t−τ ) Bu(τ )dτ (3.227)
0
Rt
where eAt x(0) is the zero input response and 0
eA(t−τ ) Bu(τ )dτ is the zero initial state
response.
Substituting (3.226) in (3.221b) gives
k−1
X
k
y(k) = CG x(0) + C Gk−j−1 Hu(j) + Du(k). (3.228)
j=0
Ψ = Gk . (3.229)
z transform approach
Rearranging,
(zI − G) X(z) = zx(0) + HU(z). (3.232)
Pre-multiplying both sides of (3.232) by (zI − G)−1 gives
and
k−1
X
Gk−j−1 Hu(j) = z−1 (zI − G)−1 HU(z) .
(3.235)
j=0
Solution
Gk = z−1 (zI − G)−1 z .
Obtaining Gk involves finding the inverse z-transform of (zI − G)−1 z. Recall that if
a function X(z) has one or more zeros at the origin, it is easier to obtain the inverse
z-transform by first expanding X(z)z
into partial fractions then multiplying the result by
z (see Section 3.2.4). For that reason, we first expand (zI − G)−1 into partial fractions.
" #−1
z −1
(zI − G)−1 =
0.16 z + 1
z+1 1
(z + 0.2)(z + 0.8) (z + 0.2)(z + 0.8)
=
−0.16 z
(z + 0.2)(z + 0.8) (z + 0.2)(z + 0.8)
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 102
4 1 5 5
3 3 3 3
z + 0.2 − z + 0.8 −
z + 0.2 z + 0.8
=
.
0.8 0.8 1 4
− 3 + 3
− 3 + 3
z + 0.2 z + 0.8 z + 0.2 z + 0.8
Multiplying the above result by z gives
4 1 5 5
3
z 3
z 3
z 3
z
−
z + 0.2 z + 0.8 −
z + 0.2 z + 0.8
−1
(zI − G) z =
,
0.8 0.8 1 4
z z z z
− 3 + 3 − 3 + 3
z + 0.2 z + 0.8 z + 0.2 z + 0.8
hence
4 1 5 5
(−0.2)k − (−0.8)k k
(−0.2) − (−0.8) k
3 3 3 3
Gk = z−1 (zI − G)−1 z =
0.8 0.8 1 4
− (−0.2)k + (−0.8)k − (−0.2)k + (−0.8)k
3 3 3 3
which is the state transition matrix.
Recall that the transfer function relates the inputs to the outputs for zero initial condi-
tions hence we can ignore x(0) in equation (3.237a) to write
The input to the plant comes from the ZOH. Recall that the output of the ZOH is constant
between sampling instants. Between the sampling instant kT and kT + T , the output of
the ZOH is u(kT ). Substituting this in (3.246) above gives
Z kT +T
AT
x(kT + T ) = e x(kT ) + eA(kT +T −τ ) Bu(kT )dτ. (3.247)
kT
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 104
or Z T
AT
x(kT + T ) = e x(kT ) + eAλ Bu(kT )dλ. (3.249)
0
G = eAT (3.251)
and Z T
H= eAλ Bdλ. (3.252)
0
Equation (3.251) is used to compute the state matrix of the discrete-time system shown
in Figure 3.52, while equation (3.252) is used to determine the input matrix. As can
be seen from equation (3.250), the output and the direct transmission matrices C and
D respectively remain unchanged. This model is known as the ZOH equivalent of the
continuous-time plant. The discrete-time model given by equations (3.249) and (3.250)
is also the step-invariant discrete-time equivalent of the continuous model (3.241).
In MATLAB, the ZOH equivalent of a continuous-time plant is obtained using the com-
1
mand c2d. Figure 3.53 shows the step response of the continuous-time system s2 +0.4s+1
(plotted using the continuous line), and the step response for the ZOH equivalent dis-
crete model (plotted with the asterisks). As can be seen from the figure, the two plots
agree at the sampling instants.
ELC 514E: STATE SPACE DESIGN AND DIGITAL CONTROL by Dr.-Ing. S. I. Kamau 105
1.6
1.4
1.2
0.8
0.6
0.4
0.2
0
0 2 4 6 8 10 12 14 16 18 20
Ts = 1; %sampling time
sysd = c2d(sys1,Ts); %discretization - ZOH equivalent
td =0:Ts:20; %sampling times for discrete model
N = length(td); %total number of sample points
AD = sysd.A; %state matrix for discrete model
BD = sysd.B; %input matrix for discrete model
CD = sysd.C; %output matrix for discrete model
DD = sysd.D; %feedthrough matrix for discrete model
yd = dstep(AD,BD,CD,DD,1,N); %step response for discrete model
figure
plot(t,y,td,yd,’*’),grid, %plotting the results
=== × === × === × === × === × === × === × === × === × === × ===
Bibliography
[4] N. S. Nise, Control Systems Engineering, Fourth Edition, John Wiley & Sons
Inc., 2004.
SK ’07
T HIS DOCU MEN T WAS PREPARED U SIN G T HE LATEX
WORD − PROCESSOR RU N N IN G ON A LIN U X PLAT F ORM.
PICT U RES GEN ERAT ED U SIN G T HE Xfig PROGRAM.
106