Description
1 Introduction
At the end of the previous lab, you arrived at various representations of a linearized model of a control
system, including transfer function form, zero-pole form, and state space form. In this course, our
primary focus is on the analysis of linear systems in state space form, which strongly relies on the
machinery of linear algebra. In this lab you will review the basic techniques of numerical linear algebra,
and then apply those techniques to the analysis of a simple linear system.
Throughout the lab, you will be guided through a number of steps which will require you to write
Matlab code. You will write your code in a Matlab script called labx.m, where x is the lab number. You
will submit this code as a group. Your code should provide certain outputs (figures, numbers, and the
like). A request for output is highlighted with a shaded area of text, such as the following.
1
Output 1. Print the eigenvalues of the matrix.
Parts of the text containing directions for the writing of your Matlab code will be highlighted in a
different colour, such as this:
Create a function pendulum.m. The first line of the function will be
function xdot = pendulum(t,x,parameters)
2 Submission guidelines and mark breakdown
Marks are assigned to groups. Members of each group get identical marks, unless special circumstances
occur. The group lab mark will be made up of three components.
Matlab code 6 pts
Lab report 4 pts
Total 10 pts
Matlab code. Your code should be clean and readable, and it should contain abundant commentary,
so that the instructors may follow your development and check its correctness. This component of the
mark will be based on the correctness of the code and its readability, and it will be assigned as follows:
0 out of 6 the code is absent
1 out of 6 the code is largely incomplete
2 out of 6 the code is largely complete, but there are parts missing and the outputs are incorrect
3 out of 6 the code is complete, but it does not produce correct outputs
4 out of 6 the code is complete, it produces correct outputs, but there is no commentary in the
code, and/or the code is poorly organized
5 out of 6 the code is complete, correct, and contains some but insufficient commentary,
and/or its organization is below par
6 out of 6 the code is complete, correct, and the commentary and organzation are adequate so
that it is easy to read
Lab report. Write a concise report containing the requested outputs, but don’t just print out a bunch
of matrices and figures. Add some flesh to the outpus that are requested within the lab document so
that one can follow the logical development of the lab. Aim for a style like this:
Output 1. The following are the bases of the subspaces V and W that were obtained:
(…)
We observe that the columns of V were already linearly independent, while those of W weren’t.
We further note that (…).
2
Output 2. Below the solution of the differential equation. There are three figures. The first two
figures depict x1(t) and x2(t), while the third figure depicts the orbit.
Do not screenshot Matlab figures. Rather, save them as jpeg files (or other formats) and include then
in the report in the appropriate order.
When the lab document asks you to comment on something, strive to provide meaningful, insightful
commentary, that demonstrates you have understood your own work and how it connects to the course
concepts. This portion of the mark will be assigned as follows:
0 out of 4 the report is either absent, or a simple printout of the Matlab outputs without
commentary
1 out of 4 the report is incomplete and the commentary is inadequate
2 out of 4 the report is incomplete xor the commentary is inadequate
3 out of 4 the report is complete and the commentary is adequate, but lacks insight/clear
understanding
4 out of 4 the report is complete and the commentary is adequate and demonstrates clear
understanding
3 Numerical linear algebra
In this lab, it is assumed that you have read the linear algebra chapter in the textbook, and have reviewed
the lecture notes on this subject. If you haven’t, this is a good time to do it.
3.1 Basic operations on subspaces
In this section, you will use the following Matlab commands.
Matlab commands
orth(A) Find an orthogonal basis for the column space of A
null(A) Find an orthogonal basis for the null space of A
inv(A) Find the inverse of A
Let X be an n-dimensional vector space, and let V, W be two subspaces of X . Two key operations
involving these subspaces are the subspace sum V + W and subspace intersection V ∩ W. In this section
we will explore how to represent these subspaces numerically, and algorithms for performing these two
operations.
Let A ∈ Rn×m be a matrix. Recall that the image of A is simply the span of the columns of A
Im(A) = span{a
1
, . . . , a
n
}
where a
i
is the ith column of A. This fact motivates us to use matrices to represent subspaces. In order
to computationally convenient and remove redundancy, we will need a procedure to find a basis for the
image of a matrix. Matlab’s orth function does exactly this.
3
Consider R4 and the subspaces
V = span
1
−1
0
1
,
1
1
0
0
,
3
1
0
1
and W = span
1
0
2
1
,
1
0
−2
0
Using the Matlab orth function, begin your lab2.m script by finding bases for these subspaces.
Output 1. Print the matrices containing the bases of these subspaces. Comment on any interesting
details. Are the given vectors linearly independent? What has happened to the basis for W?
We now want to develop procedures for finding the subspace sum and intersection. First we address
the subspace sum. Recall the definition:
V + W = {x ∈ X : x = v + w, v ∈ V, w ∈ W} .
If we have V = Im(V) and W = Im(W) for suitable matrices V and W, then
V + W = Im V W.
Using this fact, do the following.
Create a function sub_sum.m. The first line of the function will be
function C = sub_sum(A, B)
The input matrices A and B contain bases for the subspaces we wish to sum. The matrix C that is to be
returned should contain a basis for the sum of the subspaces.
The procedure for finding the intersection is slightly more involved. Suppose that V and W are land m-dimensional subspaces respectively, so that V ∈ Rn×l and W ∈ Rn×m. For a vector x to be in
V ∩ W, it must be the case that there exist vectors λ ∈ Rl
, µ ∈ Rm such that
x = Vλ = Wµ .
Rearranging this, we get
0 =
V − W
”
λ
µ
#
.
For any vector ”
λ
µ
#
in the null space of
V − W
, the vector x = Vλ = Wµ will be in Im(V) ∩ Im(W).
Create a function sub_intersect.m. The first line of the function will be
function C = sub_intersect(A, B)
The input matrices A and B contain bases for the subspaces we wish to intersect. The matrix C that is to
be returned should contain a basis for the intersection of the subspaces. You will need to use the Matlab
command null, which returns a basis for the null space (or kernel) of a given matrix.
Output 2. Return to your lab2.m script. Using the functions you developed above, find matrices spanning the sum and intersection of the vector spaces from earlier. Print those matrices.
3.2 Linear transformations and linear equations
Matlab commands
rank(A) Find the rank of A
A\b Find a particular solution of Ax = b
We now turn our attention to the basic operations required for linear transformations and systems
of linear equations. The first basic operation we consider is how to perform a change of basis. Let
P =
x
1
· · · x
n
∈ Rn×n be a matrix whose columns constitute a basis for Rn
. We want to express a vector
x ∈ Rn
in terms of this basis. That is, we want coefficients z1, . . . , zn such that x = z1x
1 + · · · + znx
n
, or
x = Pz.
The vector z ∈ Rn
represents the coordinates of x in the basis {x
1
, . . . , x
n}, and we see from the above
that z is given by z = P
−1x.
Output 3. In your lab2.m script, consider the basis {x
1
, x
2} of R2 given by
x
1 =
”
1
1
#
, x
2 =
”
2
1
#
.
Find the coordinates z of the vector x =
2 1⊤
in this basis and print it. Then, numerically verify that
the identity x = z1x
1 + z2x
2 holds.
Next we turn to matrix representations of linear transformations. For this, review the Procedure
to find the matrix representation in Section 2.4 of the textbook. Consider a matrix A ∈ Rm×n and
5
the associated linear transformation A : Rn → Rm given by A(x) = Ax. Suppose we have matrices
P =
x
1
· · · x
n
∈ Rn×n and Q =
y
1
· · · y
m
∈ Rm×m whose columns constitute bases for Rn and
Rm, respectively. To find the matrix representation Aˆ of the transformation A in the given bases, we
perform the following steps
1. For each j ∈ {1, . . . , n}, we compute A(x
j
), where x
j
is the j-th column of P.
2. We express A(x
j
) in the coordinates of the basis {y
1
, . . . , y
m}, where y
i
is the i-th column of Q.
3. The m-vector of coordinates found in point 2 is the j-th column of Aˆ.
Output 4. • Using the conceptual procedure above, find the expression of the matrix Aˆ in terms of
the matrices A, P, Q. Note that the textbook solves this problem in the special case when P = Q
and m = n.
• Using your expression for Aˆ, in your lab2.m script find the matrix representation of the transformation R4 → R5 defined as
A(x) =
1 2 0 −1
0 1 −1 −2
1 0 3 4
0 −1 2 3
0 0 2 2
x (1)
in the bases
R
4 = span
1
0
0
0
,
1
1
0
0
,
1
0
1
0
,
1
0
0
1
, R
5 = span
1
1
0
0
0
,
1
−1
0
0
0
,
0
0
1
1
0
,
0
0
1
−1
0
,
1
0
0
0
1
.
For this part, you might want to review Sections 1.7-1.9 of the textbook. In order to determine the
key properties of injectivity and surjectivity of a transformation, we first recall the rank-nullity theorem
from linear algebra. Recall that this theorem states that for a linear transformation A : Rn → Rm,
rank(A) + nullity(A) = n ,
where the rank of A is the dimension of the image of A, and the nullity is the dimension of the kernel
of A. Since the rank is easy to check computationally, we can use it for simple tests of the injectivity and
surjectivity of a transformation.
Output 5. In your lab2.m script, consider the matrix A given in (1).
6
• Print the rank of this matrix.
• Using the rank-nullity theorem, print the dimension of Ker(A).
• Display simple messages saying whether the matrix is injective, surjective, both, or neither.
N.B.: please do this programmatically. That is, we should be able to replace the above matrix in
your code with any other matrix, and the script would still work appropriately.
Finally, consider a system of linear equations Ax = b, A ∈ Rm×n
, x ∈ Rn
, b ∈ Rm. We want to
determine if this system has any solutions, and if so if it has just one or many. To do this, recall that the
rank of a matrix is by definition the dimension of the column space. In order for the equation Ax = b
to have a solution, it is necessary and sufficient that the vector b lie in the column space of A. So by
constructing the augmented matrix
h
A bi
and comparing its rank to the rank of A, we can determine if Ax = b has any solutions1
.
The question of whether Ax = b has only one or infinitely many solutions, provided it has at least
one, comes down to determining the kernel of A, since if x is a particular solution of Ax = b and y is in
the kernel of A, then A(x + y) = Ax + Ay = b. So the general solution of Ax = b is x + y, where x is a
particular solution and y is an arbitrary vector in Ker(A)
2
. We have already seen the Matlab command
to find a basis for the kernel of a matrix. If a particular solution of Ax = b exists, in Matlab you can find
it using the backslash command, A\b.
Output 6. Return to your lab2.m script, and consider again the matrix A in (1).
• For each of the following vectors b, determine if Ax = b has no solutions, one solution, or infinitely
many solutions.
b =
1
0
0
0
0
and b =
1
1
−2
−2
−2
• Print each of the preceding vectors, and indicate the number of solutions. If a solution exists, print
it. If infinite solutions exist, provide two different solutions to the equation.
1See Theorem 1.34 in the textbook.
2See Theorem 1.36 in the textbook
7
4 A-invariance and Representation Theorem
In this course we make heavy use of the notion of A-invariant subspace. We use it, for instance, to compute Kalman decompositions. In this section you will learn (a) how to check whether a given subspace
is A-invariant, and (b) how to numerically compute the coordinate transformation in the Representation
Theorem for A-invariant subspace3
, as well as the block upper triangular matrix representation of A in
new coordinates.
Testing for A-invariance. Suppose we have a subspace V ⊂ Rn
represented by a matrix V ∈ Rn×k
whose columns are a basis for V, and that we are given a matrix A ∈ Rn×n
. To check if V is A-invariant,
we need to check whether Avi ∈ V for each of the columns v
i of V, or equivalently whether Im(AV) ⊂
Im(V). Numerically, this condition can be easily checked by testing whether rank
AV V
= rank V
(see Section 2.9 of the textbook).
Finding the matrix P of the Representation Theorem. We already have that the columns of V form
a basis for the subspace V. In order to find the coordinate transformation T(x) = P
−1x of the Representation Theorem, it only remains to find a basis for a subspace W that is an independent complement of
V, i.e., such that V ⊕ W = Rn
. Equivalently, we need to find a matrix W ∈ Rn×(n−k)
such that the matrix
P :=
V W
is invertible. Once this is done, the linear transformation T(x) = P
−1x will guarantee that
the matrix representation of A in new coordinates, Aˆ = P
−1AP, is block upper triangular.
To find the matrix W ∈ Rn×(n−k)
, we rely on the following observation. Let the columns of W be a
basis for Ker(V
⊤), so that V
⊤W = 0. Then, the columns of V and W are mutually orthogonal vectors in
Rn
, and one can show (we will not do it) that Im(W) is an independent complement of Im(V) (in fact,
the orthogonal complement). In conclusion, we numerically find W by finding a basis for Ker(V
⊤), which
we already know how to do in Matlab.
Output 7. Return to your script lab2.m, and consider the matrix
A =
1 2 2
1 −3 −5
−1 2 4
,
and the subspace
V = span
0
1
−1
,
1
−2
1
.
• Numerically show that V is A-invariant.
• Find the matrix P of the Representation Theorem, print it out, and verify that P
−1AP is block
upper triangular.
3Theorem 2.46 in the textbook.
8
5 Controllability and Kalman decomposition
Matlab commands
ctrb(A,B) Construct the controllability matrix of the pair (A, B)
For a control system (A, B), the Kalman decomposition for controllability is a straightforward application
of the Representation Theorem in the special case when the A-invariant subspace of interest is Im(Qc),
the image of the controllability matrix
Qc
:=
h
B AB · · · A
n−1B
i
.
The tools we have developed throughout this lab will now allow you to determine if the following
system is controllable, and if not, what its Kalman decomposition for controllability is.
Return to your lab2.m script, and consider the control system
x˙ =
5 −11 5
0 −6 0
−5 5 −1
x +
1 −2
0 0
1 2
u .
• Use the Matlab command ctrb to construct the controllability matrix for this system and check
if (A, B) is controllable. You will find that the system isn’t controllable. We have seen in class4
that the subspace V = Im(Qc) is A-invariant, and that Im(B) ⊂ Im(Qc). The application of the
Representation Theorem gives the Kalman decomposition for controllability, and for this you will
use the concepts presented in Section 4.
• Find a full-rank matrix V whose columns form a basis for Im(Qc).
• Using the concepts of Section 4, find the matrix P of the Representation Theorem.
• Using the coordinate transformation z = T(x) = P
−1x, find the Kalman decomposition (Aˆ, Bˆ) =
(P
−1AP, P
−1B), and verify that it has the form predicted in Section 5.7 of the textbook.
Output 8. • Print (Aˆ, Bˆ), the Kalman decomposition for controllability of the above system.
• Print the controllable and uncontrollable subsystems.
