『Python数値計算ノート』ではアフィリエイトプログラムを利用して商品を紹介しています。

ガウス・ジョルダンの消去法(掃き出し法)

≪【前の記事】基本行列と基本変形

ガウス・ジョルダンの消去法

以前の記事で連立方程式 $A\boldsymbol{x}=\boldsymbol{b}$ を $\boldsymbol{x}=A^{-1}\boldsymbol{b}$ の形で解く方法を解説しました。しかし、連立方程式を解くために明示的に逆行列を求める必要はありません。今回はガウス・ジョルダンの消去法(掃き出し法)とよばれる手法を解説します。

拡大係数行列

例として、
 \[A=\begin{bmatrix}4&1\\2&-1\end{bmatrix},\quad \boldsymbol{b}=\begin{bmatrix}9\\3\end{bmatrix}\tag{1}\]
とおいて、
 \[\begin{bmatrix}4&1\\2&-1\end{bmatrix}\begin{bmatrix}x\\y\end{bmatrix}=\begin{bmatrix}9\\3\end{bmatrix}\tag{2}\]
という連立方程式を考えます。最初に左辺の係数行列と右辺のベクトルをまとめて
 \[[A\ \ \boldsymbol{b}]=\begin{bmatrix}4&1&9\\2&-1&3\end{bmatrix}\tag{3}\]
という行列をつくります。これを拡大係数行列または拡大行列とよびます。ガウス・ジョルダンの消去法は、一般的な連立方程式を解く手順、すなわち、ある方程式に別の方程式の定数倍を足して (または引いて) 変数を消去していく手法と本質的に同じものです。つまり拡大行列 (3) に対して

・ある行に別の行の定数倍を加える
・ある行全体を定数倍する
・ある行と別の行を交換する
 
を繰り返して最終的には
 \[\begin{bmatrix}1&0&\alpha\\0&1&\beta\end{bmatrix}\tag{4}\]
という形を作り出します。行列形式の方程式で表すと
 \[\begin{bmatrix}1&0\\0&1\end{bmatrix}\begin{bmatrix}x\\y\end{bmatrix}=\begin{bmatrix}\alpha\\\beta\end{bmatrix}\tag{5}\]
となります。すなわち
 \[\begin{bmatrix}x\\y\end{bmatrix}=\begin{bmatrix}\alpha\\\beta\end{bmatrix}\tag{6}\]
なので、拡大行列の右端の列 $\alpha$ と $\beta$ は解となります。それでは実際に拡大行列
 \[\begin{bmatrix}4&1&9\\2&-1&3\end{bmatrix}\tag{7}\]
を操作して連立方程式を解いてみましょう。

【操作①】$1$行目に $2$行目を加える
 \[\begin{bmatrix}6&0&12\\2&-1&3\end{bmatrix}\tag{8}\]
【操作②】$2$行目から $1$ 行目の$\cfrac{2}{6}$ 倍を引く
 \[\begin{bmatrix}6&0&12\\2&-1&-1\end{bmatrix}\tag{9}\]
【操作③】$1$行目を $\cfrac{1}{6}$ 倍
 \[\begin{bmatrix}1&0&2\\0&-1&-1\end{bmatrix}\tag{10}\]
【操作④】$2$行目を $(-1)$ 倍
 \[\begin{bmatrix}1&0&2\\0&1&1\end{bmatrix}\tag{11}\]
左側の正方行列の部分は単位行列となりました。
 
右端のベクトル $\begin{bmatrix}x\\y\end{bmatrix}=\begin{bmatrix}2\\1\end{bmatrix}$ が方程式 (2) の解となります。

以上の手順は前回記事で学んだ 基本行列 によって再現できます。
たとえば、1 行目に 2 行目を加えることは
 \[E_{1,2,1}=\begin{bmatrix}1&1\\0&1\end{bmatrix}\tag{12}\]
を乗じるに等しい操作です。実際、
 \[\begin{bmatrix}1&1\\0&1\end{bmatrix}\begin{bmatrix}4&1&9\\2&-1&3\end{bmatrix}=\begin{bmatrix}6&0&12\\2&-1&3\end{bmatrix}\tag{13}\]
となって、これは ① と等しい操作になっています。同様に方程式 (9) の右辺の拡大行列に左側から
 \[E_{2,1,-2/6}=\begin{bmatrix}1&0\\-\frac{2}{6}&1\end{bmatrix}\tag{14}\]
を乗じると
 \[\begin{bmatrix}1&0\\-\frac{2}{6}&1\end{bmatrix}\begin{bmatrix}6&0&12\\2&-1&3\end{bmatrix}=\begin{bmatrix}6&0&12\\0&-1&-1\end{bmatrix}\tag{15}\]
となって、これは ② と同じ操作です。さらに
 \[C_{1,1/6}=\begin{bmatrix}\frac{1}{6}&0\\0&1\end{bmatrix}\tag{16}\]
を左から掛ければ
 \[\begin{bmatrix}\frac{1}{6}&0\\0&1\end{bmatrix}\begin{bmatrix}6&0&12\\0&-1&-1\end{bmatrix}\begin{bmatrix}1&0&2\\0&-1&-1\end{bmatrix}\tag{17}\]
となって、操作 ③ と同じ結果を得ます。最後に
 \[C_{2,-1}=\begin{bmatrix}1&0\\0&-1\end{bmatrix}\tag{18}\]
を掛ければ $2$ 行目が $-1$ 倍されます (操作④):
 \[\begin{bmatrix}1&0\\0&-1\end{bmatrix}\begin{bmatrix}1&0&2\\0&-1&-1\end{bmatrix}=\begin{bmatrix}1&0&2\\0&1&1\end{bmatrix}\tag{19}\]
以上より、操作 ① ~ ④ は拡大係数行列 $[A\ \ \boldsymbol{b}]$ の左側に
 \[C_{2,-1}\ C_{1,1/6}\ E_{2,1,-2/6}\ E_{1,2,1}=\begin{bmatrix}\cfrac{1}{6}&\cfrac{1}{6}\\\cfrac{1}{3}&-\cfrac{2}{3}\end{bmatrix}\tag{20}\]
を掛けることで再現できることがわかります:

ガウス・ジョルダン消去法による逆行列の計算

ガウス・ジョルダンの消去法を使って逆行列を求めることもできます。
定義より、$A$ と逆行列 $A^{-1}$ の積は単位行列となります。
 \[AA^{-1}=I\tag{21}\]
$A^{-1}$ が求めるべき未知行列です。いま、$A$ が $3\times 3$ の正方行列であるとして、$A^{-1}$ の各列を $\boldsymbol{x}_1,\ \boldsymbol{x}_2,\ \boldsymbol{x}_3$ とおくと、
 \[A\ [\boldsymbol{x}_1\ \ \boldsymbol{x}_2\ \ \boldsymbol{x}_3]=[\boldsymbol{e}_1\ \ \boldsymbol{e}_2\ \ \boldsymbol{e}_3]\tag{22}\]
と表せます。$\boldsymbol{e}_1\ \boldsymbol{e}_2\ \boldsymbol{e}_3$ は単位行列 $I$ の各列を成分にもつ単位ベクトルです:
 \[\boldsymbol{e}_1=\begin{bmatrix}1\\ 0\\ 0\end{bmatrix},\quad \boldsymbol{e}_2=\begin{bmatrix}0\\ 1\\ 0\end{bmatrix},\quad \boldsymbol{e}_3=\begin{bmatrix}0\\ 0\\ 1\end{bmatrix}\tag{23}\]
(22) は 3 個の連立方程式
 \[\begin{align*}&A\boldsymbol{x}_1=\boldsymbol{e}_1\\[6pt]&A\boldsymbol{x}_2=\boldsymbol{e}_2\\[6pt]&A\boldsymbol{x}_3=\boldsymbol{e}_3\end{align*}\tag{24}\]
を表しています。拡大係数行列
 \[[A\ \ I]=[A\ \ \boldsymbol{e}_1\ \ \boldsymbol{e}_2\ \ \boldsymbol{e}_3]\tag{25}\]
に操作をほどこして、最終的に $[I\ \ B]$ の形を作り出せば (ガウス・ジョルダンの消去法)、
 \[I\ [\boldsymbol{x}_1\ \ \boldsymbol{x}_2\ \ \boldsymbol{x}_3]=B\tag{26}\]
となるので、$B$ が求める逆行列 $A^{-1}$ となります。

ガウス・ジョルダン法を使って $3\times 3$ の正方行列
 \[A=\begin{bmatrix}1&1&-2\\2&1&1\\1&-1&2\end{bmatrix}\tag{27}\]
の逆行列を求めてみます。拡大行列 $[A\ \ I]$ をつくると
 \[[A\ \ I]=\begin{bmatrix}1&1&-2&1&0&0\\2&1&1&0&1&0\\1&-1&2&0&0&1\end{bmatrix}\tag{28}\]
【操作①】$1$ 行目に $3$ 行目を加える
 \[\begin{bmatrix}2&0&0&1&0&1\\2&1&1&0&1&0\\1&-1&2&0&0&1\end{bmatrix}\tag{29}\]
【操作②】$3$ 行目に $2$ 行目を加える
 \[\begin{bmatrix}2&0&0&1&0&1\\2&1&1&0&1&0\\3&0&3&0&1&1\end{bmatrix}\tag{30}\]
【操作③】$2$ 行目から $1$ 行目を引く
 \[\begin{bmatrix}2&0&0&1&0&1\\0&1&1&-1&1&-1\\3&0&3&0&1&1\end{bmatrix}\tag{31}\]
【操作④】$3$ 行目から $1$ 行目の $\cfrac{3}{2}$ 倍を引く
 \[\begin{bmatrix}2&0&0&1&0&1\\0&1&1&-1&1&-1\\0&0&3&-\frac{3}{2}&1&-\frac{1}{2}\end{bmatrix}\tag{32}\]
【操作⑤】$2$ 行目から $3$ 行目の $\cfrac{1}{3}$ 倍を引く
 \[\begin{bmatrix}2&0&0&1&0&1\\0&1&0&-\frac{1}{2}&\frac{2}{3}&-\frac{5}{6}\\0&0&3&-\frac{3}{2}&1&-\frac{1}{2}\end{bmatrix}\tag{33}\]
【操作⑥】$1$ 行目を $1/2$ 倍する
 \[\begin{bmatrix}1&0&0&\frac{1}{2}&0&\frac{1}{2}\\0&1&0&-\frac{1}{2}&\frac{2}{3}&-\frac{5}{6}\\0&0&3&-\frac{3}{2}&1&-\frac{1}{2}\end{bmatrix}\tag{34}\]
【操作⑦】$3$ 行目を $1/3$ 倍する
 \[\begin{bmatrix}1&0&0&\frac{1}{2}&0&\frac{1}{2}\\0&1&0&-\frac{1}{2}&\frac{2}{3}&-\frac{5}{6}\\0&0&1&-\frac{1}{2}&\frac{1}{3}&-\frac{1}{6}\end{bmatrix}\tag{35}\]
したがって、$A$ の逆行列は
 \[A^{-1}=\begin{bmatrix}\frac{1}{2}&0&\frac{1}{2}\\-\frac{1}{2}&\frac{2}{3}&-\frac{5}{6}\\-\frac{1}{2}&\frac{1}{3}&-\frac{1}{6}\end{bmatrix}\tag{36}\]
となります。

操作 ① ~ ⑦ に対応する行列を並べてみます:
 \[\begin{align*}E_{1,3,1}&=\begin{bmatrix}1&0&1\\0&1&0\\0&0&1\end{bmatrix}\\[6pt]E_{3,2,1}&=\begin{bmatrix}1&0&0\\0&1&0\\0&1&1\end{bmatrix}\\[6pt]E_{2,1,-1}&=\begin{bmatrix}1&0&0\\-1&1&0\\0&0&1\end{bmatrix}\\[6pt]E_{3,1,-3/2}&=\begin{bmatrix}1&0&0\\0&1&0\\-\frac{3}{2}&0&1\end{bmatrix}\\[6pt]E_{2,3,-1/3}&=\begin{bmatrix}1&0&0\\0&1&-\frac{1}{3}\\0&1&0\end{bmatrix}\\[6pt]C_{1,1/2}&=\begin{bmatrix}\frac{1}{2}&0&0\\0&1&0\\0&0&1\end{bmatrix}\\[6pt]C_{3,1/3}&=\begin{bmatrix}1&0&0\\0&1&0\\0&0&\frac{1}{3}\end{bmatrix}\end{align*}\]
行列積をとると、逆行列 $A^{-1}$ が得られます。
 \[\begin{align*}&C_{3,1/3}\ C_{1,1/2}\ E_{2,3,-1/3}\ E_{3,1,-3/2}\ E_{2,1,-1}\ E_{3,2,1}\ E_{1,3,1}\\[6pt]&=\begin{bmatrix}\frac{1}{2}&0&\frac{1}{2}\\-\frac{1}{2}&\frac{2}{3}&-\frac{5}{6}\\-\frac{1}{2}&\frac{1}{3}&-\frac{1}{6}\end{bmatrix}\end{align*}\]
計算が面倒であれば、Python で確認しましょう。前回記事で作成した matrix_c() と matrix_e() も再掲します。

# python_gauss–jordan_elimination

# In[1]

import numpy as np

# i行を定数倍する関数
def matrix_c(n, i=1, k=1):
    arr = np.eye(n)
    arr[i-1, i-1] = k
    return arr

# i行にj行の定数倍を加える関数
def matrix_e(n, i=1, j=1, k=1):
    arr = np.eye(n)
    arr[i-1, j-1] = k
    return arr

np.set_printoptions(precision=5)

# 操作①~⑦に対応する基本行列
m1 = matrix_e(3, 1, 3, 1)
m2 = matrix_e(3, 3, 2, 1)
m3 = matrix_e(3, 2, 1, -1)
m4 = matrix_e(3, 3, 1, -3/2)
m5 = matrix_e(3, 2, 3, -1/3)
m6 = matrix_c(3, 1, 1/2)
m7 = matrix_c(3, 3, 1/3)

# 行列積を計算
inv_A = np.linalg.multi_dot([m7, m6, m5, m4, m3, m2, m1])

print(inv_A)
# [[ 0.5      0.       0.5    ]
#  [-0.5      0.66667 -0.83333]
#  [-0.5      0.33333 -0.16667]]

≫【次の記事】下三角行列と上三角行列
≫ Pythonで学ぶ線形代数トップページ

コメント

  1. HNaito より:

    下記は誤植と思われますので、ご確認ください。

    (14)式で、E_2,1,1 → E_2,1,-2/6
    (19)式の下の文章で、係数行列 → 拡大係数行列
    (20)式で、E_2,1,1 → E_2,1,-2/6
    (23)式で、e1=[1 1 0] → [1 0 0]
    (25)式の次の式番号で、(25) → (26)

    以降の式番号は 1 ずつ増えるとことになりますが、(33)式以降はそのままでOK
    (36)式のa_32要素で、1 → 1/3
    E_1,3,1の基本行列で、[1 0 0] →[1 0 1]
               [0 1 0] [0 1 0]
               [1 0 1] [0 0 1]
    その4行下の基本行列名で、E_2,3,-1 → E_2,3,-1/3
    さらに4行下の基本行列名で、E_2,3,-1 → E_2,3,-1/3