Assignment 1

Due on September 26 at 11:59pm electronically via email.

How to submit the assignment

Use a word-processor like Microsoft Word or LibreOffice (free), or \( \LaTeX \) (also free) to create a single PDF file with all of your answers. If you want to include code, do not put them in this file. This file should include your graphs, equations, tables, and explanations. The name of this file should have the following format: LASTNAME_STUDENTNUMBER_A1.pdf. You can hand-write your answers, but you will have to scan your papers, or take a picture of them. In any case, send me one PDF file. I will not accept paper submissions.

If you want to include multiple code files (.m or .py files), put them all in one folder, zip the folder, and send me only one zip file. The name of this file should look like this: LASTNAME_STUDENTNUMBER_A1.zip

Attach both files to an email and send it to msamani[at]fields.utoronto.ca. Make sure the title of the email starts with MATH3090.

Question 1: Consider the number \(x=1.43\).

1. Represent \(x\) in binary and hexadecimal formats.
2. Did you have to round the last bit up or down?
3. What are the absolute and relative errors associated with rounding?

Questions from the textbook: 2.3, 2.5, 2.9, 2.10, 2.19.

Notes for Python users

Question 2.5

The fucntions in part (a) are as follows:

import numpy as np
import scipy.linalg

# I could not find an implementation of magic in python, so we have to define it ourselves.
# Code was adopted from https://stackoverflow.com/questions/47834140/numpy-equivalent-of-matlabs-magic and modified on Sep 10, 2018.
def magic(n):
	n = int(n)
	if n < 3:
		raise ValueError("Size must be at least 3")
	if n % 2 == 1:
		p = np.arange(1, n + 1)
		return n * np.mod(p[:, None] + p - (n + 3) // 2, n) + np.mod(p[:, None] + 2 * p - 2, n) + 1
	elif n % 4 == 0:
		J = np.mod(np.arange(1, n + 1), 4) // 2
		K = J[:, None] == J
		M = np.arange(1, n * n + 1, n)[:, None] + np.arange(n)
		M[K] = n * n + 1 - M[K]
	else:
		p = n // 2
		M = magic(p)
		M = np.block([[M, M + 2 * p * p], [M + 3 * p * p, M + p * p]])
		i = np.arange(p)
		k = (n - 2) // 4
		j = np.concatenate((np.arange(k), np.arange(n - k + 1, n)))
		M[np.ix_(np.concatenate((i, i + p)), j)
			] = M[np.ix_(np.concatenate((i + p, i)), j)]
		M[np.ix_([k, k + p], [0, k])] = M[np.ix_([k + p, k], [0, k])]
	return M

# scipy.linalg has a built-in command for Hilbert matrices
H = scipy.linalg.hilbert(n) # this will make a Hilbert matrix of size 3

# scipy.linalg also has a built in package for Pascal matrices
P = scipy.linalg.pascal(n)

# numpy has an equivalent eye command
E = np.eye(n, n)
		

I have translated the functions defined in section 2.7 into Python:

# Assumes L is a lower-triangular square matrix
def forward(L, b):
	N = len(L[0])
	x = np.zeros(N)
	for k in range(N):
		j = slice(k)
		x[k] = (b[k] - L[k,j].dot(x[j]) ) / L[k,k]
	return x
	
# Assumes U is an upper-triangular square matrix
def backsubs(U, b):
	N = len(U[0])
	x = np.zeros(N)
	for k in range(N-1, -1, -1):
		j = slice(k,N)
		x[k] = (b[k] - U[k,j].dot(x[j])) / U[k,k]
	return x


# LU factorization
def LUtx(AA):
	A = np.copy(AA)
	N = len(A[0])
	p = np.array(range(N)).T
	for k in range(N-1):
		# Fint the largest element below diagonal in k-th column
		m = np.argmax(np.abs(A[k:N,k])) + k
		
		# Skip elimination if column is zero
		if A[m,k] != 0:
			# Swap pivot now
			if m != k:
				A[[k,m],:] = A[[m,k],:]
				p[[k,m]] = p[[m,k]]
			
			# Compute multipliers
			i = slice(k+1,N)
			A[i,k] = A[i,k]/A[k,k]
			
			# Update the remainder of the matrix
			j = slice(k+1,N)
			A[i,j] = A[i,j] - A[i,k:k+1] * A[k:k+1,j]
	L = np.tril(A, k=-1) + np.eye(N)
	U = np.triu(A, k=0)
	return L, U, p

# It solves A.dot(x) = b and returns x
def bslashtx(A,b):
	N = len(A[0])
	if np.all(np.triu(A, k=1) == np.zeros((N,N))):
		# Lower triangular
		return forward(A, b)
	elif np.all(np.tril(A, k=-1) == np.zeros((N,N))):
		# Upper triangular
		return backsubs(A, b)
	elif np.all(A == A.T):
		# symmetric
		L = np.linalg.cholesky(C)
		y = forward(L, b)
		x = backsubs(L.T,y)	
	else:
		L,U,p = LUtx(A)
		y = forward(L, b[p])
		x = backsubs(U,y)
		return x		
		

Question 2.19 The sparse matrix libraries are in scipy.sparse