Advanced Matrix Algorithms for Efficient Calculations in CS542G

Notes • Assignment 1 is out (due October 5) • Matrix storage: usually column-major cs542g-term1-2006

Block Approach to LU • Rather than get bogged down in details of GE (hard to see forest for trees) • Partition the equation A=LU • Gives natural formulas for algorithms • Extends to block algorithms cs542g-term1-2006

Cholesky Factorization • If A is symmetric positive definite, can cut work in half: A=LLT • L is lower triangular • If A is symmetric but indefinite, possibly still have the Modified Cholesky factorization: A=LDLT • L is unit lower triangular • D is diagonal cs542g-term1-2006

Pivoting • LU and Modified Cholesky can fail • Example: if A11=0 • Go back to Gaussian Elimination ideas: reorder the equations (rows) to get a nonzero entry • In fact, nearly zero entries still a problem • Possibly due to cancellation error => few significant digits • Dividing through will taint rest of calculation • Pivoting strategy: reorder to get the biggest entry on the diagonal • Partial pivoting: just reorder rows • Complete pivoting: reorder rows and columns (expensive) cs542g-term1-2006

Pivoting in LU • Can express it as a factorization:A=PLU • P is a permutation matrix: just the identity with its rows (or columns) permuted • Store the permutation, not P! cs542g-term1-2006

Symmetric Pivoting • Problem: partial (or complete) pivoting destroys symmetry • How can we factor a symmetric indefinite matrix reliably but twice as fast as unsymmetric matrices? • One idea: symmetric pivoting PAPT=LDLT • Swap the rows the same as the columns • But let D have 2x2 as well as 1x1 blocks on the diagonal • Partial pivoting: Bunch-Kaufman (LAPACK) • Complete pivoting: Bunch-Parlett (safer) cs542g-term1-2006

Reconsidering RBF • RBF interpolation has advantages: • Mesh-free • Optimal in some sense • Exponential convergence (each point extra data point improves fit everywhere) • Defined everywhere • But some disadvantages: • It’s a global calculation(even with compactly supported functions) • Big dense matrix to form and solve(though later we’ll revisit that… cs542g-term1-2006

Gibbs • Globally smooth calculation also makes for overshoot/undershoot(Gibbs phenomena) around discontinuities • Can’t easily control effect cs542g-term1-2006

Noise • If data contains noise (errors), RBF strictly interpolates them • If the errors aren’t spatially correlated, lots of discontinuities: RBF interpolant becomes wiggly cs542g-term1-2006

Linear Least Squares • Idea: instead of interpolating data + noise, approximate • Pick our approximation from a space of functions we expect (e.g. not wiggly -- maybe low degree polynomials) to filter out the noise • Standard way of defining it: cs542g-term1-2006

Rewriting • Write it in matrix-vector form: cs542g-term1-2006

Normal Equations • First attempt at finding minimum:set the gradient equal to zero(called “the normal equations”) cs542g-term1-2006

Good Normal Equations • ATA is a square kk matrix(k probably much smaller than n) • Symmetric positive (semi-)definite cs542g-term1-2006

Bad Normal Equations • What if k=n?At least for 2-norm condition number, k(ATA)=k(A)2 • Accuracy could be a problem… • In general, can we avoid squaring the errors? cs542g-term1-2006

Advanced Matrix Algorithms for Efficient Calculations in CS542G