Programming Assigment #2 – Linear Algebra Solved

(a)   Please submit a PDF file that includes all answers and ALL THE CODE you used (even inquestions that do not require to hand in an *.py file - your code should be in the PDF). The code in the PDF file must be IDENTICAL to the code in the *.py files you submit (below in (b)). This PDF file should be readable so that you could get constructive remarks. Make sure you edit your code so that it looks nice, and in particular, DO NOT USE Print Screens! The code should appear in the PDF file as text. You should submit a printed copy to the checkers box and and also a digital copy as instructed below in (b).

(b)   In addition you should submit four *.py files (mylu.py; myeig.py; myGD.py; and Krank.py), and one *.npy file (’cities coor.npy’).

(c)    Compress the PDF file of (a) and the files mylu.py, myeig.py, myGD.py, Krank.py and ’cities coor.npy’ of (b), into a SINGLE *.zip (you can also use *.rar, but no other format!) named YourID.zip (for example: 201920209.zip) and upload the file through the ”matalot lehagasha” section in the moodle page of the course. Make sure that the file names, function signatures and input and output data structures are EXACTLY as the ones specified in the questions below. Files with names or signatures which are different in any way from what is instructed WILL NOT BE CHECKED.

(d)   The functions mylu, myeig, myGD and Krank must NOT print to the terminal, display plots or images, read or write to files or in general, do anything except for executing the algorithm for which they are intended to according to the instructions in the questions below.

(e)   If you call functions from imported Python libraries (e.g numpy and etc.) all import statements MUST be included in the *.py file. For example, if you used a function from numpy in the implementation of myeig, you must include the import statement inside myeig.py.

Some general (and hopefully helpful) advice:

(a)   There’s no need to check the validity of the input inside the functions, this is not a programming course. For instance, the function mylu.py below requires a symmetric matrix A, there’s no need to check if the inputed A is indeed symmetric.

(b)   You may and are encouraged to use built-in library functions of Python whenever possible,except for (obviously) implementations of Python for an algorithm which you are required to implement in this exercise. For example, you can and should use some Python implementation to sort the eigenvalues in question 2 and not write your own Sort function (but obviously cannot use Python’s LU decomposition).

(c)    For students who are having trouble with coding, take a look athttps://docs.scipy.org/doc/numpy/user/quickstart.html. It includes a nice quick tutorial for Python’s numpy library which should have most of what you need for this course.



1.   Implement LU factorization.

Your function should have the signature mylu(A,pivot), and return a Python list [L,U,P,Q] where

            A            Matrix to factorize

pivot 0 - no pivoting, 1 - partial pivoting, 2 - complete pivoting.

            L,U           L and U as defined by the LU algorithm.

P          Permutation matrix as defined for partial and complete pivoting. Q      Permutation matrix as defined for complete pivoting.

The algorithm should satisfy

pivot = 0
A        =LU
pivot = 1
PA =LU
pivot = 2
PAQ=LU
Important: The argument and output matrices A,L,U,P and Q must all be numpy arrays, that is, of type ndarray.

(See https://docs.scipy.org/doc/numpy/reference/generated/numpy.array.html)



2.   Implement Jacobi’s iterations to find the eigenvalues and eigenvectors of a symmetric matrix.

Your function should have the signature myeig(A), and return a Python list [D,V] where

A Matrix to diagonalize.

D Diagonal matrix with SORTED eigenvalues (largest to smallest) on the diagonal. V Matrix of eigenvectors.

The matrices V and D should satisfy V DV T = A. You may assume that the entries of A are real (i.e., ai,j ∈ R).

Important: The argument and output matrices A,D and V must all be numpy arrays, that is, of type ndarray.

(See https://docs.scipy.org/doc/numpy/reference/generated/numpy.array.html)

3.   Download the file cities.npy from the moodle site and use the load function of numpy to load the file into an ndarray variable named ’dists’. The file contains the pairwise distances between 312 cities, specifically, dists is a matrix of size 312 × 312, with

dists(i,j) = kxi − xjk

where xi,xj ∈ R2 are the coordinates in the plane of cities i and j respectively, and k·k is the Euclidean norm.

(a)                    Let X ∈ R312×2 be a matrix whose ith row is xi. Find an expression for XXT in terms of the matrix dists.

(b)                    The coordinates x1,...,x312 cannot be recovered uniquely from XXT . Explain why.

(c)                     The coordinates of the first 5 cities are

 



Find x1,...,x312, and write them as a 312 × 2 ndarray into a .npy file named ’cities coor.npy’, using the save function of numpy. In addition, write the coordinates of the first 10 cities to the PDF file.



4.   Implement Gradient Descent iterations ~xk+1 = ~xk − a · ∇f(~xk) for a general function f(x) : Rn → R. Your function should have the signature myGD(f,gradf,x0,a,tol,maxiter), and return a Python list [y, iter] where

gradf    Function handle to ∇f(x) x0  Initial vector to start the iterations from

(a)                    Step size

             tol             Tolerance

maxiter Maximal number of iterations

y           Final approximation of the vector that minimizes f iter    Number of iterations until convergence

where the tolerance tol is an upper bound on kxi+1−xik to be used as the stopping criteria. Use myGD to minimize the Rosenbrock function f(x,y) = (1 − x)2 + 100(y − x2)2 with tol = 10−6, the parameter a = 0.001 and the starting point x0 = (2,2). How many iterations did you need in order to achieve the approximation?



5.   In this question we would like to compress an image using the SVD decomposition.

•   Download the file Madril.png from the moodle site.

•   Read the image into an ndarray variable in Python named ’bwmandril’ using the imread function of matplotlib library.

•   Display the image using the imshow function of matplotlib library, with the flag cmap=’gray’ (you should see a gray-scale image of a mandrill monkey displayed on the screen now).

Answer the following questions:

(i)      Plot the singular values of the decomposition and decide according to their graph howmany singular values are needed.

(ii)     Create a function which gives a rank-k approximation of the image. Your function should have the signature Krank(A,k), and return a Python list [Uk, Sk, Vk] where

A A matrix of size m × n k The rank of the approximation

Uk m × k matrix

Sk k × k matrix of the singular values

Vk n × k matrix

The rank-k approximation of A should be Ak = Uk · Sk · VkT

Important: The argument and output matrices A,Uk,Sk and Vk must all be numpy arrays, that is, of type ndarray.

(See https://docs.scipy.org/doc/numpy/reference/generated/numpy.array.html)

(iii)   Use the Krank function to reconstruct an approximation of the image using just k singular values (the number k is for you to decide upon). Add your approximated image to the PDF file (use the matplotlib command imsave, NO Print Screens).