数値解析 第11回 (3) LU 分解法 Step 1 の手順
Step 1 の原理
$A=LU$ の成分を、左上から
$
\dps{\mat{c}{
\quad\overrightarrow{\qquad\qquad}\quad \\[-0.4em]
\overrightarrow{\qquad\qquad} \\[-0.7em]
\vdots \\
\overrightarrow{\qquad\qquad} \\
}}$
の順に見てゆくと、$i$ 行目では
$\ell_{i1}$, $\ell_{i2}$, $\cdots$, $\ell_{ii}$,
$u_{i,i+1}$, $u_{i,i+2}$, $\cdots$, $u_{in}$
の順に $L$ や $U$ の成分が求まります。
実際、未知の数を
赤字で表すと、
各 $a_{ij}$ の式にひとつずつ未知の数が入っていていることがわかります。
$\require{color}$
- 第 $1$ 行：
\begin{align}
a_{11} &= \textcolor{red}{\ell_{11}} \times 1 \\
a_{12} &= \ell_{11} \times \textcolor{red}{u_{12}} \\
a_{13} &= \ell_{11} \times \textcolor{red}{u_{13}} \\
\qquad \vdots \\
a_{1n} &= \ell_{11} \times \textcolor{red}{u_{1n}} \\
\end{align}
- 第 $2$ 行：
\begin{align}
a_{21} &= \textcolor{red}{\ell_{21}} \times 1 \\
a_{22} &= \ell_{21} \times u_{12} + \textcolor{red}{\ell_{22}} \times 1 \\
a_{23} &= \ell_{21} \times u_{13} + \ell_{22} \times \textcolor{red}{u_{23}} \\
&\quad \vdots \\
a_{2n} &= \ell_{21} \times u_{1n} + \ell_{22} \times \textcolor{red}{u_{2n}}\\
\end{align}
- 第 $i$ 行では：
\begin{align}
a_{i1} &= \textcolor{red}{\ell_{i1}} \times 1 \\
a_{i2} &= \ell_{i1} \times u_{12} + \textcolor{red}{\ell_{i2}} \times 1 \\
a_{i3} &= \ell_{i1} \times u_{13} + \ell_{i2} \times u_{23} + \textcolor{red}{\ell_{i3}} \times 1 \\
&\quad \vdots \\
a_{ii} &= \ell_{i1} \times u_{1i} + \,\cdots\cdots\cdot\cdot + \textcolor{red}{\ell_{ii}} \times 1 \\
a_{i,i+1} &= \ell_{i1} \times u_{1,i+1} + \cdots\cdots + \ell_{ii} \times \textcolor{red}{u_{i,i+1}} \\
&\quad \vdots \\
a_{in} &= \ell_{i1} \times u_{1n} + \cdots\cdots\cdot\cdot + \ell_{ii} \times \textcolor{red}{u_{in}}\\
\end{align}
といった具合です。$L$, $U$ の積は
$\require{color}
\dps{
\mat{c}{
\vdots \\[-1em]
\vdots \\
\fbox{$\ \ell_{i1}\ \ell_{i2}\ \cdots \ell_{ii}\ 0\ \cdots\ 0\ $} \\
\vdots \\[-1em]
\vdots \\
}
\mat{ccc}{
\cdots\cdots &
\fbox{$\begin{array}{c}u_{1j}\\[-1em] u_{2j}\\[-1em] \vdots\\[-1em] 1\\[-1em] 0\\[-1em] \vdots \\ 0\end{array}$} &
\cdots\cdots \\
}
}$
のように計算するので、「うしろの方が $0$ 」なのが効いています。
Step 1 の手順 ver.1
具体的には
Algorithm ( Step 1 ) ver.1
for $i=1$ to $n$
for $j=1$ to $i$
$\ell_{ij}=a_{ij}-\dps{\sum_{k=1}^{j-1} \ell_{ik}\,u_{kj}}$
for $j=i+1$ to $n$
$\dps{u_{ij}=\left(a_{ij}-\sum_{k=1}^{i-1} \ell_{ik}\,u_{kj}\right)\,/\,\ell_{ii}}$
となり、プログラムではこの $\sum$ を更に for 文で書き下します。
Step 1 の手順 ver.2
次のようにすれば $L$ や $U$ の成分を $A$ の成分に上書きしてメモリを節約することができます：
Algorithm ( Step 1 ) ver.2
for $i=1$ to $n-1$
$c=a_{ii}$
for $j=i+1$ to $n$
$a_{ij}=a_{ij}\,/\,c$
for $k=i+1$ to $n$
for $j=i+1$ to $n$
$a_{kj}=a_{kj}-a_{ki}\,a_{ij}$
このアルゴリズム終了時には
- $i \geqq j$ のとき $a_{ij}$ のメモリには $\ell_{ij}$ が
- $i \lt j$ のとき $a_{ij}$ のメモリには $u_{ij}$ が
それぞれ上書きされています。
※ $a_{ij}$ は、$\ell_{ij}$ あるいは $u_{ij}$ を求めた時点で不要になるので上書きしても構わない、という訳です。
ただし元の $A$ が残っていませんので検算ができません。上級者向けのアルゴリズムです。