313x Filetype PDF File size 0.14 MB Source: www.math.usm.edu
Jim Lambers
MAT461/561
Spring Semester 2009-10
Lecture 15 Notes
These notes correspond to Section 6.6 in the text.
Special Types of Matrices
The solution of a system of linear equations Ax = b can be obtained using Gaussian elimination
with pivoting in conjunction with back substitution for any nonsingular matrix A. However, there
are several classes of matrices for which modifications of this process are more appropriate. We
will now examine some of these classes.
Diagonally Dominant Matrices
Amatrix A is strictly diagonally dominant if each diagonal element aii is larger in magnitude than
the sum of the magnitudes of all other elements in the ith row. That is,
n
X ∣a ∣<∣a ∣, i=1,2,...,n.
ij ii
j=1,j∕=i
Wealso say that A is strictly row diagonally dominant.
It can be shown that if A is strictly diagonally dominant, then
∙ A is nonsingular.
∙ Gaussian elimination can be applied to A without performing any row interchanges.
This is useful for ensuring that there is a unique solution of the system of linear equations that
must be solved to obtain the coefficients of a cubic spline, as such a system is strictly diagonally
dominant.
Symmetric Positive Definite Matrices
Areal, n×n matrix A is positive definite if, for any nonzero vector x,
xTAx>0.
Manyapplications involve symmetric positive definite matrices, which are positive definite and also
T
symmetric; that is, A = A . A positive definite matrix is the generalization to n×n matrices of a
positive number.
If A is symmetric positive definite, then it has the following properties:
1
∙ A is nonsingular.
∙ All of the diagonal elements of A are positive.
∙ The largest element of the matrix lies on the diagonal.
∙ All of the eigenvalues of A are positive.
∙ det(A) > 0.
It is worth noting that the first two properties hold even if A is not symmetric.
In general it is not easy to determine whether a given n×n symmetric matrix A is also positive
definite. One approach is to check the matrices
⎡ a a ⋅ ⋅ ⋅ a ⎤
11 12 1k
⎢ a a ⋅ ⋅ ⋅ a ⎥
A =⎢ 21 22 2k ⎥, k = 1,2,...,n.
k ⎢ . . . ⎥
⎣ . . . ⎦
. . .
a a ⋅ ⋅ ⋅ a
k1 k2 kk
These matrices are called the leading principal submatrices of A, or the principal minors. It can be
shown that A is positive definite if and only if det(A ) > 0 for k = 1,2,...,n.
k
One desirable property of symmetric positive definite matrices is that Gaussian elimination
can be performed on them without pivoting, and all pivot elements are positive. Furthermore,
Gaussianelimination applied to such matrices is robust with respect to the accumulation of roundoff
error. However, Gaussian elimination is not the most practical approach to solving systems of
linear equations involving symmetric positive definite matrices, because it is not the most efficient
approach in terms of the number of floating-point operations that are required.
Instead, it is preferable to compute the Cholesky factorization of A,
A=LLT,
where L is a lower triangular matrix with positive diagonal entries. Because A is factored into two
matrices that are the transpose of one another, the process of computing the Cholesky factorization
requires about half as many operations as the LU decomposition.
The algorithm for computing the Cholesky factorization can be derived by matching entries of
the product LLT with those of A. This yields the following relation between the entries of L and
A,
k
a =Xℓijℓ , i,k=1,2,...,n, i≥k.
ik kj
j=1
From this relation, we obtain the following algorithm.
2
for j = 1,2,...,n do
√
ℓ = a
jj jj
for i = j + 1,j + 2,...,n do
ℓ =a /ℓ
ij ij jj
for k = j,j +1,...,i do
a =a −ℓijℓ
ik ik kj
end
end
end
If A is not symmetric positive definite, then the algorithm will break down, because it will attempt
to compute ℓjj, for some j, by taking the square root of a negative number, or divide by a zero ℓjj.
In fact, due to the expense involved in computing determinants, the Cholesky factorization is also
an efficient method for checking whether a symmetric matrix is also positive definite.
Once the Cholesky factor L of A is computed, a system Ax = b can be solved by first solving
Ly=bbyforward substitution, and then solving LTx = y by back substitution.
T
An alternative factorization of a symmetric positive definite matrix is the LDL factorization
T
A=LDL ,
where L is unit lower triangular, and D is a diagonal matrix whose diagonal entries are positive.
Thecomputation of the factors L and D is very similar to the algorithm for the LU decomposition,
with the following exceptions:
∙ During elimination, only elements in the lower triangle of A, that is, entries a where i ≥ j,
ij
need to be updated, due to the symmetry of A.
∙ Before eliminating elements in column j, the jth diagonal element, d , is assigned the value
jj
(j) (j)
of ajj ; that is, djj is set equal to the (j,j) element of A .
T
The LDL factorization is also known as the “square-root-free Cholesky factorization”, since it
computes factors that are similar in structure to the Cholesky factors, but without computing any
T 1/2
square roots. In fact, if A = GG is the Cholesky factorization of A, then G = LD .
Banded Matrices
An n×n matrix A is said to have upper bandwidth p if a =0whenever j ≥ i+p. Similarly, A
ij
has lower bandwidth q if a =0whenever i ≥ j +q. A matrix that has upper bandwidth p and
ij
lower bandwidth q is said to have bandwidth w = p+q +1.
Any n×n matrix A has a bandwidth w ≤ 2n−1. If w < 2n−1, then A is said to be banded.
However, cases in which the bandwidth is O(1), such as when A is a tridiagonal matrix for which
p = q = 1, are of particular interest because for such matrices, Gaussian elimination, forward
substitution and back substitution are much more efficient. This is because
3
∙ If A has lower bandwidth q, and A = LU is the LU decomposition of A (without pivoting),
then L has lower bandwidth q, because at most q elements per column need to be eliminated.
∙ If A has upper bandwidth p, and A = LU is the LU decomposition of A (without pivoting),
then U has upper bandwidth p, because at most p elements per row need to be updated by
row operations.
It follows that if A has O(1) bandwidth, then Gaussian elimination, forward substitution and back
substitution all require only O(n) operations each, provided no pivoting is required.
When a matrix A is banded with bandwidth w, it is wasteful to store it in the traditional
2-dimensional array. Instead, it is much more efficient to store the elements of A in w vectors of
length at most n. Then, the algorithms for Gaussian elimination, forward substitution and back
substitution can be modified appropriately to work with components of these vectors instead of
entries of a matrix. For example, to perform Gaussian elimination on a tridigonal matrix, we can
proceed as in the following algorithm. We assume that the main diagonal of A is stored in the
vector a, the subdiagonal (entries a ) is stored in the vector l, and the superdiagonal (entries
j+1,j
a ) is stored in the vector u.
j,j+1
for j = 1,2,...,n−1 do
lj = lj/aj
a =a −l u
j+1 j+1 j j
end
Then, we can use these updated vectors to solve the system Ax = b using forward and back
subtitution as follows:
y =b
1 1
for i = 2,3,...,n do
y =b −l y
end i i i−1 i−1
x =y /a
n n n
for i = n−1,n−2,...,1 do
x =(y −ux )/a
i i i i+1 i
end
Note that after Gaussian elimination, the components of the vector l are the subdiagonal entries of
Linthe LU decomposition of A, and the components of the vector u are the superdiagonal entries
of U.
4
no reviews yet
Please Login to review.