T4: Matrix operations#
Note
Click the and open this notebook in Colab to enable interactivity.
Note
To save your progress, make a copy of this notebook in Colab File > Save a copy in Drive
and you’ll find it in My Drive > Colab Notebooks
.
From this point in the course onward, the tutorials are mostly theoretical, but we’ll demonstrate a few numerical solutions here to verify the results.
Determinants#
Determinants of a matrix can be computed using the np.linalg.det(A)
function.
We’ll demonstrate this by computing the determinant of \(A = \begin{bmatrix} 1 & 1 & 1 & 1 \\ 1 & 1 & 1 & 2 \\ 1 & 1 & 2 & 3 \\ 1 & 2 & 3 & 4 \end{bmatrix}\).
import numpy as np
A = np.array([[1, 1, 1, 1], [1, 1, 1, 2], [1, 1, 2, 3], [1, 2, 3, 4]])
np.linalg.det(A)
np.float64(-1.0)
Matrix inversion#
If \(A\) is a square matrix with \(\det(A) \neq 0\), its inverse \(A^{-1}\) exists. A system of linear equation \(A \vec{x} = \vec{b}\) can be solved by multiplying \(A^{-1}\) on the left. \(A^{-1}\) can be computed using Gaussian elimination.
We’ll demonstrate this by computing the inverse of \(A = \begin{bmatrix} 1 & 1 & 1 & 1 \\ 1 & 1 & 1 & 2 \\ 1 & 1 & 2 & 3 \\ 1 & 2 & 3 & 4 \end{bmatrix}\), and then solve for \(A \vec{x} = \begin{bmatrix} 1 \\ 1 \\ 2 \\ 2 \end{bmatrix}\). We will rely on several functions here:
np.linalg.inv(A)
- Calculate the inverse of matrixA
.np.linalg.solve(A, b)
- Solve the equation \(A \vec{x} = \vec{b}\) for \(\vec{x}\).sympy.Matrix.rref()
- Calculate the reduced row echelon form.
A = np.array([[1, 1, 1, 1], [1, 1, 1, 2], [1, 1, 2, 3], [1, 2, 3, 4]])
b = np.array([[1], [1], [2], [2]])
# use rref
from sympy import Matrix
aug = Matrix(np.concatenate([A, np.eye(4)], axis=1))
rref = aug.rref()[0]
Ainv = rref[:, 4:] # stays as sympy Matrix, can convert with sympy.matrix2numpy(Ainv)
display(Ainv @ A) # verify
# Multiply by b to get the solution
x = Ainv @ b
display(x)
print()
# Directly reduce [A, b]
aug = Matrix(np.concatenate([A, b], axis=1))
rref = aug.rref()[0]
x = rref[:, -1]
display(x)
print()
# Direct solve with NumPy inverse
Ainv = np.linalg.inv(A)
x = Ainv @ b
print(x)
print()
# Direct solve with NumPy
x = np.linalg.solve(A, b)
print(x)
print(A @ x - b) # verify
[[ 1.]
[-1.]
[ 1.]
[ 0.]]
[[ 1.]
[-1.]
[ 1.]
[ 0.]]
[[0.]
[0.]
[0.]
[0.]]