連立1次方程式(Linear Equations)は,次のような形をしている.
式(7)は行列とベクトルで書くと,式がすっきりして
考えやすくなる.書き直すと,
である.それぞれの行列とベクトルは,
を表す.
通常,連立1次方程式(7)は
と書き表せる.このようにすると,見通しがかなり良くなる.
ガウス・ジョルダン(Gauss-Jordan)法というのは,連立方程式
(
10)を次にように変形させて,解く方法である.
この式から明らかに,求める解
![$ x_{i}=b_{i}^\prime$](img49.png)
となる.これをどうやって求めるか?.
コンピューターで実際に計算する前に,人力でガウス・ジョルダン法で計算してみる.例
として,
をガウス・ジョルダン法で解を求める.
解くべき,方程式
2行
![$ -2\times$](img52.png)
1行
3行
![$ -2\times$](img54.png)
1行
![$ -\frac{1}{2}\times$](img56.png)
2行
1行
![$ -2\times$](img58.png)
2行
3行
![$ +2\times$](img60.png)
2行
![$ -\frac{1}{2}\times$](img62.png)
3行
![$ 2行-\frac{3}{2}\times$](img64.png)
3行
これで,ガウス・ジョルダン法による対角化の作業は完了である.これか
ら,
![$ (x_1,x_2,x_3)=(-1,0,1)$](img66.png)
と分かる.
ガウス・ジョルダンを使って,逆行列が求められる.以下のようにする.解くべき,方程式
とする.
2行
1行
3行
![$ -2\times$](img70.png)
1行
![$ -\frac{1}{2}\times$](img72.png)
2行
1行
![$ -2\times$](img74.png)
2行
3行
![$ +2\times$](img76.png)
2行
![$ -\frac{1}{2}\times$](img78.png)
3行
![$ 2行-\frac{3}{2}\times$](img80.png)
3行
これで,ガウス・ジョルダン法による対角化の作業は完了である.これか
ら,
![$ (x_1,x_2,x_3)=(-1,0,1)$](img82.png)
と分かる.さらに,逆行列が
と分かる.
ピボット選択は行わないで,逆行列も求めないのガウス・ジョルダン法で連立方程式を計
算するプログラムを示す.このプログラムの動作は,次の通りである.
- 仮引数「n」は,解くべき連立方程式の未知数の数である.
- 仮引数の配列「a」と「b」は,係数行列
と非同次項
である.
値は,呼び出し元からにアドレス渡しで送られる.
- 係数行列は,配列「a[1][1]」〜「a[n][n]」に格納されている.
- 非同次項は,配列「b[1]」〜「b[n]」に格納されている.
- 連立方程式の解
は,配列「b[1]」〜「b[n]」に格納される.
- このプログラムでの処理が終了すると,配列「a[1][1]」〜「a[n][n]」は単位
行列になる..
/* ========== ガウスジョルダン法の関数 =================*/
void gauss_jordan(int n, double a[][100], double b[]){
int ipv, i, j;
double inv_pivot, temp;
for(ipv=1 ; ipv <= n ; ipv++){
/* ---- 対角成分=1(ピボット行の処理) ---- */
inv_pivot = 1.0/a[ipv][ipv];
for(j=1 ; j <= n ; j++){
a[ipv][j] *= inv_pivot;
}
b[ipv] *= inv_pivot;
/* ---- ピボット列=0(ピボット行以外の処理) ---- */
for(i=1 ; i<=n ; i++){
if(i != ipv){
temp = a[i][ipv];
for(j=1 ; j<=n ; j++){
a[i][j] -= temp*a[ipv][j];
}
b[i] -= temp*b[ipv];
}
}
}
}
ホームページ:
Yamamoto's laboratory著者:
山本昌志
Yamamoto Masashi
平成18年11月27日