Pivoting for LU Factorization - UPS

1 Pivoting forLUFactorizationMatthew W. ReidApril 21, 2014 University of Puget SoundE-mail: (C) 2014 Matthew W. Reid. Permission is granted to copy, distribute and/or modifythis document under the terms of the GNU Free Documentation License, Version or any laterversion published by the Free Software Foundation; with no Invariant Sections, no Front-CoverTexts, and no Back-Cover Texts. A copy of the license is included in the section entitled GNUFree Documentation License .1 INTRODUCTION11 IntroductionPivoting forLUfactorization is the process of systematically selecting pivots for Gaussian elimina-tion during theLUfactorization of a matrix. TheLUfactorization is closely related to Gaussianelimination, which is unstable in its pure form. To guarantee the elimination process goes to com-pletion, we must ensure that there is a nonzero pivot at every step of the elimination process.

2 Thisis the reason we need Pivoting when computingLUfactorizations. But we can do more with piv-oting than just making sure Gaussian elimination completes. We can reduce roundoff errors duringcomputation and make our algorithm backward stable by implementing the right Pivoting on the matrixA, someLUdecompositions can become numerically unstable if relativelysmall pivots are used. Relatively small pivots cause instability because they operate very similarto zeros during Gaussian elimination. Through the process of Pivoting , we can greatly reduce thisinstability by ensuring that we use relatively large entries as our pivot elements. This prevents largefactors from appearing in the computedLandU, which reduces roundoff errors during computa-tion. The three Pivoting strategies I am going to discuss are partial Pivoting , complete Pivoting ,and rook Pivoting , along with an explanation of the benefits obtained from using each Permutation MatricesThe process of swapping rows and columns is crucial to Pivoting .

3 Permutation matrices are usedto achieve this. A permutation matrix is the identity matrix with interchanged rows. When thesematrices multiply another matrix they swap the rows or columns of the matrix. Left multiplicationby a permutation matrix will result in the swapping of rows while right multiplication will swapcolumns. For example, in order to swap rows 1 and 3 of a matrixA, we right multiply by a permu-tation matrixP, which is the identity matrix with rows 1 and 3 swapped. To interchange columns1 and 3 of the matrixA, we left multiply by the same permutation matrixP. At each stepkofGaussian elimination we encounter an optional row interchange, represented byPk, before eliminat-ing entries via the lower triangular elementary matrixMk. If there is no row interchange required,Pkis simply the identity matrix. Think of this optional row interchange as switching the order ofthe equations when presented with a system of equations.

4 For Gaussian elimination to terminatesuccessfully, this process is necessary when encountering a zero pivot. When describing the differentpivoting strategies, I will denote a permutation matrix that swaps rows withPkand will denote apermutation matrix that swaps columns by refering to the matrix asQk. When computing theLUfactorizations of matrices, we will routinely pack the permutation matrices together into a singlepermutation matrix. They are simply a matrix product of all the permutation matrices used toachieve the Factorization . I will define these matrices computingPA=LU,P=PkPk (1)wherekis the index of the computingPAQ=LU,Q= 1Qk(2)3 ROLE OF PIVOTING2wherekis the index of the Role of PivotingThe role of Pivoting is to reduce instability that is inherent in Gaussian elimination. Instabilityarises only if the factorLorUis relatively large in size compared toA. Pivoting reduces instabliltyby ensuring thatLandUare not too big relative toA.

5 By keeping all of the intermediate quantitesthat appear in the elimination process a managable size, we can minimize rounding errors. As aconsequence of Pivoting , the algorithm for computing theLUfactorization is backward stable. Iwill define backward stability in the upcoming Zero PivotsThe first cause of instability is the situation in which there is a zero in the pivot position. With azero as a pivot element, the calculation cannot be completed because we encounter the error of di-vision by zero. In this case, Pivoting is essential to finishing the computation. Consider the [0 11 2]When trying to compute the factorsLandU, we fail in the first elimination step because ofdivision by Small PivotsAnother cause of instability in Gaussian elimination is when encountering relatively small elimination is guaranteed to succeed if row or column interchanges are used in order toavoid zero pivots when using exact calculations.

6 But, when computing theLUfactorization numer-ically, this is not necessarily true. When a calculation is undefined because of division by zero, thesame calculation will suffer numerical difficulties when there is division by a nonzero number thatis relatively [10 20112]When computing the factorsLandU, the process does not fail in this case because there is nodivision by [1010201],U=[10 20102 10 20]3 ROLE OF PIVOTING3 When these computations are performed in floating point arithmetic, the number 2 10 20isnot represented exactly but will be rounded to the nearest floating point number which we will sayis 10 20. This means that our matrices are now floating point matricesL andU .L =[1010201],U =[10 2010 1020]The small change we made inUto getU shows its significance when we computeL U .L U =[10 20110]The matrixL U should be equal toA, but obviously that is not the case. The 2 in position(2,2) of matrixAis now 0.

7 Also, when trying to solve a system such asAx=busing theLUfactorization, the factorsL U would not give you a correct answer. TheLUfactorization was astable computation but not backward times we computeLUfactorizations in order to solve systems of equations. In that case,we want to find a solution to the vector equationAx=b, whereAis the coefficient matrix ofinterest andbis the vector of constants. Therefore, for solving this vector equation specifically, thedefinition of backward stability is as algorithm is stable for a class of matricesCif for every matrixA C, thecomputed solution by the algorithm is the exact solution to a nearby problem. Thus, for a linearsystem problemAx=ban algorithm is stable for a class of matricesCif for everyA Cand for eachb, it produces acomputed solution xthat satisfies(A+E) x=b+ bfor someEand b, where (A+E) is close toAandb+ bis close without Pivoting is not backward stable because the computed solution xfor the problemL U x=bis not the exact solution to a nearby problemAx=b.

8 The computedsolution toL U x=b, withb=[13], is x=[01]. Yet the correct solution isx [11].4 PARTIAL PIVOTING44 Partial PivotingThe goal of partial Pivoting is to use a permutation matrix to place the largest entry of the firstcolumn of the matrix at the top of that first column. For ann nmatrixB, we scannrows of thefirst column for the largest value. At stepkof the elimination, the pivot we choose is the largest ofthen (k+ 1) subdiagonal entries of columnk, which costsO(nk) operations for each step of theelimination. So for an nmatrix, there is a total ofO(n2) comparisons. Once located, this entry isthen moved into the pivot positionAkkon the diagonal of the matrix. So in the first step the entryis moved into the (1,1) position of matrixB. We interchange rows by multiplyingBon the leftwith a permutation matrixP. After we multiply matrixBbyPwe continue theLUfactorizationand use our new pivot to clear out the entries below it in its column in order to obtain the uppertriangular matrixU.

9 The general Factorization of an nmatrixAto a lower triangular matrixLand an upper triangular matrixUusing partial Pivoting is represented by the following k= (Pn 1 Pk+1)Mk(Pk+1 Pn 1)(3)(Mn 1Mn 2 M2M1) 1=L(4)Mn 1Pn 1Mn 2Pn 2 M2P2M1P1A=U(5)In order to get a visual sense of how pivots are selected, let represent the largest value in thematrixBand the bold represent entries of the matrix being scanned as potential xx x xxx x x 1x x xxx x x 1x x x0xx x0xx x0 2x x 1x x x0 2x x00xx00 3x 1x x x0 2x x00 3x000x Aftern 1 permutations, the last permutation becomes trivial as the column of interest hasonly a single entry. As a numerical example, consider the 1 2 42 1 33 2 4 In this case, we want to use a permutation matrix to place the largest entry of the first columninto the positionA11. For this matrix, this means we want 3 to be the first entry of the first we use the permutation matrixP1,P1= 0 0 10 1 01 0 0 ,P1A= 3 2 42 1 31 2 4 5 COMPLETE PIVOTING5 Then we use the matrixM1to eliminate the bottom two entries of the first 10 0 2/3 1 0 1/3 0 1 ,M1P1A= 3240 1/3 1/304/38/3 Then, we use the permutation matrixP2to move the larger of the bottom two entries in columntwo ofAinto the 1 0 00 0 10 1 0 ,P2M1P1A= 32404/38/30 1/3 1/3 As we did before, we are going to use matrix multiplication to eliminate entries at the bottomof column 2.

10 We use matrixM2to achieve 1000100 1/4 1 ,M2P2M1P1A= 3240 4/3 8/3001 This upper triangular matrix product is ourUmatrix. TheLmatrix comes from the matrixproduct (M2P2M1P1) (M2P2M1P1) 1= 1001/3102/3 1/4 1 To check, we can show thatPA=LU, whereP= 0 0 11 0 00 1 0 1 2 42 1 33 2 4 = 1001/3102/3 1/4 1 3240 4/3 8/3001 =LUThrough the process of partial Pivoting , we get the smallest possible quantites for the eliminationprocess. As a consequence, we will have no quantity larger than 1 in magnitude which minimizesthe chance of large errors. Partial Pivoting is the most commmon method of Pivoting as it is thefastest and least costly of all the Pivoting Complete PivotingWhen employing the complete Pivoting strategy we scan for the largest value in the entire subma-trixAk:n,k:n, wherekis the number of the elimination, instead of just the next subcolumn.