Diagonal Dominance and Pivoting - Executable Engineering: Interactive Numerical Methods in Python

import numpy as np
import matplotlib.pyplot as plt

try:
    import executable_engineering as exe
except ImportError:
    %pip install -q executable_engineering
    import executable_engineering as exe

The Need for Pivoting¶

The order of equations (and therefore rows in the matrix) is mathematically arbitrary, but the physical order in memory can break Gaussian Elimination.

For example, a system with the following augmented matrix works perfectly with our standard GE algorithm:

\begin{bmatrix} 2 & -1 & 0 & 1\\ -1 & 2 & -1 & 0\\ 0 & -1 & 1 & 0 \end{bmatrix}

(1)

But the exact same system, simply reordered, will fail:

\begin{bmatrix} 0 & -1 & 1 & 0\\ -1 & 2 & -1 & 0\\ 2 & -1 & 0 & 1 \end{bmatrix}

(2)

Because the pivot element in the first row is 0, calculating the multiplier $m_{2,1} = -1 / 0$ results in a division by zero!

The remedy is to pivot (swap) row 1 with another row to put a non-zero element on the diagonal.

The Permutation Matrix¶

Swapping rows in matrix algebra is achieved via a permutation matrix $\mathbf{P}$ (which is simply a row-swapped identity matrix $\mathbf{I}$ ).

To mathematically swap rows 3 and 1 of a matrix, we multiply it by:

\mathbf{P} = \begin{bmatrix} 0 & 0 & 1 \\ 0 & 1 & 0 \\ 1 & 0 & 0 \end{bmatrix}

(3)

Roundoff Error and Small Pivots¶

Even if the pivot isn’t exactly 0, a very small pivot element causes catastrophic roundoff error.

Consider a matrix where the top-left element is some incredibly small number $\epsilon$ :

\begin{bmatrix} \epsilon & -1 & 1 & 0\\ -1 & 2 & -1 & 0\\ 2 & -1 & 0 & 1 \end{bmatrix}

(4)

After the first pass of Gaussian Elimination, the matrix becomes:

\begin{bmatrix} \epsilon & -1 & 1 & 0\\ -1 & 2 - 1/\epsilon & -1 + 1/\epsilon & 0\\ 2 & -1 - 1/\epsilon & 1/\epsilon & 1 \end{bmatrix}

(5)

We can immediately see the danger as $\epsilon \rightarrow 0$ . The terms $1/\epsilon$ become infinitely massive, completely dominating and erasing the true values of the other coefficients due to floating-point truncation (roundoff error).

This motivates our desire to have the largest possible numbers on the diagonal!

Diagonal Dominance¶

A matrix is said to be diagonally dominant if the absolute value of each diagonal element is strictly larger than the sum of the absolute values of the other elements in that row.

|A_{ii}| \ge \sum_{j=1,j \ne i}^{n} |A_{ij}|

(6)

Which of these is diagonally dominant?

A) $$

\begin{bmatrix} -2 & 4 & -1 \\ -1 & -1 & 3 \\ 4 & -2 & 1 \end{bmatrix}

(7)

$$

B) $$

\begin{bmatrix} 4 & -2 & 1 \\ -2 & 4 & -1 \\ -1 & -1 & 3 \end{bmatrix}

(8)

$$

(Answer: Matrix B is diagonally dominant! Notice how Matrix A is just Matrix B without row-pivoting. By swapping the rows of A, we can achieve diagonal dominance!)

The Importance of Diagonal Dominance¶

It can be mathematically proven that if a matrix is diagonally dominant, it will never require pivoting during Gaussian Elimination, and it is inherently immune to the catastrophic $\epsilon \rightarrow 0$ roundoff errors shown above.

Diagonal dominance contributes broadly to numerical stability. It is a critical criteria for matrix factorization methods and ensures the rapid convergence of iterative methods!