CSC338. Homework 4
Due Date: Wednesday Feburary 5, 9pm
Please see the guidelines at https://www.cs.toronto.edu/~lczhang/338/homework.html
What to Hand In
Please hand in 2 files:
• Python File containing all your code, named hw4.py.
• PDF file named hw4_written.pdf containing your solutions to the written parts of the assignment. Your solution can be hand-written, but must be legible. Graders may deduct marks for illegible or poorly presented solutions.
If you are using Jupyter Notebook to complete the work, your notebook can be exported as a .py file (File -> Download As -> Python). Your code will be auto-graded using Python 3.6, so please make sure that your code runs. There will be a 20% penalty if you need a remark due to small issues that renders your code untestable.
Make sure to remove or comment out all matplotlib or other expensive code before submitting your homework!
Submit the assignment on MarkUs by 9pm on the due date. See the syllabus for the course policy regarding late assignments. All assignments must be done individually.
In [ ]:
import math
import numpy as np
Question 1
For this question, we will again start from code from tutorial 3.
In [ ]:
# Code from tutorial 3
def backward_substitution(A, b):
“””Return a vector x with np.matmul(A, x) == b, where
* A is an nxn numpy matrix that is upper-triangular and non-singular
* b is an nx1 numpy vector
“””
n = A.shape[0]
x = np.zeros_like(b, dtype=np.float)
for i in range(n-1, -1, -1):
s = 0
for j in range(n-1, i, -1):
s += A[i,j] * x[j]
x[i] = (b[i] – s) / A[i,i]
return x
def eliminate(A, b, k):
“””Eliminate the k-th row of A, in the system np.matmul(A, x) == b,
so that A[i, k] = 0 for i < k. The elimination is done in place."""
n = A.shape[0]
for i in range(k + 1, n):
m = A[i, k] / A[k, k]
for j in range(k, n):
A[i, j] = A[i, j] - m * A[k, j]
b[i] = b[i] - m * b[k]
def gauss_elimination(A, b):
"""Return a vector x with np.matmul(A, x) == b using
the Gauss Elimination algorithm, without partial pivoting."""
for k in range(A.shape[0] - 1):
eliminate(A, b, k)
x = backward_substitution(A, b)
return x
Part (a) -- 1 pt
Solve the system $Ax = b$ for the values of $A$ and $b$ below using the function gauss_elimination from last time. Save the solution you obtain in the variable soln_nopivot.
In [ ]:
e = pow(2, -100)
A = np.array([[e, 1],
[1, 1]])
b = np.array([1 + e, 2])
soln_nopivot = None
Part (b) -- 2 pt
Write a helper function partial_pivot that performs partial pivoting on A at column k, so that the function gauss_elimination_partial_pivot performs Gauss Elimination with Partial Pivoting.
In [ ]:
def partial_pivot(A, b, k):
"""Perform partial pivoting for column k. That is, swap row k
with row j > k so that the new element at A[k,k] is the largest
amongst all other values in column k below the diagonal.
This function should modify A and b in place.
“””
# TODO
def gauss_elimination_partial_pivot(A, b):
“””Return a vector x with np.matmul(A, x) == b using
the Gauss Elimination algorithm, with partial pivoting.”””
for k in range(A.shape[0] – 1):
partial_pivot(A, b, k)
eliminate(A, b, k)
x = backward_substitution(A, b)
return x
Part (c) — 1 pt
Solve the system $Ax = b$ for the values of $A$ and $b$ below using gauss_elimination_partial_pivot. Save the solution you obtain in the variable soln_pivot.
In [ ]:
e = pow(2, -100)
A = np.array([[e, 1],
[1, 1]])
b = np.array([1 + e, 2])
soln_pivot = None
Part (d) — 2 pt
Do your answers in parts (a) and (d) match? If not, which is the correct answer? Include your explanation in the PDF write-up.
Question 2
Part (a) — 3 pt
Consider the following matrices M1, M2, and M3. Compute each of their $L_1$, $L_2$, and $L_\infty$ norms. Save the results in the variables below.
For the $L_2$ norm you may find the function np.linalg.norm helpful. You can compute the $L_1$ and $L_\infty$ norms either by hand or write a function.
In [ ]:
M1 = np.array([[3., 0.],
[-4., 2.]])
M2 = np.array([[2., -2., 0.3],
[0.5, 1., 0.9],
[-4., -2., 5]])
M3 = np.array([[0.2, -0.2],
[1.0, 0.2]])
# fill in these answers
M1_l_1 = 0
M1_l_2 = 0
M1_l_infty = 0
M2_l_1 = 0
M2_l_2 = 0
M2_l_infty = 0
M3_l_1 = 0
M3_l_2 = 0
M3_l_infty = 0
Part (b) — 4 pt
Show that an induced matrix norm $|| \cdot ||$ has the property $$||A + B|| \leq ||A|| + ||B||$$
Part (c) — 4 pt
Is it true that for a vector, $||v||_\infty \le ||v||_1$? What about for a matrix: is it true that for a matrix, $||M||_\infty \le ||M||_1$? Include your solution and justification in your pdf writeup.
Question 3
Part (a) [3 pt]
Write a function matrix_condition_number that computes the condition number of a $2 \times 2$ matrix. Use the $$L_1$$ matrix norm.
In [ ]:
def matrix_condition_number(M):
“””
Returns the condition number of the 2×2 matrix M.
Use the $L_1$ matrix norm.
Precondition: M.shape == [2, 2]
M is non-singular
>>> matrix_condition_number(np.array([[1., 0.], [0., 1.]]))
1
“””
Part (b) [2 pt]
Classify each of the following matrices A1, A2, A3 and A4, as well-conditioned or ill-conditioned.
You may do this question either by hand, or by using the function above.
Save the classifications in a Python array called conditioning. Each item of the array should be either the string “well” or the string “ill”.
In [ ]:
A1 = np.array([[3, 0],
[0, pow(3, -30)]])
A2 = np.array([[pow(3, 20), 0],
[0, pow(5, 50)]])
A3 = np.array([[pow(4, -40), 0],
[0, pow(2, -80)]])
A4 = np.array([[pow(5, -17), pow(5, -16)],
[pow(5, -18), pow(5, -17)]])
# whether A1, A2, A3, A4 is “well” or “ill” conditioned
conditioning = [None, None, None, None]
Part (c) [2 pt]
It should be immediate obvious that the matrix
$$\begin{bmatrix}100 & 2 \\ 201 & 4\end{bmatrix}$$
is ill-conditioned. Explain why. Include your answer in your pdf write-up.
Part (d) [2 pt]
Suppose that $A$ and $B$ are two $n \times n$ matrices, and both are well-conditioned. Is $A(B^{-1})$ also well-conditioned? Why or why not? Include your answer and justifications in your pdf write-up. Be specific.
Part (e) [4 pt]
Describe an efficient algorithm to compute $d^T B^T A^{-1} B d$
Where:
• $A$ is an invertible $n \times n$ matrix,
• $B$ is an $n \times n$ matrix, and
• $d$ is an $n \times 1$ vectors
Be clear and specific. Include your strategy in your pdf write-up.