DOC PREVIEW
UMD ASTR 415 - Numerical Linear Algebra

This preview shows page 1 out of 4 pages.

Save
View full document
View full document
Premium Document
Do you want full access? Go Premium and unlock all 4 pages.
Access to all documents
Download any document
Ad free experience
Premium Document
Do you want full access? Go Premium and unlock all 4 pages.
Access to all documents
Download any document
Ad free experience

Unformatted text preview:

Class 7. Numerical Linear Algebra, Part 2LU Decomposition• Suppose we can write A as a product of two matrices: A = LU, where L is lower triangularand U is upper triangular:L =× 0 0× × 0× × ×U =× × ×0 × ×0 0 ×• Then Ax = (LU) x = L (Ux) = b, i.e., must solve,(1) Ly = b; (2) Ux = y.• Can reuseL and U for subsequent calculations.• Why is this better?– Solving triangular matrices is easy: just use forward substitution for (1), back-substitution for (2).• Problem is, how to decompose A into L and U?– Expand matrix multiplication LU to get n2equations for n2+ n unknowns (ele-ments of L and U plus n extras because diagonal elements counted twice).– Get an extra n equations by choosing Lii= 1 (i = 1 , n).– Then use Crout’s algorithm for finding solution to these n2+ n equations “triv-ially” (NRiC §2.3).LU decomposition in NRiC• The routines ludcmp() and lubksb() perform LU decomposition and backsubstituion,respectively.• Can easily compute A−1(solve for the identity matrix column by column) and det(A)(find the product of the diagonal elements of the LU decomposed matrix)—see NRiC§2.3.• Warning: for large matrices, computing det(A) can overflow or underflow the com-puter’s floating-point dynamic range (there are ways around this).Iterative Improvement• For large sets of linear equations Ax = b, roundoff error may become a problem.• We want to know x but we only have x+δx, which is an exact solution to A(x + δx) =b + δb.1• Subtract the first equation from the second, and use the second to eliminate δb:Aδx = A(x + δx) − b.• The RHS is known, hence can solve for δx. Subtract this from the wrong solution toget an improved solution (make sure to use doubles!). See mprove() in NRiC.Tridiagonal Matrices• Many systems can be written as (or reduced to):aixi−1+ bixi+ cixi+1= dii = 1, ni.e., a tridiagonal matrix:b1c100sa2b2c2a3b3c3.........an−1bn−1cn−100s anbnx1x2x3...xn−1xn=d1d2d3...dn−1dn.Here a1and cnare associated with “boundary conditions” (i.e., x0and xn+1).• LU decomposition and backsubstitution is very efficient for tri-di systems: O(n) oper-ations as opposed to O( n3) in general case.Sparse Matrices• Operations on many sparse systems in general can be optimized, e.g.,tridiagonal;band diagonal with bandwidth M;block diagonal;banded.• See NRiC §2.7 for various systems and techniques.Iterative methods• For very large systems, direct solution methods (e.g., LU decomposition) are slow a ndRE prone.• Often iterative methods much more efficient:1. Guess a trial solution x0.2. Compute a correction x1= x0+ δx.3. Iterate pro cedure until convergence, i.e., |δx| < ∆.• E.g., congugate gradient method for sparse systems (NRiC §2.7).2Singular Value Decomposition• Can diagnose or (nearly) solve singular or near-singular systems.• Used for solving linear least-squares problems.• Theorem: any m × n matrix A (with m rows and n columns) can be written:A = UWVT,where U (m × n) and V (n × n) are orthogonal1and W (n × n) is a diagonal matrix.• Proof: buy a good linear algebra textbook...• The n diagonal values wiof W are zero or positive and are called the “singular values.”• The NRiC routine svdcmp() returns U, V, and W given A. You have to trust it (ortest it yourself!).– Uses Householder reduction, QR diagonalization, etc.• If A is square, then we know2A−1= V [diag(1/wi)] UT.– This is fine so long as no wiis too small (or zero). Otherwise, the presence ofsmall or zero witell you how singular your system is...Definitions• Condition numbercond(A) = (max wi)/(min wi).– If cond(A) = ∞, A is singular.– If cond(A) very large (∼ e−1m), A is ill-conditioned.• Consider Ax = b. If A is singular, there is some subspace of x (the nullspace) suchthat Ax = 0.• The nullity of A is the dimension of the nullspace (the number of linearly independentvectors x that can be found in it).• The subspace of b such that Ax = b is the range of A.• The rank of A is the dimension of the range.1U has orthonormal columns while V, being square, has both orthonormal rows and columns.2Since U and V are square and orthogonal, their inverses are equal to their transposes, and since W isdiagonal, its inverse is a diagona l matrix whose elements are 1/wi.3The homogeneous equation• SVD constructs orthonormal bases for the nullspace and range of a matrix.• Columns of U with corresponding non-zero wiare an orthonormal basis for the range.• Columns of V with corresponding zero wiare an orthonormal basis for the nullspace.• Hence immediately have solution for Ax = 0, i.e., the columns of V with correspondingzero wi.Residuals• If b (6= 0) lies in the range of A, then the singular equations do in fact have a solution.• Even if b is outside the range of A, can get solution which minimizes residualr =|Ax − b|.– Trick: replace 1/wiby 0 if wi= 0 and computex = V [diag(1/wi)] (UTb).• Similarly, can set 1/wi= 0 if wivery small.Approximation of matrices• Can write A = UWVTasAij= ΣNk=1wkUikVjk.• If most of the singular values wkare small, then A is well-approximated by only a fewterms in the sum (strategy: sort wk’s in descending order).• For large memory savings, just store the columns of U and V corresponding to non-negligible wk’s.• Useful technique for digital image


View Full Document

UMD ASTR 415 - Numerical Linear Algebra

Download Numerical Linear Algebra
Our administrator received your request to download this document. We will send you the file to your email shortly.
Loading Unlocking...
Login

Join to view Numerical Linear Algebra and access 3M+ class-specific study document.

or
We will never post anything without your permission.
Don't have an account?
Sign Up

Join to view Numerical Linear Algebra 2 2 and access 3M+ class-specific study document.

or

By creating an account you agree to our Privacy Policy and Terms Of Use

Already a member?