MATH 417 Lecture 20

« previous | Thursday, April 3, 2014 | next »

Gaussian Elimination with Partial Pivoting

Step 0. System of equations Ax=b

$\left[{\begin{array}{ccc|c}0&2&2&2\\1&2&2&1\\1&1&2&3\end{array}}\right]$

Step 1. cell [1,1] is 0. That's unfortunate, so Switch rows 1 & 2

$\left[{\begin{array}{ccc|c}1&2&2&1\\0&2&2&2\\1&1&2&3\end{array}}\right]$

Step 2. subtract 1*(top row) from (bottom row)

$\left[{\begin{array}{ccc|c}1&2&2&1\\0&2&2&2\\0&-1&0&2\end{array}}\right]$

Step 3. multiply (second row) by 1/2

$\left[{\begin{array}{ccc|c}1&2&2&1\\0&1&1&1\\0&-1&0&2\end{array}}\right]$

Step 4. multiply (bottom row) by -1

$\left[{\begin{array}{ccc|c}1&2&2&1\\0&1&1&1\\0&1&0&-2\end{array}}\right]$

Step 5. subtract (second row) from (bottom row)

$\left[{\begin{array}{ccc|c}1&2&2&1\\0&1&1&1\\0&0&-1&-3\end{array}}\right]$

Step 5. multiply (bottom row) by -1

$\left[{\begin{array}{ccc|c}1&2&2&1\\0&1&1&1\\0&0&1&3\end{array}}\right]=U$

Step 7. Perform back substitution to get

${\begin{aligned}x_{3}&=3\\x_{2}&=-2\\x_{1}&=-1\end{aligned}}$

At each step, we left multiply by an elementary matrix. This elementary matrix is the identity matrix with the same operation performed on it:

1. $I={\begin{bmatrix}1&0&0\\0&1&0\\0&0&1\end{bmatrix}}\quad \longrightarrow \quad E_{1}={\begin{bmatrix}0&1&0\\1&0&0\\0&0&1\end{bmatrix}}$
2. $I={\begin{bmatrix}1&0&0\\0&1&0\\0&0&1\end{bmatrix}}\quad \longrightarrow \quad E_{2}={\begin{bmatrix}1&0&0\\0&1&0\\-1&0&1\end{bmatrix}}$
3. $I={\begin{bmatrix}1&0&0\\0&1&0\\0&0&1\end{bmatrix}}\quad \longrightarrow \quad E_{3}={\begin{bmatrix}1&0&0\\0&{\frac {1}{2}}&0\\0&0&1\end{bmatrix}}$

and so on...

Observation: Each $E_{k}$ has two possibilities:

1. row permutations
2. lower triangular matrix

Note, if we perform row permutations $P_{1}\cdot P_{2}\dots P_{n}$ at the start, then we can proceed with Gaussian elimination without pivoting.

Let ${\tilde {A}}:=P_{n}\,\dots \,P_{2}\,P_{1}\,A$

Now since $E_{k}\dots E_{3}\,E_{2}\,{\tilde {A}}=U$, we have $A=(E_{k}\,\dots \,E_{3}\,E_{2})^{-1}\,U$, so $(E_{k}\dots E_{3}\,E_{2})^{-1}=L=E_{2}^{-1}\,E_{3}^{-1}\,\dots \,E_{k}^{-1}$

What is $E_{i}^{-1}$?

1. row-scaling: $E={\begin{bmatrix}\lambda &0&0&0\\0&1&0\\0&0&1\end{bmatrix}}\quad \longrightarrow \quad E^{-1}={\begin{bmatrix}{\frac {1}{\lambda }}&0&0\\0&1&0\\0&0&1\end{bmatrix}}$
2. row-addition: $E={\begin{bmatrix}1&0&0\\0&1&0\\\alpha &0&1\end{bmatrix}}\quad \longrightarrow \quad E^{-1}={\begin{bmatrix}1&0&0\\0&1&0\\-\alpha &0&1\end{bmatrix}}$

To multiply these inverses, take superposition of matrices

• combine nonzero entries not along diagonal
• multiply corresponding entries on diagonals

${\begin{bmatrix}1&0&0\\0&1&0\\-1&0&1\end{bmatrix}}\cdot {\begin{bmatrix}1&0&0\\0&{\frac {1}{2}}&0\\0&0&1\end{bmatrix}}={\begin{bmatrix}1&0&0\\0&{\frac {1}{2}}&0\\-1&0&1\end{bmatrix}}$

Therefore, given $A$, if we can run gaussian elimination without pivoting, then $A=L\,U$ and $L={\begin{bmatrix}\ell _{11}&0\\\ell _{ij}&l_{nn}\end{bmatrix}}$; $U={\begin{bmatrix}1&u_{ij}\\0&1\end{bmatrix}}$

In practice, we can store $L$ and $U$ in the same matrix: ${\begin{bmatrix}\ell _{11}&u_{ij}\\\ell _{ij}&\ell _{nn}\end{bmatrix}}$

In contrast, we can do Gaussian elimination without pivoting (i.e. leave diagonal entries ≠ 1). In this case, our LU decomposition is

$L={\begin{bmatrix}1&0\\\ell _{ij}&1\end{bmatrix}}$ $U={\begin{bmatrix}u_{11}&u_{ij}\\0&u_{nn}\end{bmatrix}}$

If we have $A=A^{T}$, can we get $A={\tilde {L}}={\tilde {L}}^{T}$? we should be able to (see later)

Diagonalization

${\begin{aligned}L\,U&={\begin{bmatrix}1&0&0\\\ell _{21}&1&0\\\ell _{31}&\ell _{32}&1\end{bmatrix}}*{\begin{bmatrix}u_{11}&u_{12}&u_{13}\\0&u_{22}&u_{23}\\0&0&u_{33}\end{bmatrix}}\\&={\begin{bmatrix}1&0&0\\\ell _{21}&1&0\\\ell _{31}&\ell _{32}&1\end{bmatrix}}*{\begin{bmatrix}u_{11}&0&0\\0&u_{22}&0\\0&0&u_{33}\end{bmatrix}}*{\begin{bmatrix}1&{\frac {u_{12}}{u_{11}}}&{\frac {u_{13}}{u_{11}}}\\0&1&{\frac {u_{23}}{u_{22}}}\\0&0&1\end{bmatrix}}\\&=L\,\Delta \,{\tilde {U}}\end{aligned}}$

In this case, the diagonal entries of $\Delta$ are the eigenvalues of the matrix.

Theorem: If there are $n$ eigenvalues, then there will be $n$ independent eigenvectors.

LDL Decomposition / Choleski Factorization

Above we saw a decomposition of the form $A=L\,D\,L^{T}$. This is called LDL decomposition.

We can alternatively write this as $L\,{\sqrt {D}}\,{\sqrt {D}}\,L^{T}=(L\,{\sqrt {D}})\,(L\,{\sqrt {D}})^{T}={\tilde {L}}\,{\tilde {L}}^{T}$

So in the example above, we have $\ell _{21}={\frac {u_{12}}{u_{11}}}$

Special Matrices

Diagonally Dominant

Given $A=(a_{ij})$, for every row, we have $\left|a_{ii}\right|\geq \sum _{j\neq i}\left|a_{ij}\right|\geq 0$

"strictly" variant: $\left|a_{ii}\right|>\sum _{j\neq i}\left|a_{ij}\right|\geq 0$

We will consider only this "strict" case

Positive Definite

For any ${\vec {x}}\neq {\vec {0}}$, we have ${\vec {x}}^{T}\,A\,{\vec {x}}\geq 0$

also called non-negative definite

"strictly" variant: ${\vec {x}}^{T}\,A\,{\vec {x}}>0$

We will consider only this "strict" case

Theorem

Theorem. If $A$ is strictly diagonally dominant or strictly positive definite, then $A$ is non-singular.

Proof. Assume to the contrary that $A$ is singular. Then $A\,{\vec {x}}={\vec {0}}$ has a non-trivial solution (${\vec {x}}\neq {\vec {0}}$)

Case 1: $A$ is diagonally dominant. Take ${\vec {x}}={\begin{bmatrix}x_{1}\\x_{2}\\\vdots \\x_{n}\end{bmatrix}}$. Find $x_{j}$ such that $\left|x_{j}\right|=\max _{1\leq i\leq n}{\left|x_{i}\right|}$

Solve the system for $a_{jj}\,x_{j}=-\sum _{k=1;k\neq j}^{n}a_{jk}\,x_{k}$. Hence

${\begin{aligned}\left|a_{jj}\right|\,\left|x_{j}\right|&=\left|\sum _{k\neq j}a_{jk}\,x_{k}\right|\leq \sum _{k\neq j}\left|a_{jk}\right|\,\left|x_{k}\right|\\\left|a_{jj}\right|\leq \sum _{k\neq j}\left|a_{kj}\right|\,\left|{\frac {x_{k}}{x_{j}}}\right|\end{aligned}}$

Observe that $\left|{\frac {x_{k}}{x_{j}}}\right|\leq 1$, so

$\left|a_{jj}\right|\leq \sum _{k\neq j}\left|a_{kj}\right|\,\left|{\frac {x_{k}}{x_{j}}}\right|\leq \sum _{k\neq j}\left|a_{jk}\right|<a_{jj}$

→CONTRADICTION

Case 2: $A$ is positive definite and $A\,{\vec {x}}={\vec {0}}$ for some ${\vec {x}}\neq {\vec {0}}$.

$0<{\vec {x}}^{T}\,A\,{\vec {x}}={\vec {x}}^{T}\,{\vec {0}}=0$

If $\lambda _{1}$ is an eigenvalue of $A$, then $A\,{\vec {x_{1}}}=\lambda _{1}\,x_{1}$ for some $x_{1}$.

${\vec {x}}_{1}\,A{\vec {x}}_{1}={\vec {x}}_{1}^{T}\,\lambda _{1}\,{\vec {x}}_{1}=\lambda _{1}\,{\vec {x}}_{1}^{T}\,{\vec {x}}_{1}>0$

Hence $\lambda _{1}>0$.

...

quod erat demonstrandum

Going back to diagonal dominance. If we do Gaussian elimination, we are guaranteed that we will not need any pivoting