3.3: Partial Pivoting
( \newcommand{\kernel}{\mathrm{null}\,}\)
When performing Gaussian elimination, the diagonal element that one uses during the elimination procedure is called the pivot. To obtain the correct multiple, one uses the pivot as the divisor to the elements below the pivot. Gaussian elimination in this form will fail if the pivot is zero. In this situation, a row interchange must be performed.
Even if the pivot is not identically zero, a small value can result in big roundoff errors. For very large matrices, one can easily lose all accuracy in the solution. To avoid these round-off errors arising from small pivots, row interchanges are made, and this technique is called partial pivoting (partial pivoting is in contrast to complete pivoting, where both rows and columns are interchanged). We will illustrate by example the LU decomposition using partial pivoting. Consider
A=\left(\begin{array}{rrr} -2 & 2 & -1 \\ 6 & -6 & 7 \\ 3 & -8 & 4 \end{array}\right) \nonumber
We interchange rows to place the largest element (in absolute value) in the pivot, or a_{11}, position. That is,
\mathrm{A} \rightarrow\left(\begin{array}{rrr} 6 & -6 & 7 \\ -2 & 2 & -1 \\ 3 & -8 & 4 \end{array}\right)=\mathrm{P}_{12} \mathrm{~A} \nonumber
where
P_{12}=\left(\begin{array}{lll} 0 & 1 & 0 \\ 1 & 0 & 0 \\ 0 & 0 & 1 \end{array}\right) \nonumber
is a permutation matrix that when multiplied on the left interchanges the first and second rows of a matrix. Note that P_{12}^{-1}=P_{12}. The elimination step is then
\mathrm{P}_{12} \mathrm{~A} \rightarrow\left(\begin{array}{rrc} 6 & -6 & 7 \\ 0 & 0 & 4 / 3 \\ 0 & -5 & 1 / 2 \end{array}\right)=\mathrm{M}_{1} \mathrm{P}_{12} \mathrm{~A} \nonumber
where
M_{1}=\left(\begin{array}{ccc} 1 & 0 & 0 \\ 1 / 3 & 1 & 0 \\ -1 / 2 & 0 & 1 \end{array}\right) \nonumber
The final step requires one more row interchange:
\mathrm{M}_{1} \mathrm{P}_{12} \mathrm{~A} \rightarrow\left(\begin{array}{rrc} 6 & -6 & 7 \\ 0 & -5 & 1 / 2 \\ 0 & 0 & 4 / 3 \end{array}\right)=\mathrm{P}_{23} \mathrm{M}_{1} \mathrm{P}_{12} \mathrm{~A}=\mathrm{U} \nonumber
Since the permutation matrices given by P are their own inverses, we can write our result as
\left(P_{23} M_{1} P_{23}\right) P_{23} P_{12} A=U \nonumber
Multiplication of M on the left by P interchanges rows while multiplication on the right by P interchanges columns. That is,
P_{23}\left(\begin{array}{cll} 1 & 0 & 0 \\ 1 / 3 & 1 & 0 \\ -1 / 2 & 0 & 1 \end{array}\right) P_{23}=\left(\begin{array}{ccc} 1 & 0 & 0 \\ -1 / 2 & 0 & 1 \\ 1 / 3 & 1 & 0 \end{array}\right) P_{23}=\left(\begin{array}{ccc} 1 & 0 & 0 \\ -1 / 2 & 1 & 0 \\ 1 / 3 & 0 & 1 \end{array}\right) \nonumber
The net result on M_{1} is an interchange of the nondiagonal elements 1 / 3 and -1 / 2.
We can then multiply by the inverse of \left(\mathrm{P}_{23} \mathrm{M}_{1} \mathrm{P}_{23}\right) to obtain
P_{23} P_{12} A=\left(P_{23} M_{1} P_{23}\right)^{-1} U \nonumber
which we write as
\mathrm{PA}=\mathrm{LU} \nonumber
Instead of L, MATLAB will write this as
\mathrm{A}=\left(\mathrm{P}^{-1} \mathrm{~L}\right) \mathrm{U} \nonumber
For convenience, we will just denote \left(\mathrm{P}^{-1} \mathrm{~L}\right) by \mathrm{L}, but by \mathrm{L} here we mean a permutated lower triangular matrix.
For example, in MATLAB, to solve A \mathbf{x}=\mathbf{b} for \mathbf{x} using Gaussian elimination, one types
\mathrm{x}=\mathrm{A} \backslash \mathrm{b} \nonumber
where \solves for \mathbf{x} using the most efficient algorithm available, depending on the form of \mathrm{A}. If \mathrm{A} is a general n \times n matrix, then first the LU decomposition of \mathrm{A} is found using partial pivoting, and then x is determined from permuted forward and backward substitution. If A is upper or lower triangular, then forward or backward substitution (or their permuted version) is used directly.
If there are many different right-hand-sides, one would first directly find the LU decomposition of A using a function call, and then solve using . That is, one would iterate for different \mathbf{b}^{\prime} s the following expressions:
\begin{aligned} &{[\mathrm{LU}]=l u(\mathrm{~A})} \\ &\mathrm{y}=\mathrm{L} \backslash \mathrm{b} \\ &\mathrm{x}=\mathrm{U} \backslash y \end{aligned} \nonumber
where the second and third lines can be shortened to
\mathrm{x}=\mathrm{U} \backslash(\mathrm{L} \backslash \mathrm{b}) \nonumber
where the parenthesis are required. In lecture, I will demonstrate these solutions in MATLAB using the matrix A=[-2,2,-1 ; 6,-6,7 ; 3,-8,4] ; which is the example in the notes.