Solve $ A \cdot x = b$ exactly:
$
\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.
Advantage of LU decomposition:
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}$.
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}
$
the system of equations reads
$
\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