3 ガウス・ジョルダン法

逆行列が不要であれば、ガウス・ジョルダン法よりも、後で述べるLU分解の法 が計算速度は速い。しかし、教育的効果を考えると、両方の方法を知っておく のは良いことです。

3.1 基本的な考え方

ガウス・ジョルダン(Gauss-Jordan)法というのは、連立方程式 (4)を次にように変形させて、解く方法です。

$$\begin{aligned}\begin{pmatrix}1 & 0 & 0 & \cdots & 0 \\ 0 & 1 &... ... \\ b_3^\prime \\ \vdots\\ b_N^\prime \end{pmatrix} \end{aligned}$$

この式から明らかに、求める解 $x_{i}=b_{i}^\prime$となります。これをどう やって求めるか?。コンピューターで実際に計算する前に、人力でガウス・ジョ ルダン法で計算してみましょう。例として、

$$\left\{ \begin{aligned}3x+2y+z&=10\\ x+y+z&=6\\ x+2y+2z&=11 \end{aligned} \right.$$

をガウス・ジョルダン法で計算してみましょう。

まずは、1行目の$x$の係数を1に、2と3行目のそれは0にします。そのために、 1行目は$x$の係数の値で割ります。2行目と3行目は、1行目に適当な係数を書 けて引きます。次の通りです。

$$\left\{ \begin{aligned}x+\frac{2}{3}y+z&=\frac{10}{3}\\ x+y+z&=... ...\ 0+\frac{4}{3}y+\frac{5}{3}z&=\frac{22}{3} \end{aligned} \right.$$

つぎに、2行目のyの係数を1にして、1と3行目のそれを0にします。

$$\left\{ \begin{aligned}x+\frac{2}{3}y+\frac{1}{3}z&=\frac{10}{3... ...n{aligned}x+0-z&=-2\\ 0+y+2z&=8\\ 0+0-z&=-3 \end{aligned} \right.$$

同じことを$z$についても繰り返します。すると、

$$\left\{ \begin{aligned}x+0-z&=-2\\ 0+y+2z&=8\\ 0+0+z&=3 \end{al... ...egin{aligned}x+0+0&=1\\ 0+y+0&=2\\ 0+0+z&=3 \end{aligned} \right.$$

となります。従って、連立方程式(10)の解は、

$$\left\{ \begin{aligned}x=1 \\ y=2 \\ z=3 \end{aligned} \right.$$

となります。これがガウス・ジョルダン法です。

ガウス・ジョルダン法のイメージが湧いたと思いますので、もう少し数学的に その内容を説明します²。 そのために、次の線形行列方程式を考えます。ここでは、紙面の関係で係数行 列が $4 \times 4$について述べますが、一般的に$N$への拡張は容易です。

$$\begin{bmatrix}a_{11} & a_{12} & a_{13} & a_{14}\\ a_{21} & a_{22... ...& 0\\ 0 & 1 & 0 & 0\\ 0 & 0 & 1 & 0\\ 0 & 0 & 0 & 1 \end{pmatrix} \end{bmatrix}$$ (15)

ここで、$\sqcup$は列拡大、つまりこの両側の括弧を取り去って行列をくっつ けて幅を広くすること意味します。この式から、容易に

$$\boldsymbol{A}\boldsymbol{x} =\boldsymbol{b}$$ (16)

$$\boldsymbol{A}\boldsymbol{Y} =\boldsymbol{E}$$ (17)

が分かります。もちろん、 $\boldsymbol{E}$は単位行列です。要するに行列方程式 (15)を解くということは、この2つの連立方 程式を解くことに他なりません。 $\boldsymbol{x}$はもとの方程式 (1)の解となり、 $\boldsymbol{Y}$は $\boldsymbol{A}$の逆行列になっ ています。

式(15)について、以下のことが容易に分か ります。

• $\boldsymbol{A}$の任意の2行を入れ替えて、それに対応する $\boldsymbol{b}$と $\boldsymbol{E}$を入れ替えた場合、 $\boldsymbol{x}$や $\boldsymbol{Y}$の値と順序は変わらない。た だし、 $\boldsymbol{E}$はもはや単位行列ではない。これは連立方程式の順序を入れ替 えて書いていることに相当している。
• $\boldsymbol{A}$の任意の行を、その行と別の行との線形結合に置き換え、同 時に対応する $\boldsymbol{b}$と $\boldsymbol{E}$も同様に置き換える。この場合、 $\boldsymbol{x}$ と $\boldsymbol{Y}$の値と順序は変わらない。当然、この場合も $\boldsymbol{E}$はもはや単 位行列ではない。
• $\boldsymbol{A}$の任意の2列を入れ替えて、それに対応する $\boldsymbol{x}$と $\boldsymbol{Y}$の行を入れ替えれば、 $\boldsymbol{b}$と $\boldsymbol{E}$の順序は入れ替える必要 は無い。この場合、解の行の順序が変わるので、最後に元に戻す操作が必要 になる。

この3つの操作を組み合わせて、係数行列 $\boldsymbol{A}$を単位行列に変換するのが ガウス・ジョルダン法です。 $\boldsymbol{A}$が単位行列に変換されれば、右辺に $\boldsymbol{x}$と $\boldsymbol{Y}$が表れます。したがって、解と逆行列が求められたことに なります。もし、逆行列が不要であれば $\boldsymbol{x}$だけ計算し、逆行列のみ必要 であれば $\boldsymbol{Y}$のみ計算することになります。

3.2 ピボット選択

先に示した、ガウス・ジョルダン法の3つの基本操作のうち、2番目しか使わな い方法を「ピボット選択なしのガウス・ジョルダン法」といいます。最初、人 力で連立方程式を解いた方法です。この方法の明らかにまずい点は、もし1に したい対角要素がゼロの場合、計算ができなくなってしまいます。この割る要 素をピボット(pivot)と言います。ゼロでないにしても、そのピボットが非常 に小さい値の場合、丸め誤差が大きくなります³。このようなことから、普通は ピボット選択なしのガウス・ジョルダン法というものは考えられません。

この問題を避けるためにどうするかというと、ピボット選択という方法を使い ます。方法は簡単で、先に示した3つの基本操作のうち、1番目と3番目を使っ て、対角に素性の良い要素をもってきます。1番目の操作のみを用いて行を入 れ替える方法を、部分ピボット選択(partial pivoting)といいます。1番目の 操作と3番目の操作を使って、行と列を入れ替えるのを完全ピボット選択(full pivoting)といいます。すでにある程度出来上がっている単位行列を壊したく ないので、ピボットの選択は操作している行のしたの行から選ばなくてはなり ません。

部分ピボット選択の方が明らかに簡単です。解の行列を入れ替える必要が無い からです。その場合、行の入れ替えしかしないので、ピボットはその列から選 ばなくてはなりません。完全ピボット選択の方が選べる要素が多いのですが、 最終的な解の精度はあまり変わらないようです。したがって、ここではプログ ラムの簡単な部分ピボット選択で計算しましょう。

次に考えなくてはならないのは、ピボットを選択する基準です。簡単に言えば、 大きな要素選択すれば大体よいです。しかし、ある行を100万倍して、それに 対応する右辺の行も100万倍することもできますので、ただ大きいというだけ では問題がありそうです。どうするかと言うと、各方程式の最大係数を1に規 格化して、最大のものをピボットに選ぶことが行われています。この方法を陰 的ピボット選択(implicit pivoting)と呼ばれています。

これで、ピボットの問題も片付いたので、フローチャートを書いてみます。

ホームページ: Yamamoto's laboratory
著者: 山本昌志