This week there are no supplied files; you will be writing everything from scratch. You should make a directory for this lab and work in it:
mkdir ~/cs371/lab7 cd ~/cs371/lab7and when you are ready to submit files make sure they are in your lab7 directory.
Note that this lab is tougher than the previous ones have been; that's why you're getting two weeks to do it. I will be available to answer questions about this lab during my office hours and will take e-mail questions and random visits throughout the week.
IMPORTANT: It is OK for you to work in small groups on this lab, but if you do so, you must let me know who you're working with in order to get credit. Put your partners' names in a comment header in your code or e-mail me the names.
[L,U] = lu(A)or you can see the permutation matrix (a representation of the pivoting) with
[L,U,P] = lu(A)The permutation matrix is another way to represent the information in the pivot vector. One way to interpret the pivot vector is to consider the ith element l(i) and say "the row that is now the ith element of the matrix A used to be the l(i)th element." So for example if you have the pivot vector
l = [ 4 1 3 2 5 ]then you know that the row that used to be the fourth row of A is now the first row, the row that used to be the first row is now the second, and the row that used to be the second row is now the fourth.
The permutation matrix P has the same information as the pivot vector l. P(i,j) == 1 exactly when l(i) = j. In other words, the above pivot vector l corresponds to the permutation matrix
| 0 0 0 1 0 | | 1 0 0 0 0 | P = | 0 0 1 0 0 | | 0 1 0 0 0 | | 0 0 0 0 1 |The advantage of using the permutation matrix is that you can swap the rows of A by premultiplying the permutation matrix -- in other words the matrix product PA is A with the rows swapped around. You should try a few examples just to be sure you understand what lu is doing. Try it also with a bad matrix (for example, one with two repeated rows or with a row of zeros).
The backslash function calculates the answer to a system of equations. When you use your LU factorization to compute the answer to a system of equations, you can use backslash to check your answer. If you have a matrix A and a column vector b and you want to solve the system of equations Ax = b, you can do
x = A\bwhere that slash is the 'backslash' and not the regular forward slash.
Your subroutine should differ from the pseudocode in that it should not overwrite the input matrix. Your first line should be
function [U,l] = gausselim(A);where U is the output matrix and l is the pivot vector. (It is traditional, but not required, MATLAB usage for matrices to be capitalized and vectors to be in lower case.)
Hint: One way to avoid overwriting A is to create the matrix U right away -- maybe make it a matrix full of zeros the same size as A by doing
U = zeros(size(A));and then build U up as the algorithm progresses.
If your code runs and gets a correct answer for square input matrices, you will get most of the credit for this problem. For full credit your program must reject non-square input matrices and contain at least one fewer for loop than the pseudocode on pp236-237.
function [L,U,l] = lufact(A);If your code runs and gets a correct answer for square input matrices, you will get most of the credit for this problem. For full credit your program must reject non-square input matrices.
function x = solve(L,U,l,b)where L, U, and l are the output of your lufact program and b is the right hand side of the original linear system Ax = b. As before, you will get most of the credit for a program which runs correctly on square input matrices, but for full credit you will need to reject non-square inputs. The output vector x should be permuted so that it corresponds to the original order of the rows in A -- that's the purpose of passing the pivot vector l into the function.
To test your functions, start with simple examples and work up to more complicated ones. For example, try doing the Gaussian elimination of the following matrices to start with:
| 1 0 0 | | 1 0 0 | | 0 1 0 | | 0 1 1 | | 0 0 1 | | 0 1 2 |and work up to more complicated examples. Checking your answers for lufact and solve will be easier because you will be able to check them against MATLAB's lu and backslash functions, respectively. I can suggest some good test matrices in lab if you wish.
~cs371/submit lab7to submit your gausselim.m, lufact.m, and solve.m programs. Please be sure to use those exact names or else submit won't find them.
These programs will be due by 11:59pm on Thursday, March 14. Please budget your time appropriately, especially if you don't feel like you got halfway there by the end of lab on the 7th.
There are no questions for the lab period on March 7th -- I will write up some questions for the 14th and link them to this page a few days before then.