Square CoCalc Logo
FeaturesSoftwarePricingInfoPoliciesShareSupport Try Sign InSign Up
Project: admcycles
Views: 77
Visibility: Unlisted (only visible to those who know the link)
Image: ubuntu2004
| Embed | Download | Raw
Kernel: SageMath 9.1

Lecture 2 : Linear Algebra

References:

Summary:
We start with a reminder on lists in Python, and functions for creating and manipulating them. Then we discuss matrices in SageMath: how to construct them, how to solve linear systems of equations, computer kernels and eigenvalues etc. Finally, we discuss linear spaces and how to intersect or add them.

Python warmup : lists

Recall that in Python, a list L is an ordered collection of objects. They can be created by hand, using square brackets [,]:

In [ ]:
L = [3, sin, pi, -12/5]; L

You can access element number i in the list (with count starting from 0) using L[i]:

In [ ]:
L[0]
In [ ]:
L[1]

This can also be used to change an element of the list:

In [ ]:
L[1] = cos L

Lists have many useful functions:

  • you can get the number of elements in L using len(L)

  • given a second list M you can compute the concatenation of L and M by L+M

  • you can check if a is contained in L by the expression a in L

Exercise

  • Create a list M with 3 elements of your choice.

In [ ]:
  • Define N to be the list obtained by concatenating L (from above) and M.

In [ ]:
  • Check that the length of N is now equal to 7.

In [ ]:

Solution (uncomment to see)

A list L also has many useful methods, like L.pop(i) which returns L[i] and removes this entry from the list. As we saw last week, you can see all these methods by typing L.+Tab to access the Tab completion.

Exercise

Consider the following list:

In [ ]:
L = [77, 23, 49, 66, 39, 26, 48, 23, 35, 70, 64, 9, 21, 15, 78, 21, 2, 56, 42, 53]
  • What is the number i such that L[i] equals 9?
    Hint: Use Tab completion and the documentation to find the right function.

In [ ]:
  • Remove the corresponding entry from the list.

In [ ]:

Solution (uncomment to see)

Finally, a nice way to create bigger lists is using list comprehension involving the for constructor. We will learn more about this and for-loops later, but for now we just mention that

[f(i) for i in range(n)]

creates the list [f(0), f(1), ..., f(n-1)], and similarly

[f(i) for i in range(a,b)]

creates the list [f(a), f(1), ..., f(b-1)].

In [ ]:
[3*i+2 for i in range(7)]

Exercise

You get an email from a local school, asking for help since they misplaced some of their math education materials. In particular, they request whether you can send them:

  • a list of the square numbers 1^2, ..., 20^2

In [ ]:
  • a multiplication table for the numbers from 1 to 10

In [ ]:

Solution (uncomment to see)

Matrices

At the center of linear algebra are matrices, which describe linear maps between finite-dimensional vector spaces.

Constructing matrices

From last time, we recall that we can create matrices using the function matrix. We saw that it can handle at least two different types of input:

  • matrix([v1, v2, ..., vm]) takes a list of vectors, forming the rows of the desired matrix,

  • matrix(m,n,f) takes the numbers m,n of rows and columns of the matrix MM, and a function f sending pairs (i,j)(i,j) to the desired entry Mi,jM_{i,j} of the matrix

Exercise

The Vandermonde matrix associated to a vector a=(a1,,an)a=(a_1, \ldots, a_n) of numbers is the matrix V(a1,,an)=(1a1a1n11a2a2n11anann1). V(a_1, \ldots, a_n)= \begin{pmatrix} 1 & a_1 & \ldots & a_1^{n-1}\\ 1 & a_2 & \ldots & a_2^{n-1}\\ \vdots & & & \vdots\\ 1 & a_n & \ldots & a_n^{n-1} \end{pmatrix}\,. Below is the beginning of an implementation of the function taking a list a of numbers and returning the corresponding Vandermonde matrix. Complete the program and test it in some examples.

In [ ]:
def vandermonde_matrix(a): n = len(a) return ???
In [ ]:

Solution (uncomment to see)

One further possibility to construct a matrix by hand is to call matrix with a dictionary. We are going to discuss these later in more detail, but below is a hopefully self-explaining example:

In [ ]:
matrix(2,3,{(0,1):17, (1,2):-3})

For many standard types of matrices, there are also pre-defined functions:

In [ ]:
matrix.zero(2,3)
In [ ]:
matrix.diagonal([5,-1,0,-1])
In [ ]:
A = matrix([[3,9],[6,10]]) block_matrix([ [A, -A], [A.transpose(), 100*A] ], subdivide=False)

On the other hand, if someone hands you a matrix M, you can find out its number of rows and columns using M.nrows() and M.ncols():

In [ ]:
M = matrix([[10,9,8],[7,6,5]]) M
In [ ]:
M.nrows()
In [ ]:
M.ncols()

You can access the entry Mi,jM_{i,j} using M[i,j]:

In [ ]:
M[0,1]

Sometimes it can be useful to explicitly specify that a matrix should have entries in a certain ring. This ring can then be given as the first argument to the matrix function.

Some common rings:

  • ZZ the integers

  • QQ the rational numbers

  • RR the real numbers (floating point operations!)

  • CC the complex numbers (floating point operations!)

  • GF(q) the finite field Fq\mathbb{F}_q with qq elements

Example

The matrix M=(1135) M = \begin{pmatrix} 1 & 1 \\ 3 & 5 \end{pmatrix} has rank 22 seen as a matrix of rational numbers.

In [ ]:
M = matrix(ZZ,[[1,1],[3,5]])
In [ ]:
M.rank()

But over the finite field F2=Z/2Z\mathbb{F}_2 = \mathbb{Z}/2\mathbb{Z} we have M=(1135)=(1111) M = \begin{pmatrix} 1 & 1 \\ 3 & 5 \end{pmatrix} = \begin{pmatrix} 1 & 1 \\ 1 & 1 \end{pmatrix} and so the matrix only has rank 11.

In [ ]:
M = matrix(GF(2),[[1,1],[3,5]])
In [ ]:
M.rank()

Exercise

Find a formula for the determinant det(V(a1,,an))\det(V(a_1, \ldots, a_n)) of the Vandermonde matrix. Bonus points if you can get the sign correct ...
Hint: Another useful ring is a polynomial ring, which you can construct as follows:

In [ ]:
R.<a,b,c> = PolynomialRing(QQ,3) # here 3 is the number of variables R
In [ ]:
f = a^2 + 2*a*b + b^2 factor(f)
In [ ]:

Solution (uncomment to see)

Linear algebra with matrices

Now that we can input matrices, we can start doing more interesting computations.

Matrix arithmetic and vector products

We already saw that we can perform matrix arithmetic with +,*,^:

In [ ]:
M = matrix([[2,6],[0,-7]]) N = matrix([[1,0],[2,-9]]) [M+N, M*N, M^2]

In particular, when M is invertible, we can compute its inverse by M^(-1) (or M.inverse()):

In [ ]:
M^(-1)
In [ ]:
M.inverse()

Note that you can also compute the expression MxM^x for a formal variable xx. Below we use that in SageMath x is by default a variable in the symbolic ring, which we will see later.

In [ ]:
M^x

We can also create vectors and multiply them with matrices.

In [ ]:
v = vector([2,-5]) M*v

Row echelon form

Another interesting computation for a matrix MM is its row echelon form (German: Zeilenstufenform), which is the simplified form of MM obtained by row-operations.

In [ ]:
M = matrix(QQ,[[2,4,-1,7],[-3,-6,3,-2],[1,2,5,-8]]) M
In [ ]:
M.echelon_form()

Note that if we work over the integers ZZ, the algorithm only does integer operations on the rows (it is not allowed to divide a row by a number), so the output will be slightly different:

In [ ]:
M = matrix(ZZ,[[2,4,-1,7],[-3,-6,3,-2],[1,2,5,-8]]) M.echelon_form()

We can also get the matrix UU describing the row-operations to bring MM into echelon form EE, so that UM=EU \cdot M = E:

In [ ]:
E, U = M.echelon_form(transformation=True)
In [ ]:
U
In [ ]:
U*M

Exercise

Does there exist a subset of the vectors v1=(253),v2=(4106),v3=(107),v4=(154),v5=(51013),v6=(5211),v7=(91530)R3 v_1 = \begin{pmatrix}-2\\ 5\\ 3 \end{pmatrix}, v_2 = \begin{pmatrix}4\\ -10\\ -6 \end{pmatrix}, v_3 = \begin{pmatrix}1 \\0 \\ -7 \end{pmatrix}, v_4 = \begin{pmatrix} -1\\ 5\\ -4 \end{pmatrix}, v_5 = \begin{pmatrix} -5\\ 10\\ 13 \end{pmatrix}, v_6 = \begin{pmatrix} 5\\ -2\\ -11 \end{pmatrix}, v_7 = \begin{pmatrix} 9\\ -15\\ -30 \end{pmatrix} \in \mathbb{R}^3 which forms a basis of R3\mathbb{R}^3? If so, specify such a subset.
Hint: To save you time, here is the list of the vectors above:

In [ ]:
V = [(-2, 5, 3), (4, -10, -6), (1, 0, -7), (-1, 5, -4), (-5, 10, 13), (5, -2, -11), (9, -15, -30)]
In [ ]:
In [ ]:

Solution (uncomment to see)

Solving linear equations

Of course, we can also solve linear systems of equations using SageMath: given a matrix M and a vector b the solution of the linear system Mx=b M \cdot x = b is computed using M.solve_right(b).

In [ ]:
M = matrix([[1,2],[3,4]]) b = vector([1,2]) [M, b]
In [ ]:
v=M.solve_right(b) v
In [ ]:
M*v

As the name suggests, there is also the methods M.solve_left(b), which solves the equation xTM=bT, x^T \cdot M = b^T, where xT,bTx^T, b^T are the transposed (or row) vectors of x,bx,b.

In [ ]:
w = M.solve_left(b) w

Exercise

Find out what M.solve_right(b) does when the system has

  • more than one solution, or

  • no solution,

either by reading the documentation or by experimenting.

Solution (uncomment to see)

Further matrix functions

As you can imagine, there are many more things that SageMath can compute about matrices. Instead of going through all of them, let's have an exercise where you can practice finding out about some of these functions yourself. In this regard, the following advice has proven useful in the past:

SageMath-Tip No. 4
Just google it.

Exercise Consider the following matrix:

In [ ]:
M = matrix(QQ,[(29, -3, -3, -27, -24),(30, -1, 6, -12, -54),(-30, 6, -4, 6, 60),(18, -3, 0, -10, -24),(3, 0, -3, -9, 8)]) M
  • What are the eigenvalues of M?

In [ ]:
  • Compute a matrix SS such that SMS1S M S^{-1} is a diagonal matrix.

In [ ]:
  • Compute a permutation matrix PP, a lower-triangular matrix LL and an upper triangular matrix UU such that A=PLUA = PLU.
    (Hint: This is often called a LU-decomposition)

In [ ]:
  • Is MM a normal matrix, i.e. does it satisfy MM=MMM \cdot M^\top = M^\top \cdot M?

In [ ]:

Solution (uncomment to see)

Linear spaces

In SageMath we can also directly create vector spaces and perform computations with subspaces (like direct sums, intersections etc.).

In [ ]:
V = VectorSpace(QQ,3) V

Then we can create subspaces of V:

In [ ]:
W1 = V.subspace([(1,2,3), (1,1,1)]) W2 = V.subspace([(1,-1,1), (0,0,1)]) W1

Here degree 3 means that the ambient space V has is of rank 3, and dimension 2 means that W1 itself is of dimension 2.

Now we can check if a vector is contained in one of those subspaces, or compute their intersection or sum in V.

In [ ]:
vector(QQ,[1,2,3]) in V
In [ ]:
W1.intersection(W2)
In [ ]:
W1+W2

An alternative way to create subspaces is using the function span:

In [ ]:
W1alt = span(QQ,[(1,2,3), (1,1,1)]) W1 == W1alt

Exercise

Consider the following vectors in Q4\mathbb{Q}^4:

In [ ]:
v = vector(QQ,[4,-3,2,6]) w1 = vector(QQ,[1,-7,-3,0]) w2 = vector(QQ,[9,-4,-4,1]) w3 = vector(QQ,[-15,21,8,-7])

Use the methods for linear spaces we saw above, to answer the following questions:

  • Are w1,w2,w3w_1, w_2, w_3 linearly independent?

In [ ]:
  • Is vv contained in the subspace of Q4\mathbb{Q}^4 generated by w1,w2,w3w_1, w_2, w_3?

In [ ]:

Solution (uncomment to see)

Of course one could also find the answer to these questions using matrices, and in the background this is what SageMath does. Given a space like W1, you can get access to the matrix whose rows form a basis of W1 using W1.basis_matrix():

In [ ]:
W1.basis_matrix()

Exercise

Let's apply the concepts from above to write a nice function ourselves, which computes orthogonal projections. Recall that given a sub-vector space VRnV \subseteq \mathbb{R}^n and a point pRnp \in \mathbb{R}^n, the orthogonal projection of pp onto VV is the unique point p0Vp_0 \in V such that the vector pp0p-p_0 is orthogonal to VV.

In [ ]:
f = (lambda u,v: u+v, lambda u,v: u, lambda u,v: v) plane = parametric_plot3d(f, (-2,2), (-2,2)) arrow = arrow3d((5/3, 4/3, 1/3), (1,2,1), 1, color='red') plane + arrow
  • Complete the definition of the function orthogonal_projection below, which takes as input the vector space V and the vector p and returns the orthogonal projection p_0 as above.

In [ ]:
def orthogonal_projection(V, p): return ???

Hint: If v1,,vmVv_1, \ldots, v_m \in V is a basis of VV, then there exist numbers x1,,xmx_1, \ldots, x_m such that p0=x1v1++xmvmp_0 = x_1 v_1 + \ldots + x_m v_m. Reformulate the condition that pp0p - p_0 is orthogonal to VV as a linear system Ax=bA \cdot x = b.

  • Test your function in some non-trivial example.

In [ ]:

Solution (uncomment to see)

Assignments

Exercise

For which values of xQx \in \mathbb{Q} is the matrix M(x)=(23x5x6110) M(x) = \begin{pmatrix} -2 & 3 & x \\ -5 & x & 6 \\ 1 & -1 & 0 \end{pmatrix} invertible?

Solution (uncomment to see)

Exercise

Write a function jordan_form(blocks) taking as input a list blocks of the form [[e1,m1], ..., [en,mn]] of pairs [ei,mi] of eigenvalues and multiplicities, and producing the corresponding Jordan block diagonal matrix. Here is an example how it would work:

jordan_form([[3,2],[4,1],[-1,3]]) > [ 3 1 0 0 0 0] [ 0 3 0 0 0 0] [ 0 0 4 0 0 0] [ 0 0 0 -1 1 0] [ 0 0 0 0 -1 1] [ 0 0 0 0 0 -1]

For an easier variant of the exercise: try to create the particular 6×66 \times 6-matrix above with as little code as possible.

Solution (uncomment to see)

Exercise

Consider the following matrix:

In [74]:
rows = [(1812, 2239, -3, 1, 17711), (-158013, -195140, 2067, -689, -1543587), (-189, -234, 153, -50, -1851), (-567, -702, 465, -152, -5553), (19791, 24441, -261, 87, 193332)] M = matrix(rows); M
Out[74]:
[ 1812 2239 -3 1 17711] [ -158013 -195140 2067 -689 -1543587] [ -189 -234 153 -50 -1851] [ -567 -702 465 -152 -5553] [ 19791 24441 -261 87 193332]

Compute its characteristic polynomial and minimal polynomial. Why are they different?
Bonus: See the solution for an explanation how the hell I came up with this matrix.

Solution (uncomment to see)

Exercise

Inside R3\mathbb{R}^3 with coordinates x1,x2,x3x_1, x_2, x_3 consider the subspace V={(x1,x2,x3):2x1x2+3x3=0}. V = \{(x_1, x_2, x_3) : 2 x_1 - x_2 + 3 x_3 = 0\}\,.

  • Compute a basis B\mathcal{B} of VV.

  • For the linear map f:R3R3f:\mathbb{R}^3 \to \mathbb{R}^3 induced by the matrix A=(350268117) A = \begin{pmatrix} 3 & 5 & 0 \\ -2 & 6 & 8 \\ -1 & 1 & 7 \end{pmatrix} compute the matrix AVA_V for the restriction fV:VR3f|_V : V \to \mathbb{R}^3 with respect to your chosen basis B\mathcal{B} of VV and the standard basis of R3\mathbb{R}^3.

Solution (uncomment to see)