embryo

エンジニアの備忘録

LU分解を実装してみた

線形代数 アルゴリズム Javascript

今日は連立方程式の解法などに用いられるLU分解についてです。

LU分解とは

とある正方行列Aに対して以下のような等式が成り立つ場合を考えます。

 Ax = B

このとき A を L(下三角行列)とU(上三角行列)に分解すると三角化された行列の計算を2回行えば 解を得られることになるので計算が簡単になるというもの。

 A = LU

{ 
L = \begin{bmatrix} 
    1 & 0 & 0 & 0 \\\
    \vdots & 1 & 0 & 0 \\\ 
    \vdots & \vdots & 1 & 0 \\\ 
    l_{n1} & l_{n2} & \ldots & 1
\end{bmatrix}
}

{ 
U = \begin{bmatrix} 
    u_{11} & u_{12} & \ldots & u_{1n} \\\
    0 & u_{22} & \ldots & u_{2n} \\\ 
    0 & 0 & \ddots & \vdots \\\ 
    0 & 0 & 0 & u_{nn}
\end{bmatrix}
}

LU分解のメリット

では一体どう簡単になるかという話ですが

 Ax = B
 LUx = B
 Ly = B
 Ux = y

と変形すると結局三角行列の方程式を２つ計算すればxが求められることになります。例えば二元連立方程式の場合



\begin{bmatrix} 
    1 &  0 \\\ 
    l_{21} & 1 
\end{bmatrix}
\begin{bmatrix} 
    y_1 \\\ 
    y_2
\end{bmatrix}
=
\begin{bmatrix} 
    b_1 \\\ 
    b_2
\end{bmatrix}


 
\begin{bmatrix} 
    u_{11} & u_{12} \\\ 
    0 & u_{22}
\end{bmatrix}
\begin{bmatrix} 
    x_1 \\\ 
    x_2
\end{bmatrix}
=
\begin{bmatrix} 
    y_1 \\\ 
    y_2
\end{bmatrix}

となるため、自明なものから順に計算すれば簡単に求めることができます。

LU分解の過程について

 A = LU



\left(
    \begin{array}
      a_{11} & a_{12} & \ldots & a_{1n} \\
      a_{21} & a_{22} & \ldots & a_{2n} \\
      \vdots & \vdots & \ddots & \vdots \\
      a_{n1} & a_{n2} & \ldots & a_{nn}
    \end{array}
\right)

=

\left(
    \begin{array}
      l_{11} & l_{12} & \ldots & l_{1n} \\
      l_{21} & l_{22} & \ldots & l_{2n} \\
      \vdots & \vdots & \ddots & \vdots \\
      l_{n1} & l_{n2} & \ldots & l_{nn}
    \end{array}
\right)
\left(
    \begin{array}
      u_{11} & u_{12} & \ldots & u_{1n} \\
      u_{21} & u_{22} & \ldots & u_{2n} \\
      \vdots & \vdots & \ddots & \vdots \\
      u_{n1} & u_{n2} & \ldots & u_{nn}
    \end{array}
\right)




\left(
    \begin{array}{c|ccc}
      a_{11} & a_{12} & \ldots & a_{1n} \\
\hline
      a_{21} & a_{22} & \ldots & a_{2n} \\
      \vdots & \vdots & \ddots & \vdots \\
      a_{n1} & a_{n2} & \ldots & a_{nn}
    \end{array}
\right)
=
\left(
    \begin{array}{c|ccc}
      l_{11} & l_{12} & \ldots & l_{1n} \\
\hline
      l_{21} & l_{22} & \ldots & l_{2n} \\
      \vdots & \vdots & \ddots & \vdots \\
      l_{n1} & l_{n2} & \ldots & l_{nn}
    \end{array}
\right)
\left(
    \begin{array}{c|ccc}
      u_{11} & u_{12} & \ldots & u_{1n} \\
\hline
      u_{21} & u_{22} & \ldots & u_{2n} \\
      \vdots & \vdots & \ddots & \vdots \\
      u_{n1} & u_{n2} & \ldots & u_{nn}
    \end{array}
\right)


 
\begin{eqnarray}

\left( \begin{array}{c|c}
   a_{11} & a_1^T \\
   \hline
   a_1 & A_1 \\
\end{array}\right)

&=&

\left( \begin{array}{c|c}
   1 & 0 \\
   \hline
   l_{1} & L_1 \\
\end{array}\right)

\left( \begin{array}{c|c}
   u_{11} & u_1^T \\
   \hline
   0 & U_1 \\
\end{array}\right)
\\
&=&

\left( \begin{array}{c|c}
   u_{11} & u_1^T \\
   \hline
   l_1u_{11} & l_{1}u_1^T + L_1U_1 \\
\end{array}\right)

\end{eqnarray}

のように分解して計算すると、恒等式の関係より

 u_{11} = a_{11}
 u_1^T = a_1^T
 l_1 = a_1 / u_{11}
 A_1 = l_1u_1^T + L_1U_1

が成り立ちます。ここで A_1 - l_{1}u_1^T を新たに Aとすることで、同様の操作を繰り返すことができます。つまりn次正方行列を仮定するとn回繰り返すことで全てLU成分が分かることになります。

gist.github.com

計算量について

LU分解までの処理はO(n3)になるため、一度だけ方程式を解く場合はガウスの消去法などと変わりありません。 ガウスの消去法は右辺も処理する必要がありますが、LU分解は一度分解した結果を使いまわせるのがメリットです。 右辺の値を変更して何度もLUを再利用する場合に利点がありそうです。

まとめ

LU分解を自前で実装する方法とメリットについて説明しました。 次回は単体法について説明したいと思います。

参考