I.Excel技術計算サンプルファイルとして“GE.xlsm”(Excel2007以上、zip形式で圧縮)を公開しています。必要に応じ、ダウンロードしてご確認ください。逆行列を用いないn元連立一次方程式(未知数xi)の解法です。上式は行列を用いて次のように書き換えられます。この行列Aの下半分の非対角成分を“0”に変換する、つまり上三角行列化すると、xnが簡単に求まります。このような上三角行列化による方法を前進消去と呼びます。xnが求まれば、あとは下から順に追いかけて、xn-1・・・x1と定まります。この方法を後退代入と呼びます。このような流れで連立方程式の解を導き出すことができます。2.計算方法2.1.前進消去n元連立方程式において、ある行の係数行列成分aijとベクトル成分ciに定数を掛けても方程式の内容は変わらない、という性質(以下に例を示します)に着目します。ここで、計算する前の準備を行います。<準備1>ベクトル成分を係数行列のn+1列成分として扱います。<準備2>添字を次のように決めておきます。i:行番号:1≦i≦nj:列番号:1≦j≦n+1k:計算繰り返し番号:0≦k≦n-1いよいよ、行列Aを上三角化行列化する手順について見ていきます。I>1行目は何もせず、2~n行目の第一列セルが0になるよう、2行目以降の行列成分とベクトル成分を下記に従い係数変換します。係数の上添字(0)は変換前の行列/ベクトル成分を表します(aij=aij(0))。この計算により、(1-2)式は次のようになります。変換式を一般化すると、次のよう表せます。I>同様係数変換します。この計算により、(1-2)式は次のようになります。変換式を一般化すると、次のよう表せます。II>を(1-3)式が得られるまで繰り返し行います(n-1回)。この時、係数変換の一般式は次のように表せます。この式が、前進消去の一般式になります。以上で、前進消去処理は完了です。2.2.ピボッティングところで(2.1-5)式を見ると、のとき、“分母=0”となるため、前進消去では計算できないことがわかります。この絶対値が非常に小さい場合、数値計算では誤差を大きくする原因になってしまいます。この問題は、次の方法で避けられる場合があります(常に避けられるとは限りません)。連立方程式の順番を入れ替えても未知数の順番に影響を与えない、という性質に着目します。そこで、のとき、k+1≦i≦n行の範囲の中で、k列の成分が絶対値となる行を選び出し、その行mとk行の成分を丸ごと入れ替えた後、前進消去を行います。この入れ替え処理をピボッティングと言います。仮にk=1の処理でa22=0になった場合、a32~an2の中で絶対値が最大になるものを選び出します(例えばa42が最大とします)。このとき、2行目の行列/ベクトル成分と4行目の行列/ベクトル成分をそっくり入れ替えることで、(2.1-5)式における0で割る操作を避けることが可能になります。2.3.後退代入これまでの作業で(1-3)式のように行列Aを上三角行列化できました。することにより、xnが簡単に求まります。次にxn-1が求まります。xn-2が求まります。これを繰り返し行うと、次の一般式が得られます。以上、(2.3-1)(2.3-2)式により、未知数xi(i=1、2、...、n)が求まります。VBAを用いて行ってみます。“GE.xlsm”(Excel2007以上、zip形式で圧縮)をダウンロードしてご確認ください。(「ファイルが破損しているため開くことができません」コメントが出る場合はこちらで対応願います)GE.xlsmを開くと下図のようになります。GE.xlsm画面サンプルファイルでは、n=15までしか計算できませんが、マクロの制限を取り外せばいくらでも対応できます(とはいえ、マシンの能力に依存しますが...)。“計算ボタン”をクリックすると、緑色セルの隣の列に答えが表示されます。GE.xlsm計算例3.2.マクロの流れGE.xlsmで設定しているマクロのフローチャートを以下に示します。マクロ・フローチャート3.3.マクロ構成マクロの構成は、以下の2つで構成されています。ガウスの消去法計算モジュール3.3.1.計算ボタン実行モジュールここでは主に以下の内容が記述されています。(1)の計算元数制限を変更する場合は、こちらの部分を修正してください。計算元数制限3.3.2.ガウスの消去法モジュールここでは、ガウスの消去法を実行するための内容が記述されています。戻り値として未知数(xi)を返します。他のマクロへガウスの消去法マクロを移設する場合は、下図にあわせて引数設定してください。GAUSSIAN_ELIMINATION(errFlag,N,a,x)errFlag:エラー処理用のフラグ、True時はエラーとして認識N:元数(n元正方行列)a(配列):行列成分(列番号1~n)/ベクトル成分(列番号n+1)x(配列):未知数解また、ピボット処理の閾値は、下図の部分で設定しています。お好みにあわせて修正してください。クリックするとpdfが開きます。