MAT 343 Arizona State University LU Factorization Lab Questions You nuse MATLAB program to solve for this lab, see the file for more information MAT 343 La

MAT 343 Arizona State University LU Factorization Lab Questions You nuse MATLAB program to solve for this lab, see the file for more information MAT 343 Laboratory 3
The LU Factorization
In this laboratory session we will learn how to
1. Find the LU factorization of a matrix using elementary matrices
2. Use the MATLAB command lu to find the LU factorization
3. Solve linear systems using the LU factorization
Instructions: For your lab write-up follow the instructions of LAB 1.
Elementary Matrices and Row Reduction
An elementary matrix is any matrix that results from applying a single elementary row operation to the identity matrix. In this section we will discover some properties of elementary
matrices and their connection to row reduction of matrices.
For example, the following three elementary matrices were obtained from the 4 × 4 identity
matrix:
?
? ?
? ?
?
1 0 0 0
?5 0 0 0
0 0 1 0
?0 1 0 0? ? 0 1 0 0? ??3 1 0 0?
?
? ?
? ?
?
(1)
?1 0 0 0? ? 0 0 1 0? ? 0 0 1 0?
0 0 0 1
0 0 0 1
0 0 0 1
(The row operations were, respectively: interchange rows 1 and 3; multiply row 1 by ?5; multiply row 1 by -3 and add to row 2).
We can construct each of the three types of elementary matrices using basic MATLAB
commands.
Elementary Matrices of Type I:
To swap rows i and j of the n × n identity matrix, set
E1 = eye(n);
E1([i,j],:) = E1([j,i],:)
Elementary Matrices of Type II:
To multiply row i of the n × n identity matrix by the scalar c, set
E2 = eye(n);
E2(i,i) = c
Elementary Matrices of Type III:
To replace row j of the n × n identity matrix by the sum of row j plus c times row i, set
E3 = eye(n);
E3(j,i) = c
THIS CONTENT IS PROTECTED AND MAY NOT BE SHARED, UPLOADED, SOLD OR DISTRIBUTED
2019 v1 Copyright@ School of Mathematical and Statistical Sciences, Arizona State University
1
Example 1
We can generate the three matrices listed in (1) using the following MATLAB commands.
>> E1 = eye (4);
>> E1 ([1 ,3] ,:)= E1 ([3 ,1] ,:)
E1 =
0
0
1
0
0
1
0
0
1
0
0
0
0
0
0
1
>> E2 = eye (4);
>> E2 (1 ,1)= -5
% generate the 4 x4 identity matrix
% swap rows 1 and 3
% generate the 4 x4 identity matrix
% change the entry in position (1 ,1) to -5
% equivalent to multiply row 1 by -5
E2 =
-5
0
0
0
0
1
0
0
>> E3 = eye (4); %
>> E3 (2 ,1)= -3 %
%
E3 =
1
0
-3
1
0
0
0
0
0
0
1
0
0
0
0
1
generate the 4 x4 identity matrix
change the entry in position (2 ,1) to -3
equivalent to multiply row 1 by -3 and add to row 2
0
0
1
0
0
0
0
1
Important Property: Suppose that E is an n × n elementary matrix obtained from the
identity by performing one of the elementary row operations. If A is an n × r matrix, then the
matrix EA is obtained by performing the same row operation on A.
EXERCISE 1
If you haven’t already done so, enter the commands in Example 1 to generate the elementary
matrices E1, E2 and E3 (you can suppress these matrices in your lab report).
Generate a random 4 × 3 matrix A (with integer entries), by typing
A = floor(10*rand(4,3))
Compute the products E1*A, E2*A and E3*A and compare A with each of the products. Describe
the effect of each multiplication on the matrix A. Be specific. For each of the three products,
describe specifically how the rows of A have changed.
What general pattern do you see? That is, what effect does multiplying a matrix A on the left
by an elementary matrix have on the matrix A?
THIS CONTENT IS PROTECTED AND MAY NOT BE SHARED, UPLOADED, SOLD OR DISTRIBUTED
2019 v1 Copyright@ School of Mathematical and Statistical Sciences, Arizona State University
2
The LU factorization using elementary matrices
EXERCISE 2
Enter the following matrix in MATLAB.
?
?
?2 1 ?4
5 13 ?
A=? 8
6 ?9 17
(a) Determine elementary matrices E1 , E2 , E3 of Type III such that E3 E2 E1 A = U with U
an upper triangular matrix. The matrix E1 should turn the element in position (2,1) into
a 0. Enter this matrix in MATLAB as E1 using commands similar to the ones in Example
1. The matrix E2 should turn the element in position (3,1) into a zero. Enter this matrix
in MATLAB as E2. Note that to zero out the entries in column 1, you need to add or
subtract a multiple of row 1. Once you have found the matrices E1 and E2 , compute
the product E2 E1 A in MATLAB. Use format rat so that the entries will be given as
fractions. Based on the result, determine the matrix E3 that turns the element in position
(3,2) into a zero. Enter this matrix as E3 in MATLAB and compute U=E3*E2*E1*A.
(b) Compute the product L = E1?1 E2?1 E3?1 . The matrix L is lower triangular with ones on
the diagonal. Enter format short and verify that A = LU by computing A ? LU in
MATLAB.
NOTE: From part (a) we have E3 E2 E1 A = U and therefore A = (E3 E2 E1 )?1 U = (E1?1 E2?1 E3?1 )U .
Since the elementary matrices and their inverses are lower triangular and have 1’s along the
diagonal, as is always true for elementary matrices of the third type, it follows that L is also
lower triangular with 1’s along the diagonal.
Permutation matrices
A permutation matrix is a square matrix that has exactly one 1 in every row and column and
0’s elsewhere. In this section we will look at properties of permutation matrices.
One way to construct permutation matrices is to permute the rows (or columns) of the
identity matrix. For example, we can construct
?
?
0 0 0 0 1
?0 0 1 0 0?
?
?
?
(2)
E=?
?0 0 0 1 0?
?1 0 0 0 0?
0 1 0 0 0
by using the MATLAB commands
p = [5 ,3 ,4 ,1 ,2]; % permutation vector that defines
% the new order of the rows
E = eye ( length ( p ) ) ;
% define E as the identity matrix
E = E (p ,:) % permute the rows of E according to the
% permutation vector p
The second command creates the 5×5 identity matrix, and the third command uses the vector p
to permute its rows, so row p(1)= 5 becomes row 1, row p(2)= 3 becomes row 2, row p(3)= 4
becomes row 3 and so on (compare these row permutations with the vector p defined above).
THIS CONTENT IS PROTECTED AND MAY NOT BE SHARED, UPLOADED, SOLD OR DISTRIBUTED
2019 v1 Copyright@ School of Mathematical and Statistical Sciences, Arizona State University
3
EXERCISE 3
If you haven’t already done so, enter the commands in the example above to generate the
permutation matrix E defined in (2) (you can suppress this matrix).
Generate a 5 × 5 matrix A with integer entries using the command
A = floor(10*rand(5))
(a) Compute the product EA and compare the answer with the matrix A. How are the
two matrices related? Describe the effect on A of left multiplication by the permutation
matrix E. Be specific!
Compute the product AE and compare the answer with the matrix A. How are the two
matrices related? Describe the effect on A of right multiplication by the permutation
matrix E. Be specific!
(b) Compute E ?1 and E T (recall that E T is computed in MATLAB with the command E’),
and observe that they are also permutation matrices. What else do you observe about
E ?1 and E T ?
The LU factorization in MATLAB
Not all matrices have an LU factorization because not all matrices can be reduced to upper
triangular form using only row operation of Type III. However, if we allow the rows to be
interchanged, then all matrices have an LU factorization (but we have to keep track of the rows
we permute).
Example 2
The matrix
?
?
2 4 6
A=? 1 2 5 ?
?1 5 3
does not have an LU factorization. Zeroing out the entries below the pivot in the first column
yields
?
?
2 4 6
? 0 0 2 ?
0 7 6
and it is impossible to put this matrix in triangular form without permuting the rows. However,
the matrix
?
?
2 4 6
P A = ? ?1 5 3 ?
1 2 5
(obtained by interchanging row 2 and row 3 of
the permutation matrix
?
1
?
P = 0
0
A) does have an LU factorization. Here P is
?
0 0
0 1?
1 0
In general, for any matrix A, it is always possible to find the LU factorization of P A, where
P is an appropriate permutation matrix.
The MATLAB command
[L,U,P] = lu(A)
THIS CONTENT IS PROTECTED AND MAY NOT BE SHARED, UPLOADED, SOLD OR DISTRIBUTED
2019 v1 Copyright@ School of Mathematical and Statistical Sciences, Arizona State University
4
returns a lower triangular matrix L (with 1’s on the diagonal), an upper triangular matrix U ,
and a permutation matrix P so that P*A = L*U.
For this example we have:
>> A=[2,4,6;1,2,5;-1,5,3]
A =
2
4
1
2
-1
5
>> [L,U,P]=lu(A)
L =
1
0
-1/2
1
1/2
0
U =
2
4
0
7
0
0
P =
1
0
0
0
0
1
6
5
3
0
0
1
6
6
2
0
1
0
and we can easily verify that P A = LU .
Remark: Note that, even if the LU factorization of A does exist, MATLAB will most likely
use a permutation matrix anyway. This is because MATLAB permutes the rows so that the
pivot is always the largest entry in the column. This technique, called partial pivoting, reduces
round-off errors.
Solving Systems using the LU factorization
Given A = LU , we can solve Ax = b in two steps:
(a) First solve Ly = b
(b) then solve U x = y.
Both systems are triangular and therefore can be easily solved using backward and forward
substitution.
We can carry out these two steps using the “” command in MATLAB:
y = Lb
x = Uy
However, in general, as observed before, MATLAB will return the LU factorization of P A
rather than A, thus, instead of solving the system Ax = b we will solve the equivalent system
P Ax = P b
(3)
To solve (3) using steps (a) and (b) defined above, we only need to replace b with P b:
y = L(P*b)
x = Uy
THIS CONTENT IS PROTECTED AND MAY NOT BE SHARED, UPLOADED, SOLD OR DISTRIBUTED
2019 v1 Copyright@ School of Mathematical and Statistical Sciences, Arizona State University
5
EXERCISE 4
Enter the matrix A and the vector b
?
?9
? 3
A=?
? ?6
?1
in MATLAB:
?
6 6 ?8
3 4 ?9 ?
?,
5 3 ?6 ?
2 6
4
?
?
?85
? 23 ?
?
b=?
? ?49 ?
?69
The exact solution to the system Ax = b is the vector x = (7, ?2, ?7, ?4)T .
(a) Enter [L,U,P] = lu(A) to find the LU factorization of the matrix P A. Verify that
P A = LU by computing P A ? LU .
(b) Use the LU factorization you found in part (a) to solve the system Ax = b. Call the
computed solution x lu (don’t forget that you will need P*b).
(c) Enter the vector x and compare your solution x lu from part (b) with the exact solution
x by computing norm(x lu – x)
NOTE: the norm function gives the magnitude of the vector,
p that is, for a vector a =
(a1 , a2 , . . . , an )T , the norm of a is defined as: norm(a) = a21 + a22 + . . . + a2n . Thus
the command norm(x lu – x) measures how close the computed solution x lu is to the
exact solution x. You should expect a very small number.
EXERCISE 5
In this exercise we will compare in MATLAB the speed of two methods for solving the system
Ax = b when A is an invertible square matrix: computing the RREF and computing the LU
factorization. Note that, although the number of operations for obtaining the LU factorization
is the same as the gaussian elimination, the LU factorization has the advantage that once the
matrix A is decomposed, the substitution step can be carried out efficiently for different values
of b. Thus the LU factorization is certainly preferable when solving the system Ax = b with
different values of b for the same A.
We will use the MATLAB tic and toc command to measure the computation times.
Enter:
A = rand(500);
x = ones(500,1);
b=A*x;
Important: Be sure to use semicolon after each command so that matrices and vectors are
not displayed. Do not print or include these large matrices and vectors in your lab write-up.
(a) Solve Ax = b using the reduced row echelon form and store the solution in x rref:
tic; R = rref([A, b]);
x_rref = R(:,end); toc
(make sure that all the commands are on the same line).
(b) Solve Ax = b using the LU factorization as you did in EXERCISE 4(a)(b), and calculate
the elapsed time using the tic toc commands. Store the solution in x lu.
Which method is faster?
NOTE: Make sure you use semicolon; the only output should be the elapsed time. Also,
don’t forget that, when using tic toc all the commands should be in one line, so you
should find L,U,P, solve for y, and solve for x lu all on one line.
(c) Compare the solutions from parts (a) and (b) with the exact solution x by computing
norm(x rref – x) and norm(x lu – x). How accurate are the solutions from parts (a)
and (b)?
THIS CONTENT IS PROTECTED AND MAY NOT BE SHARED, UPLOADED, SOLD OR DISTRIBUTED
2019 v1 Copyright@ School of Mathematical and Statistical Sciences, Arizona State University
6

Purchase answer to see full
attachment