### Create a StudySoup account

#### Be part of our community, it's free to join!

Already have a StudySoup account? Login here

# 285 Class Note for CS 51500 at Purdue

### View Full Document

## 21

## 0

## Popular in Course

## Popular in Department

This 20 page Class Notes was uploaded by an elite notetaker on Friday February 6, 2015. The Class Notes belongs to a course at Purdue University taught by a professor in Fall. Since its upload, it has received 21 views.

## Similar to Course at Purdue

## Popular in Subject

## Reviews for 285 Class Note for CS 51500 at Purdue

### What is Karma?

#### Karma is the currency of StudySoup.

#### You can buy or earn more Karma at anytime and redeem it for class notes, study guides, flashcards, and more!

Date Created: 02/06/15

51 Banded Storage u temperature uh temperature at gridpoints Sr H O L l f 110 ap ace 5 equa 1on u0 1 grid size 1 u1 1 S H The ve point difference operator um my 7 21496711 WW 7 My WM h 7 21496711 WM 7 h 7 0 7 hz 7 hz for interior mesh points 3 natural correspondence between unknowns and equations Suppose we label the unknowns thus 10 11 12 7 8 9 4 5 6 1 2 3 123 Then7 for exarnple7 4U8U7U9U5U11 us lt gt 07 4U1 0 U2 1 U4 U1 lt gt 0 4 71 71 ul 1 71 4 71 71 uz 1 71 4 71 ug 1 41 4 41 41 w T 71 71 4 71 71 W 0 71 71 4 71 u6 0 41 4 41 41 W T 71 71 4 71 71 us 0 71 71 4 71 ug 0 41 4 41 U10 6 1 1 4 1 U11 0 1 1 4 U12 0 The upper bandwidth is 3 and the lower bandwidth is 3 The overall bandwidth 3 1 3 7 This is also a 4gtlt4 block tridiagonal matrix with 3gtlt3 blocks Other orderings of gridpoints are possible 4 8 12 9 11 12 4 12 8 11 6 12 3 7 11 6 8 10 3 11 7 4 10 5 2 6 10 3 5 7 2 10 6 8 3 9 1 5 9 1 2 4 1 9 5 1 7 2 one way red nested dissection black De nition A matrix A for which aij0 for jlti7p or jgtiq is said to have lower bandwidth p7 upper bandwidth q7 and bandwidth p 1 q p subdiagonals 1 diagonal q superdiagonals 124 i p iq Storage is row by row ail aiq ith row gt a aii aiiq 11m Pth subdiagonal diagonal q superdiagonal lf calj is stored in an n by p 1 q array A7 how does one access a For solution algorithms there are three types of matrices 1 symmetric positive de nite7 2 strictly diagonally dlominant7 3 general Symmetric positive de nite matrices 125 0 0 A 0 G store I 0 this J 1 part P P To do a Cholesky factorization A GGT7 1 change range of loops to avoid zeros nonstored zeros7 2 then adjust indexing as above Strictly diagonally dominant matrices An example is a tridiagonal matrix from a simple 2 point boundary value problem One can do Gaussian elimination without pivoting q lt gt i A NE Following is an example of symoblic77 elimination7 for which gtlt represents an element of the reduced matrix and 0 represents a multiplier XEX3X XgtltXgtlt XXX 0X 0X ZXXXXX XXX XXX XXXXXX 0XgtltgtltXX OOXXXX xxxxxx gtltgtltgtltgtltgtltgtlt 0gtltgtltgtltgtltgtlt gtltgtltgtltxgtlt XXXXXX XgtltgtltgtltgtltX General matrices Use partial pivoting Here is an example for p 27 q 1 of worst77 case interchanges 126 gtltgtlt gtltgtltXgt XXXX gtltgtltgtltgtlt Xgtltgtltgtlt XXX OXX 0gtltx XX X ooXXX X XXX Xgtltgtltgtlt OXXF OXFF gt x gtltgtltgtltx XX XXXX XXX XXXX 0gtltgtltF OXXF XXXX Here F represents ll m which is the introduction of a value not known to be zero into a nonstored zero In this example and m generaL the reduced upper triangular matrix has the same bandwidth as A storage for A 711 q 1 storage for decornposed A n2p q 1 oats n 7 1 ints worst case row 1 lt gt row p 1 row 2 lt gt row p 2 Let us do an operation count for partial pivoting 127 kth stage pp q rnults Ignoring end effects 107 1 ltlt n pp 171 rnults Operation count for case q p full matrix banded matrix no pivoting partial pivoting factorization 7133 712 2712 backsolve 712 271p 371p For Laplace7s equation on square p 7112 no pivoting partial pivoting factorization n2 2712 backsolve 27132 37132 tridiagonal 01 Ci 52 02 02 bnil anil Cnil bn an 128 forwardeliminationonmatrix fork17277n71do bk1 bk1ak ak1 ak1 bk1 6 Af forwardeliminationonf fork 17277n1d0fk1fk1bk1fk baclLsubstitution zn for kn17ni27 7 1 doxk fkick k1akg Review questions 1 What does it mean for the lower bandwidth of a matrix to be p the upper bandwidth to be q If one uses banded storage to store a matrix A of lower bandwidth p and upper band width q7 in which element of the array is aij stored Assume that array indexing starts at 1 For a symmetric positive de nite matrix A of bandwidth 2p 17 which elements do we store Which elements of its Cholesky factor G do we store What is llin How much memory would you need to solve a tridiagonal system using Gaussian elim ination with partial pivoting For a matrix of lower bandwidth p and upper bandwidth q7 what llin might be intro duced by Gaussian elimination without pivoting by Gaussian elimination with partial pivoting For an n by 71 matrix of lower bandwidth p and upper bandwidth q7 approximately how many multiplications are required for the factorization stage of Gaussian elimination without pivoting by its backsolve stage by the factorization stage of Gaussian elimination with partial pivoting by its backsolve stage Construct the algorithm for Gaussian elimination without pivoting for a tridiagonal matrix 129 Exercises 1 Complete the following algorithm for performing Gaussian elimination with partial pivoting on a tridiagonal matrix The original matrix and the reduced matrix are stored as 01 C1 01 C1 d1 52 02 02 02 02 d2 and 39 39 39 7 respectively bnil anil Cnil anil Cnil bn an an with multipliers in bg7 b37 7b and row interchange information in 717737 7rn1 fork17277n72do rk k dk 0 iflbk1lgtlakl interchange ak and bk interchange ck and 1H1 interchange dk and CH1 rk k 1 M71 71 i 1 iflbnlgtlan71l interchange owl and b interchange an and an Mei n bn bnanili an an 39 bn Cnili 2 Given below is an algorithm for the Choleski factorization of a dense symmetric positive de nite matrix A forj17277ndo 3971 ajj W Zia 03912 for j17j277ndo 130 3971 CW aw 221 aikajklajj Consider now a symmetric positive de nite tridiagonal matrix whose subdiagonal el ements are b27193 71 and whose diagonal elements are 0102 0 Create an algorithm in a nonambiguous algorithmic language which overwrites the arrays b and c with the Cholesky triangle of the tridiagonal matrix Use no arrays other than b and c in your algorithm Let A be the matrix of the ve point difference operator on a uniform square mesh of N2 interior points Show that A is nonsingular This is usually done as follows Let u 1111711217u11u127u2277u12u1Ju2J7uU be a vec tor satisfying Au 0 Let um be component of u such that luul Hulloo lgagNluul If I or J 1 or N7 then it follows that uu 0 Show this Otherwise7 it follows that the neighboring values luIiLJl luLJill this Clearly we can repeat this argument as many times as necessary to get that llulloo luIJl luI1Jl lUNJl 0 Hence we have shown that the only solution of Au 0 is u 07 so A is nonsingular 52 Pro le1 Storage Consider only symmetric positive de nite matrices An example from ODEs with periodic BCs d 71 71 71 d 39 39 71 71 71 d 1also called envelope7 skyline7 variable bandwidth 131 where d gt 2 For pro le storage scheme in each row store from 1st nonzero through diagonal d 71 d 71 d 71 d 71 0 0 71 d Another less special example 6 8 9 3 5 7 1 2 4 4 gtlt gtlt 71 4 gtlt gtlt gtlt 71 0 4 gtlt gtlt gtlt 71 0 4 gtlt gtlt gtlt 71 1 0 4 gtlt gtlt gtlt 1 0 0 4 gtlt gtlt 171 0 4 gtlt gtlt 71 1 0 4 gtlt 1 1 4 Use an array of variable length arrays a 4771747717 0747717 0747717 717074 717 07 0747 717 717 0747717 717074 717 7174 How do we reference an element aij where j S 239 am is ammm 7am aka17 7ah where jmm 239 7 az length 7 1 132 The element aij is in alz39Hj 7 jmm ifj 2 jmm otherwise aij 0 In C one can calculate alzllength using pointer arithmetic For a Cholesky factorization A GGT7 the Cholesky triangle G has same pro le as A For the 5 point difference operator on W gtlt W mesh7 pro le storage reduces storage by g and computation by over banded storage Pro le storage is commonly used for stiffness matrices arising from the application of the nite element method In the nite element method the spatial domain is partitioned into triangular or rectangular elements Unknown values are sought at nodes in the domain The structure of the matrix is such that the 2397jth entry of the stiffness matrix is nonzero only if nodes 239 and j are in the same element of the spatial domain 53 Sparse Storage Saad Section 34 X X X X X X X X X X X X A possible storage scheme a X7 X7 X7X7 X7X7X7 X7X7 X7X7 X ja 17375717 27 37 1747275747 For a symmetric matrix Additional storage needed for ll in An example of no pivoting gtlt gtlt gtlt gtlt gtlt gtlt gtlt gtlt gtlt 0 gtlt F F 0 gtlt gtlt gtlt 0 gtlt gtlt gtlt gtlt gtlt gtlt 0 F gtlt F a 0 gtlt gtlt gtlt a 0 0 gtlt gtlt gtlt gtlt 0 F gtlt 0 0 gtlt l X X j l X X j l X X j et cetera Review questions 1 For general sparse storage of a matrix how would you access its 2397jth element 2 Do symbolic Gaussian elimination on a sparse matrix to determine llin 133 Exercises 1 What is the sparsity pattern for the Cholesky triangle of a symmetric matrix whose lower half has the following sparsity pattern X X gtltgtltgtlt X X gtltgtltgtlt X XX Note we leave blank those elements known to be zero we put gtlt7s elsewhere 2 For a matrix A with sparsity pattern determine the sparsity patterns of its L7 U factors Assume no permutation of the rows or columns of A 54 Reordering Saad Sections 32 332 pp 82784 Section 36 Problem how to number variablesequations so as to minimize bandwidth or pro le or ll in If a permutation matrix P re orders ie7 re numbers the equations and Q re orders vari ables7 then the new matrix is FAQ We consider symmetric positive de nite matrices rst To preserve symmetry7 choose Q PT7 ie7 maintain the same correspondence between equations and variables The sparsity structure of a symmetric matrix can be represented by an undirected graph graph nodes7 edges nodes 12771 edges 2 jlz 79739 and aij 7 0 134 09 Nodes 239 andj are adjacent if Lj is an edge The llin that occurs during the elimination of the kth unknown can be represented graphically by eliminating node k from the graph and connecting all unconnected nodes that are adjacent to node k Each new edge corresponds to a lled in matrix element PAPT ltgt relabeling the nodes reverse Cuthill McKee ordering A very successful algorithm for relabeling for pro le storage is the reverse Cuthill MeKee algorithm label a node with an n l n 7 1 forkn7n71772do While there are unlabeled nodes adjacent to k label with an l a node of lowest degree which is unlabeled amp adjacent to k ll71 The degree of a node is the number of adjacent nodes 0 o 9 9 nested dissection ordering Consider a 64 gtlt 64 matrix with the following graph gt 135 separator ktwo disconnected graphs Order the nodes in the left part 17 27 732 Order the nodes in the right part 337347 756 Order the nodes in the separator 577587 7 64 Apply this idea recursively to the left and right part hence nested 19 20 28 24 54 46 5 49 I21 22 2725 6347 525 17 18 26 23 62 45 51 48 29 30 31 32 61 5616116354439 7 8 15121 593d 440 3 414 10 58 34 4138 1 213 9 57 33 4 37 If one programs the Cholesky factorization to exploit all sparsity7 then for W gtlt W mesh 31 storage gm logz n 0717 829 rnults ns2 On log 136 minimum degree ordering This ordering is most popular for direct solvers for symmetric positive de nite systems the previous two orderings are preferred for iterative methods At each step of Gaussian elimination7 the pivot is chosen to correspond to a node of minimum degree in the elimination graph For a general matrix we may want to pivot interchange rows andor columns for nu merical stability avoid small pivots For a sparse matrix we also want to pivot to preserve sparseness eg7 X Xl llin X 0 no llin X X 72 con icting considerations Compromise threshold pivotingipivot for sparseness subject to the condition that the pivot is not too small lpivot l 2 Tl largest possible pivot l7 eg7 739 01 The choice 739 1 is the usual partial pivoting Review questions 1 What is the primary purpose of reordering a banded matrix a general sparse matrix 2 Explain how to construct nodes and edges of a graph representing the sparsity structure of a symmetric matrix 3 What happens to the graph representing the sparsity structure of a symmetric matrix A if an unknown is eliminated from the linear system Ax b lnterpret the change in the graph as a change in the matrix 4 What matrix operation corresponds to a relabeling ofthe nodes of a graph representing the sparsity structure of a symmetric matrix 5 What is the degree of a node the number of adjacent nodes 6 What is a separator of a connected graph 7 Explain the idea of nested dissection 8 lllustrate the idea of nested dissection on an example 9 If the graph of a symmetric positive de nite matrix is a W by W mesh7 the number of operations required for nested dissection to compute a Cholesky factorization is proportional to what power of n 10 For a general sparse matrix what are two reasons for reordering 11 What is threshold pivoting When does it reduce to partial pivoting 137 Exercises 1 Given below is a graph representing the sparsity structure of an 8 by 8 symmetric pos itive de nite matrix If one applies Cholesky decomposition or equivalently Gaussian elimination7 how many of the zeros in the lower triangular part of the matrix ll in Show your work 2 Consider a symmetric positive de nite matrix A having the following sparsity structure XX X X gtltgtltgtltgtlt gtltgtltgtlt gtltgtltgtltgtlt gtlt gtltgtltgtlt gtlt gtltgtltgtltgtlt gtltgtltgtlt gtltgtltgtltgtlt X X XX By using a graph to represent the sparsity structure of A7 determine a permutation matrix P such that the bandwidth of PAPT is minimized 3 What would be the sparsity pattern of the following symmetric matrix if we were to re order rows and columns using the reverse Cuthill McKee algorithm X X X XXX XXX XXX XXX XXX XXX gtltgtlt 138 55 Fast Direct Methods oddeven reduction for tridiagonal matrices augmented matrix X X XXX XXX XXX XXX XXX XXX XXX XXX gtltgtlt gtltgtltgtltgtltgtltgtltgtltgtltgtltgtlt pivots fork2747678do use equation k to eliminate M from equations k 7 1 and k 1 use equation 10 to eliminate 10 from equation 9 gtlt0gtlt gtlt xxx x gtlt0gtlt0gtlt gtlt xxx x gtlt0gtlt0gtlt gtlt x x x gtlt0gtlt0gtlt gtlt xxx x gtlt0gtlt0gtlt x x Note that equations 17 37 57 77 9 involve only 17 37 57z7797a reduced problem If we can solve this7 we can get 27 x47z67x87 10 from equations 27 47 67 87 10 The reduced problem is X X XXX XXX XXX gtltgtlt gtltgtltgtltgtltgtlt 139 Eliminate 2 and m using equation 2 and 4 gtlt 0 gtlt gtlt x x x x gtlt 0 gtlt 0 gtlt gtlt x x x x gtlt 0 gtlt gtlt If we can solve equations 17 37 57 then we can also get 27 n X X X X X X X X X X X X X X X X X X X X X X X X X Solve this block cyclic reduction This is a generalization of oddeven reduction to block tridiagonal matrices7 eg7 D F 0 F D F 0 F DJ where DF FD7 eg7 4 71 0 D 71 4 71 0 71 4 andF 7 Fourier analysis If we have a known spectral decomposition A QAQT7 we can write A71 QAilQT39 Multiplication by QT and Q each require 712 multiplies normally7 but for special Q7 the cost is only 0nlogn using the FFT assuming 71 has a factorization into small primes The FFT is based on a factorization of Q into Ologn sparse matrices 140 fast direct methods block Cyclic Reduction 0nlog2 n ops Fourier Analysis 0n logz n ops combination7 FACRU 0nloglog n ops These are applicable to discretized separable BlP7s7 ie7 separation of variables Software FISHPACK Review questions 1 What is the idea of oddeven reduction for a tridiagonal matrix 2 What do we call the generalization of oddeven reduction to a block tridiagonal matrix 3 What is the operation count of a fast direct method for solving systems of linear equations Name two examples of fast direct methods 141 142

### BOOM! Enjoy Your Free Notes!

We've added these Notes to your profile, click here to view them now.

### You're already Subscribed!

Looks like you've already subscribed to StudySoup, you won't need to purchase another subscription to get this material. To access this material simply click 'View Full Document'

## Why people love StudySoup

#### "I was shooting for a perfect 4.0 GPA this semester. Having StudySoup as a study aid was critical to helping me achieve my goal...and I nailed it!"

#### "I signed up to be an Elite Notetaker with 2 of my sorority sisters this semester. We just posted our notes weekly and were each making over $600 per month. I LOVE StudySoup!"

#### "I was shooting for a perfect 4.0 GPA this semester. Having StudySoup as a study aid was critical to helping me achieve my goal...and I nailed it!"

#### "Their 'Elite Notetakers' are making over $1,200/month in sales by creating high quality content that helps their classmates in a time of need."

### Refund Policy

#### STUDYSOUP CANCELLATION POLICY

All subscriptions to StudySoup are paid in full at the time of subscribing. To change your credit card information or to cancel your subscription, go to "Edit Settings". All credit card information will be available there. If you should decide to cancel your subscription, it will continue to be valid until the next payment period, as all payments for the current period were made in advance. For special circumstances, please email support@studysoup.com

#### STUDYSOUP REFUND POLICY

StudySoup has more than 1 million course-specific study resources to help students study smarter. If you’re having trouble finding what you’re looking for, our customer support team can help you find what you need! Feel free to contact them here: support@studysoup.com

Recurring Subscriptions: If you have canceled your recurring subscription on the day of renewal and have not downloaded any documents, you may request a refund by submitting an email to support@studysoup.com

Satisfaction Guarantee: If you’re not satisfied with your subscription, you can contact us for further help. Contact must be made within 3 business days of your subscription purchase and your refund request will be subject for review.

Please Note: Refunds can never be provided more than 30 days after the initial purchase date regardless of your activity on the site.