サンプルファイルとして “GE.xlsm” (Excel2007以上、zip形式で圧縮)を公開しています。 必要に応じ、ダウンロードしてご確認ください。

(「ファイルが破損しているため開くことができません」 コメントが出る場合は こちらで対応願います)


1.ガウスの消去法とは?

ガウスの消去法とは、 逆行列を用いないn元連立一次方程式(未知数xi)の解法です。


n元連立方程式

上式は行列を用いて次のように書き換えられます。


行列による表現

この行列Aの下半分の非対角成分を“0”に変換する、つまり上三角行列化すると、 xnが簡単に求まります。


行列の上三角行列化

このような上三角行列化による方法を前進消去と呼びます。

xnが求まれば、あとは下から順に追いかけて、 xn-1・・・x1と定まります。

この方法を後退代入と呼びます。


このような流れで連立方程式の解を導き出すことができます。


2.計算方法

2.1.前進消去

n元連立方程式において、ある行の係数行列成分aijとベクトル成分ciに 定数を掛けても方程式の内容は変わらない、という性質(以下に例を示します)に着目します。


方程式の係数倍

ここで、計算する前の準備を行います。


<準備1>

ベクトル成分を係数行列のn+1列成分として扱います。


ベクトル成分

<準備2>

添字を次のように決めておきます。

i:行番号 :1≦i≦n
j:列番号 :1≦j≦n+1
k:計算繰り返し番号 :0≦k≦n-1

いよいよ、行列Aを上三角化行列化する手順について見ていきます。


<STEP I>

1行目は何もせず、2~n行目の第一列セルが0になるよう、 2行目以降の行列成分とベクトル成分を下記に従い係数変換します。

なお、係数の上添字(0)は変換前の行列/ベクトル成分を表します (aij=aij(0))。


係数変換

この計算により、(1-2)式は次のようになります。


係数変換

また、変換式を一般化すると、次のよう表せます。


係数変換の一般式

<STEP II>

<STEP I>で変換完了した行列/ベクトル成分に対し、 3~n行目の第二列セルが0になるよう、<STEP I>同様係数変換します。


係数変換

この計算により、(1-2)式は次のようになります。


係数変換

また、変換式を一般化すると、次のよう表せます。


係数変換の一般式

<STEP III>

<STEP 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が簡単に求まります。


行列の上三角行列化
未知数n番目の解

次にxn-1が求まります。


未知数n-1番目の解

さらにxn-2が求まります。


未知数n-2番目の解

これを繰り返し行うと、次の一般式が得られます。


後退代入の一般式

以上、(2.3-1)(2.3-2)式により、未知数xi (i=1、2、…、n)が求まります。


3.マクロ(VBA)


ここでは、ガウスの消去法をExcel VBAを用いて行ってみます。


3.1.サンプルファイルのダウンロード


必要に応じてサンプルファイル “GE.xlsm” (Excel2007以上、zip形式で圧縮)をダウンロードしてご確認ください。

(「ファイルが破損しているため開くことができません」 コメントが出る場合は こちらで対応願います)


GE.xlsmを開くと下図のようになります。


図3.1-1 GE.xlsm画面
図3.1-1 GE.xlsm画面

サンプルファイルでは、n=15までしか計算できませんが、 マクロの制限を取り外せばいくらでも対応できます (とはいえ、マシンの能力に依存しますが...)。

計算するには、桃色セルに係数項である行列成分を、緑色セルに定数項であるベクトル成分を入力し、 “計算ボタン”をクリックすると、緑色セルの隣の列に答えが表示されます。


図3.1-2 GE.xlsm計算例
図3.1-2 GE.xlsm計算例

3.2.マクロの流れ

GE.xlsmで設定しているマクロのフローチャートを以下に示します。


図3.2-1 マクロ・フローチャート
図3.2-1 マクロ・フローチャート

3.3.マクロ構成

マクロの構成は、以下の2つで構成されています。


(1)計算ボタンクリック時に実行するモジュール

(2)ガウスの消去法計算モジュール


図3.3-1 計算ボタン実行モジュール
図3.3-1 計算ボタン実行モジュール

図3.3-2 ガウスの消去法計算モジュール
図3.3-2 ガウスの消去法計算モジュール

3.3.1.計算ボタン実行モジュール

ここでは主に以下の内容が記述されています。


(1)行列計算元数制限処理(n=2~15)

(2)行列/ベクトル成分(セル値)取得処理

(3)ガウスの消去法計算モジュールの呼び出し

(4)計算結果の出力


特に(1)の計算元数制限を変更する場合は、こちらの部分を修正してください。


図3.3.1-1 計算元数制限
図3.3.1-1 計算元数制限

3.3.2.ガウスの消去法モジュール

ここでは、ガウスの消去法を実行するための内容が記述されています。


(1)前進消去

(2)ピボッティング

(3)後退代入


行列/ベクトル成分(aij、ci)を引数として受け取り、 戻り値として未知数(xi)を返します。

他のマクロへガウスの消去法マクロを移設する場合は、下図にあわせて引数設定してください。


Sub GAUSSIAN_ELIMINATION(errFlag,N,a,x)
errFlag :エラー処理用のフラグ、True時はエラーとして認識
N :元数(n元正方行列)
a(配列) :行列成分(列番号1~n)/ベクトル成分(列番号n+1)
x(配列) :未知数解

また、ピボット処理の閾値は、下図の部分で設定しています。お好みにあわせて修正してください。

図3.3.2-1 ピボッティング閾値設定
図3.3.2-1 ピボッティング閾値設定

コメント


役に立った 役に立たなかった

広告

広告募集中
サイズ:150×150

広告募集中
サイズ:150×150

広告に関するお問い合せ

イベント・学会・展示会


学べる探せる設計技術-沐エンジニアリング

Ads by Google

Ads by Google

Ads by Google


技術計算製作所のソフトウエア購入は


Vectorソフトライブラリ/ビジネス

Amazon

技術計算製作所

技術計算製作所

画像をクリックするとpdfが開きます。