< >
 Ch. 2 Sec. 1

2.2 Exact Methods

Solve $A \cdot x = b$ exactly:

2.2.1 Gauss Elimination and Back Substitution

$\begin{eqnarray} \left( \begin{array}{ccccc} a_{11} & a_{12} & . & . & a_{1N} \\ a_{21} & a_{22} & & & \\ . & & . & & \\ . & & & . & \\ . & & & & a_{NN} \end{array} \right) \cdot \left( \begin{array}{c} x_{1} \\. \\.\\. \\x_{N} \end{array} \right) &=& \left( \begin{array}{c} b_{1} \\. \\.\\. \\b_{N} \end{array} \right) \end{eqnarray}$

Convert this to triangular form:

$\begin{eqnarray} \left( \begin{array}{ccccc} a_{11}'& a_{12}' & . & . & . \\ 0 & a_{22}' & & & \\ . & . & . & & \\ . & & . & . & \\ 0 & . & . & 0 & a_{NN}' \end{array} \right) \cdot \left( \begin{array}{c} x_{1} \\. \\.\\. \\x_{N} \end{array} \right) &=& \left( \begin{array}{c} b_{1}' \\. \\.\\. \\b_{N}' \end{array} \right) \end{eqnarray}$

Then solve the system by Back Substitution.

2.2.2 Householder Transformation

Task: Strip column and/or row vectors of some of their rear elements.

Given $A$ $\longrightarrow$ $A'$ triangular, tridiagonal, or otherwise simple.

Let $a$ be a vector (such as a matrix column vector), and $e_{1}$ a unit vector. Define an auxiliary vector $b \equiv a \pm |a| e_{1}$ and normalize it.

Householder matrix:

$P \equiv 1-2 b_{0} \cdot b_{0}^{T}$

See for yourself that in the transformed vector $a' \equiv P \cdot a$ all elements but the first are equal to zero.

Apply this transformation successively to simplify a matrix.

EXAMPLE: Let

$A=\left( \begin{array}{ccc}1&2&3\\ \\4&5&6 \\ \\7&8&9\end{array} \right)$

From $a \equiv (1,4,7)^{T}$ construct the auxiliary vector $b=(-7.124, 4, 7)^{T}$, normalize to $b_{0}=(-0.662, 0.372, 0.651)^{T}$. The Householder matrix is

$P_{1}=\left( \begin{array}{ccc}0.123&0.492&0.862\\ \\0.492&0.724&-0.484 \\ \\0.862&-0.484&0.153\end{array} \right)$

Check: Multiplying $A$ by $P_{1}$ we have

$P_{1} \cdot A= \left( \begin{array}{ccc}8.124&2.708&11.078\\ \\0&4.602&1.464 \\ \\0&-0.696&1.062\end{array} \right)$

Next step: treat the lower right $2 \times 2$ submatrix of $P_{1} \cdot A$. From $a = (4.602, -0.696)^{T}$ construct a $2 \times 2$ Householder matrix which is promoted to a $3 \times 3$ matrix by adding a trivial first line and column:

$P_{2} = \left( \begin{array}{ccc}1&0&0\\ \\0&0.989&-0.149 \\ \\0&-0.149&-0.989\end{array} \right)$

Result:

$P_{2} \cdot P_{1} \cdot A= \left( \begin{array}{ccc}8.124&2.708&11.078\\ \\0&4.655&1.289 \\ \\0&0&-1.269\end{array} \right)$

2.2.3 LU Decomposition

Task: Solve a linear system $A \cdot x=b$

Method (Banachiewicz, Cholesky and Crout): Factorize a matrix into two triangular matrices $L$(ower) and $U$(pper):

$A= L \cdot U$

with

 L$= \left( \begin{array}{cccc} l_{11} & 0 & . & 0 \\ l_{21} & l_{22} & . & . \\ . & & . & 0 \\ l_{N1} & . & . & l_{NN} \end{array} \right) ; \;$ U$= \left( \begin{array}{cccc} u_{11} & u_{12} & . & u_{1N} \\ 0 & u_{22} & . & . \\ . & . & . & . \\ 0 & . & 0 & u_{NN} \end{array} \right)$

When this is done, write $A \cdot x=b$ as

$L \cdot \left( U \cdot x \right) = b$

such that we have the two simpler systems of equations

$L \cdot y = b$

and

$U \cdot x = y$

Since $L$ and $U$ are triangular these equations are easy to solve.

Forward substitution:
$\begin{eqnarray} y_{1} &=& \frac{1}{l_{11}} b_{1} \\ y_{i} &=& \frac{1}{l_{ii}} \left( b_{i} - \sum_{j=1}^{i-1} l_{ij} y_{j} \right) ; \;\;\; i=2, \dots ,N \end{eqnarray}$

Back substitution:
$\begin{eqnarray} x_{N} &=& \frac{1}{u_{NN}} y_{N} \\ x_{i} &=& \frac{1}{u_{ii}} \left( y_{i} - \sum_{j=i+1}^{N} u_{ij} x_{j} \right) ; \;\;\; i=N\!-\!1, \dots ,1 \end{eqnarray}$

But first: find the matrices $L$ and $U$!

Write $L \cdot U= A$ as

$\sum_{k=1}^{N} l_{ik} u_{kj} = a_{ij} ; \; i=1, \dots N; \; j=1, \dots N$

These are $N^{2}$ equations but $N^{2}+N$ unknowns $l_{ij}, u_{ij}$.

$\Longrightarrow$ Choose $l_{ii}=1 (i=1, \dots N)$. Also, due to the triangular structure of $L$ and $U$ we may write
$\begin{eqnarray} {\rm for} \;\;\; i \leq j &:& \; \sum_{k=1}^{i} l_{ik} u_{kj} = a_{ij} \; \\ {\rm for} \;\;\; i > j &:& \; \sum_{k=1}^{j} l_{ik} u_{kj} = a_{ij} \end{eqnarray}$

$\Longrightarrow$ Crout's procedure:

LU decomposition: For $j=1,2, \dots N$ compute
$\begin{eqnarray} u_{1j} &=& a_{1j} \\ u_{ij} &=& a_{ij} - \sum_{k=1}^{i-1} l_{ik} u_{kj} ; \;\;\; i=2, \dots ,j \\ l_{ij} &=& \frac{1}{u_{jj}} \left( a_{ij} - \sum_{k=1}^{j-1} l_{ik} u_{kj} \right) ; \;\;\; i=j+1, \dots ,N \end{eqnarray}$

Side result: Determinant $\left| A \right| = \left| L \right| \cdot \left| U \right| = u_{11} u_{22} \dots u_{NN}$

EXAMPLE: Let

$A=\left( \begin{array}{cc}1&2\\ \\3&4\end{array} \right)$

$\Longrightarrow$ (Crout)
$\begin{eqnarray} j=1, \; i=1: \;\;\; u_{11}&=&a_{11}=1 \\ j=1, \; i=2: \;\;\; l_{21}&=&\frac{1}{u_{11}} a_{21}=3 \\ && \\ j=2, \; i=1: \;\;\; u_{12}&=&a_{12}=2 \\ j=2, \; i=2: \;\;\; u_{22}&=&a_{22}-l_{21} u_{12}=-2 \end{eqnarray}$

and thus $\begin{eqnarray} \left( \begin{array}{cc}1&2\\ \\3&4\end{array} \right)&=&\underbrace{\left( \begin{array}{cc}1&0\\ \\3&1\end{array} \right)} \cdot \underbrace{\left( \begin{array}{cc}1&2\\ \\0&-2\end{array} \right)} \\ && \;\;\;L \;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;\; U \end{eqnarray}$
At each step $(j,i)$ the required elements $l_{ik}$, $u_{kj}$ are already available.

The vector $b$ has not been touched. Therefore we may use the factors $L$ and $U$ of a given matrix $A$ with different vectors $b$.

To find the Inverse $A^{-1}$: Solve $A \cdot x_{j} = e_{j}$ with the $N$ unit vectors $e_{j}\;$; combine the column vectors $x_{j}$ to find $A^{-1}$.

2.2.4 Recursion Method

Find solution $x$ if $A$ is tri-diagonal (maybe after application of Householder).

With
$\begin{eqnarray} A &\equiv& \left( \begin{array}{cccccc} \beta_{1}& \gamma_{1} & 0 & . & . & 0 \\ \alpha_{2} & \beta_{2} & \gamma_{2} & 0 & . & 0 \\ 0 & \alpha_{3} & \beta_{3} & \gamma_{3} & 0 & . \\ . & . & . & . & . & . \\ . & . & . & \alpha_{N-1} & \beta_{N-1} & \gamma_{N-1} \\ . & . & . & 0 & \alpha_{N} & \beta_{N} \end{array} \right) \end{eqnarray}$

$\begin{eqnarray} \beta_{1}x_{1} + \gamma_{1}x_{2}&=&b_{1} \\ \alpha_{i}x_{i-1}+\beta_{i}x_{i} + \gamma_{i}x_{i+1}&=&b_{i} ; \;\; i=2,\dots,N\!-\!1 \\ \alpha_{N}x_{N\!-\!1} + \beta_{N}x_{N}&=&b_{N} \end{eqnarray}$

Introducing auxiliary variables $g_{i}$ and $h_{i}$ by the recursive ansatz

$\begin{eqnarray} x_{i+1}&=&g_{i} x_{i} + h_{i} ; \; i=1, \dots , N\!-\!1 \end{eqnarray}$

we find the "downward recursion formulae"

$\textstyle \begin{eqnarray} g_{N-1}=\frac{-\alpha_{N}}{\beta_{N}} \; & ; & \; h_{N-1}=\frac{b_{N}}{\beta_{N}} \\ g_{i-1}=\frac{-\alpha_{i}}{\beta_{i}+\gamma_{i} g_{i}} \; & ; & \; h_{i-1}=\frac{b_{i}-\gamma_{i} h_{i}}{\beta_{i}+\gamma_{i} g_{i}} ; \;\; i=N-1, \dots,2 \end{eqnarray}$

Having arrived at $g_{1}$ and $h_{1}$ we insert the known values of $g_{i}, h_{i}$ in the "upward recursion formulae"

$\begin{eqnarray} x_{1}&=&\frac{b_{1}-\gamma_{1} h_{1}}{\beta_{1}+\gamma_{1} g_{1}} \\ x_{i+1}&=&g_{i}x_{i}+h_{i} ; \;\; i=1, \dots,N\!-\!1 \end{eqnarray}$

(The equation for the starting value $x_{1}$ follows from $\beta_{1} x_{1} + \gamma_{1}x_{2}=b_{1}$ and $x_{2}=g_{1}x_{1}+h_{1}$.)

EXAMPLE: In $A \cdot x = b$, let

$A \equiv \left( \begin{array}{cccc} \beta_{1}& \gamma_{1} & 0 & 0 \\ \alpha_{2} & \beta_{2} & \gamma_{2} & 0 \\ 0 & \alpha_{3} & \beta_{3} & \gamma_{3} \\ 0 & 0 & \alpha_{4} & \beta_{4} \end{array} \right) = \left( \begin{array}{cccc} 2 & 1 & 0 & 0 \\ 2 & 3 & 1 & 0 \\ 0 & 1 & 4 & 2 \\ 0 & 0 & 1 & 3 \end{array} \right) \;\; \;\; {\rm and} \;\;\; b = \left( \begin{array}{c} 1 \\2 \\3 \\4 \end{array} \right)$

Downward recursion: $g_{3}=-\alpha_{4}/\beta_{4}=-1/3$, $h_{3}= b_{4}/\beta_{4}=4/3$, and

$\begin{eqnarray} i=3: \;\; g_{2}= -3/10 \; &,& \;\; h_{2}= 1/10 \\ i=2: \;\; g_{1}= -20/27 \; &,& \;\; h_{1}= 19/27 \end{eqnarray}$

Upward recursion: $x_{1} = 8/34$, and
$\begin{eqnarray} i=1: \;\; x_{2} &=& 9/17 \\ i=2: \;\; x_{3} &=& -1/17 \\ i=3: \;\; x_{4} &=& 23/17 \end{eqnarray}$

vesely 2005-10-10

 < >