Lecture 10 Cholesky Decomposition Banded Matrices Considerations for High Performance AMath 352 Mon Apr 19 1 10 Example 1 2 1 2 1 0 1 L1 2 1 0 3 2 1 0 1 0 2 L2 2 1 0 0 3 2 1 0 4 3 L1 1 2 1 0 L2 0 0 0 1 0 0 0 1 0 0 2 3 1 A LU L 1 2 U 0 1 0 0 2 3 1 2 1 0 0 3 2 1 0 0 4 3 0 2 1 0 1 0 0 1 Matrices for which Pivoting is not Required If A is symmetric and positive de nite SPD then pivoting is not required Instead of factoring such a matrix as LU one may factor it as LLT 2 10 Matrices for which Pivoting is not Required If A is symmetric and positive de nite SPD then pivoting is not required Instead of factoring such a matrix as LU one may factor it as LLT Example 2 1 1 0 1 0 2 1 2 L1 0 2 1 0 3 2 1 2 0 1 L2 0 2 1 0 3 2 1 4 3 0 0 L1 1 0 0 1 2 1 0 0 1 0 L2 0 0 1 0 0 1 0 2 3 1 A LU L 1 1 2 0 0 1 0 0 2 3 1 U 0 2 1 0 3 2 1 4 3 0 0 2 10 Can write U 0 3 2 2 0 0 0 0 0 4 3 0 0 1 0 1 1 2 0 1 2 3 DLT Thus A LDLT LD 1 2 LD 1 2 T L LT L 2 2 2 0 0 cid 112 3 2 0 0 2 3 cid 112 3 2 2 3 This is called the Cholesky factorization Cholesky Factorization A LU L 1 1 2 0 0 1 0 0 2 3 1 U 0 2 1 0 3 2 1 4 3 0 0 3 10 Cholesky Factorization A LU L 1 1 2 0 0 1 0 0 2 3 1 U 0 2 1 0 3 2 1 4 3 0 0 Can write U 0 2 0 3 2 0 0 0 0 4 3 1 1 2 0 0 1 0 0 2 3 1 DLT Thus A LDLT LD 1 2 LD 1 2 T L LT 2 2 2 0 0 cid 112 3 2 0 0 2 3 cid 112 3 2 2 3 L This is called the Cholesky factorization 3 10 Banded Matrices Many applications give rise to banded matrices a1 b1 b1 0 0 0 a2 b2 b2 0 0 0 a3 b3 a4 b3 3 ops a1 b1 0 0 0 0 a2 b2 b2 0 0 0 a3 b3 a4 b3 3 ops a1 b1 0 0 0 0 a2 b2 0 0 0 0 a3 b3 a4 b3 3 ops a1 b1 0 0 0 0 a2 b2 0 0 0 0 a3 b3 a4 0 Total 3n ops 4 10 Suppose half bandwidth is m and we don t have to pivot Use row one to eliminate in col one of rows 2 through m 1 the entry in col one of rows m 2 through n is already 0 a1 m 1 a2 m 1 a2 m 2 am 1 2m 2 a11 a21 am 1 1 an n m Banded Matrices Cont Half bandwidth is m if aij 0 for i j m Half bandwidth of a tridiagonal matrix is 1 5 10 Banded Matrices Cont Half bandwidth is m if aij 0 for i j m Half bandwidth of a tridiagonal matrix is 1 Suppose half bandwidth is m and we don t have to pivot Use row one to eliminate in col one of rows 2 through m 1 the entry in col one of rows m 2 through n is already 0 a1 m 1 a2 m 1 a2 m 2 am 1 2m 2 a11 a21 am 1 1 an n m 5 10 Banded Matrices Cont At rst stage work on m rows modify m 1 entries Work at rst stage reduced from 2n n 1 to 2m m 1 Total work n 1 cid 88 j 1 2 min m n j min m 1 n j 1 2m2n 6 10 Banded Matrices Cont At rst stage work on m rows modify m 1 entries Work at rst stage reduced from 2n n 1 to 2m m 1 Total work n 1 cid 88 j 1 2 min m n j min m 1 n j 1 2m2n 6 10 a1 b1 0 c1 a2 b2 0 0 0 0 0 c2 0 0 a3 b3 c3 0 a4 b4 c4 a5 0 0 0 P1 c1 a2 b2 a1 b1 0 0 0 c2 0 0 a3 b3 c3 0 a4 b4 c4 a5 0 0 0 L1 0 0 0 0 0 0 0 0 a4 b4 c4 a5 c1 a2 b2 0 0 0 0 c2 0 0 c3 0 a3 b3 Bandwidth at most doubles Banded Matrices with Pivoting Suppose we have a banded matrix but we do need to pivot 7 10 Bandwidth at most doubles Banded Matrices with Pivoting Suppose we have a banded matrix but we do need to pivot a1 b1 c1 0 0 0 0 a2 b2 c2 0 0 0 0 a3 b3 c3 0 0 0 0 a4 b4 a5 c4 P1 c1 a1 b1 c2 0 0 0 0 0 0 a2 b2 0 0 a3 b3 c3 0 0 0 0 a4 b4 a5 c4 L1 0 c1 a2 b2 0 0 a3 b3 c2 0 c3 0 0 0 0 0 0 0 0 a4 b4 a5 c4 7 10 Banded Matrices with Pivoting Suppose we have a banded matrix but we do need to pivot a1 b1 c1 0 0 0 0 a2 b2 c2 0 0 0 0 a3 b3 c3 0 0 0 0 a4 b4 a5 c4 P1 c1 a1 b1 c2 0 0 0 0 0 0 a2 b2 0 0 a3 b3 c3 0 0 0 0 a4 b4 a5 c4 L1 0 c1 a2 b2 0 0 a3 b3 c2 0 c3 0 0 0 0 0 0 0 0 a4 b4 a5 c4 Bandwidth at most doubles 7 10 Consider GE code To eliminate the i j entry at stage j Fetch the pivot row j and row i i j operate on row i store the result 3 n j 1 memory accesses 2 n j 1 ops Implementation Considerations for High Performance cid 73 Memory refs may be costly cid 73 Usually there is a memory hierarchy larger storage such as main memory or disk are slower smaller storage such as local cache are faster cid 73 When an arithmetic op is performed operands are fetched from storage loaded into registers where the arithmetic is performed and stored back to memory cid 73 The more work that can be done with a word …