# 2.3. Partial Pivoting#

**Last Revised September 23, 2024**

**References:**

Section 2.4.1

*Partial Pivoting*of [Sauer, 2022].Section 6.2

*Pivoting Stratgies*of [Burden*et al.*, 2016].Section 7.1 of [Chenney and Kincaid, 2012].

Some references describe the method of *scaled* partial pivoting, but here we present instead a version without the “scaling”, because not only is it simpler, but modern research shows that it is esentially always as good, once the problem is set up in a “sane” way.

## 2.3.1. Introduction#

The basic row reduction method can fail due to divisoion by zero (and to have very large rouding errors when a denominator is extremely close to zero.
A more robust modification is to swap the order of the equations to avaid these problems: *partial pivoting*.
Here we look at a particularly robust version of this strategy, *Maximal Element Partial Pivoting*.

```
# As in recent sections, we import some items from modules individually, so they can be used by "first name only".
from numpy import array
```

## 2.3.2. What can go wrong with Naive Gaussian Elimination?#

We have noted two problems with the naive algorithm for Gaussian elimination: total failure due the division be zero, and loss of precision due to dividing by very small values — or more preciselt calculations the lead to intermediate values far larger than the final results.
The culprits in all cases are the same: the denominators are first the *pivot elements* \(a_{k,k}^{(k-1)}\) in evaluation of \(l_{i,k}\) during row reduction and then the \(u_{k,k}\) in back substitution.
Further, those \(a_{k,k}^{(k-1)}\) are the final updated values at indices \((k,k)\), so are the same as \(u_{k,k}\).
Thus it is exactly these *main diagonal elements* that we must deal with.

## 2.3.3. The basic fix: Partial Pivoting#

The basic strategy is that at step \(k\), we can swap equation \(k\) with any equation \(i\), \(i > k\).
Note that this involves swapping those rows of array `A`

and also those elements of the array `b`

for the right-hand side: \(b_k \leftrightarrow b_i\).

This approach of swapping equations (swapping rows in arrays `A`

and `b`

) is called **pivoting**, or more specifically **partial pivoting**, to distinguish from the more elaborate strategy where to columns of `A`

are also reordered (which is equivalent to reordering the unknowns in the equations).
The row that is swapped with row \(k\) is sometimes called the **pivot row**, and the new denominator is the corresponding **pivot element**.

This approach is robust so long as one is using exact arithmetic: it works for any well-posed system because so long as the \(Ax = b\) has a unique solution — so that the original matrix \(A\) is non-singular — at least one of the \(a_{i,k}^{(k-1)}, i \geq k\) will be non-zero, and thus the swap will give a new element in position \((k,k)\) that is non-zero. (I will stop caring about superscripts to distinguish updates, but if you wish to, the elements of the new row \(k\) could be called either \(a_{k,j}^{(k)}\) or even \(u_{k,j}\), since those values are in their final state.)

## 2.3.4. Handling rounding error: Maximal Element Partial Pivoting#

The final refinement is to seek the smallest possible magnitudes for intermediate values, and thus the smallest absolute errors in them, by making the multipliers \(l_{i,k}\) small, in turn by making the denominator \(a_{k,k}^{(k-1)} = u_{k,k}\) as large as possible in magnitude:

*At step \(k\), choose the pivot row \(p_k \geq k\) so that \(|a_{p_k,k}^{(k-1)}| \geq |a_{i,k}^{(k-1)}|\) for all \(i \geq k\). If there is more that one such element of largest magnitude, use the lowest value:
in particular, if \(p_k = k\) works, use it and do not swap!*

## 2.3.5. Swapping rows in Python#

I will not give a detailed algorithm for this, since we will soon implement an even better variant.

However, here are some notes on swapping values and how to avoid a possible pitfall.

### Exercise 1, on Python coding#

a) Explain why we cannot just swap the relevant elements of rows \(k\) and \(p\) with:

```
for j in range(k,n):
A[k,j] = A[p,j]
A[p,j] = A[k,j]
```

or with vectorized “slicing”:

```
A[k,k:] = A[p,k:]
A[p,k:] = A[k,k:]
```

Describe what happens instead.

b) A common strategy to avoid this problem uses an intermediate temporary copy of each value being “moved”.
This can be combined with slicing, but be careful: arrays (including slices of arrays) must be copied with the method `.copy()`

:

```
temp = A[k,k:].copy()
A[k,k:] = A[p,k:]
A[p,k:] = temp
```

Python also has an elegant alternative for swapping a pair of values:

for j in range(k,n): ( A[k,j] , A[p,j] ) = ( A[p,j] , A[k,j] )

This can also be done with slicing, with care to copy those slices:

```
( A[k,k:] , A[p,k:] ) = ( A[p,k:].copy() , A[k,k:].copy() )
```

### Some demonstrations#

No row reduction is done here, so entire rows are swapped rather than just the elements from column \(k\) onward:

```
A = array([[1. , -6. , 2.],[3. , 5. , -6.],[4. , 2. , 7.]])
n = 3
print(f"Initially,\nA =\n {A}")
```

```
Initially,
A =
[[ 1. -6. 2.]
[ 3. 5. -6.]
[ 4. 2. 7.]]
```

```
k = 0
p = 2
temp = A[k,k:].copy()
A[k,k:] = A[p,k:]
A[p,k:] = temp
print("After swapping rows 1 <-> 3 (row indices 0 <-> 2) using slicing and a temporary row,")
print(f"A =\n {A}")
```

```
After swapping rows 1 <-> 3 (row indices 0 <-> 2) using slicing and a temporary row,
A =
[[ 4. 2. 7.]
[ 3. 5. -6.]
[ 1. -6. 2.]]
```

```
k = 1
p = 2
for j in range(n):
( A[k,j] , A[p,j] ) = ( A[p,j] , A[k,j] )
print(f"After swapping rows 2 <-> 3 using a loop and tuples of elements, no temp,")
print(f"A =\n {A}")
```

```
After swapping rows 2 <-> 3 using a loop and tuples of elements, no temp,
A =
[[ 4. 2. 7.]
[ 1. -6. 2.]
[ 3. 5. -6.]]
```

```
k = 0
p = 1
( A[k,k:] , A[p,k:] ) = ( A[p,k:].copy() , A[k,k:].copy() )
print(f"After swapping rows 1 <-> 2 using tuples of slices, no loop or temp,")
print(f"A =\n {A}")
```

```
After swapping rows 1 <-> 2 using tuples of slices, no loop or temp,
A =
[[ 1. -6. 2.]
[ 4. 2. 7.]
[ 3. 5. -6.]]
```