3 連立一次方程式(消去法)

3.1 連立方程式の表現方法

連立1次方程式(Linear Equations)は，次のような形をしている．

$$\begin{aligned}a_{11}x_1+a_{12}x_2+a_{13}x_3+\cdots+a_{1N}x_N&=... ...a_{M1}x_1+a_{M2}x_2+a_{M3}x_3+\cdots+a_{MN}x_N&=b_M \end{aligned}$$

式(25)は行列とベクトルで書くと，式がすっきりして 考えやすくなる．書き直すと，

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

である．それぞれの行列とベクトルは，

$$\begin{aligned}\boldsymbol{A}&= \begin{bmatrix}a_{11} & a_{12} ... ...rix}b_1 \\ b_2 \\ b_3 \\ \vdots\\ b_N \end{bmatrix} \end{aligned}$$

を表す．

通常，連立1次方程式(25)は

$$\begin{aligned}\begin{bmatrix}a_{11} & a_{12} & a_{13} & \cdots... ...rix}b_1 \\ b_2 \\ b_3 \\ \vdots\\ b_N \end{bmatrix} \end{aligned}$$

と書き表せる．このようにすると，見通しがかなり良くなる．

3.2 ガウス・ジョルダン法の基本的な考え方

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

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

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

$$\left\{ \begin{aligned}3x+6y+9z&=6\\ 2x+2y+3z&=1\\ 2x+2y+1z&=-1 \end{aligned} \right.$$

をガウス・ジョルダン法で解を求める． 解くべき，方程式

$$\begin{bmatrix}3 & 6 & 9 \\ 2 & 2 & 3 \\ 2 & 2 & 1 \end{bmatrix} ... ...rix}x_1 \\ x_2 \\ x_3 \end{bmatrix} = \begin{bmatrix}9 \\ 1 \\ -1 \end{bmatrix}$$

$\frac{1}{3}\times$1行

$$\begin{bmatrix}1 & 2 & 3 \\ 2 & 2 & 3 \\ 2 & 2 & 1 \end{bmatrix} ... ...rix}x_1 \\ x_2 \\ x_3 \end{bmatrix} = \begin{bmatrix}2 \\ 1 \\ -1 \end{bmatrix}$$

2行$-2\times$1行

$$\begin{bmatrix}1 & 2 & 3 \\ 0 & -2 & -3 \\ 2 & 2 & 1 \end{bmatrix... ...ix}x_1 \\ x_2 \\ x_3 \end{bmatrix} = \begin{bmatrix}2 \\ -3 \\ -1 \end{bmatrix}$$

3行$-2\times$1行

$$\begin{bmatrix}1 & 2 & 3 \\ 0 & -2 & -3 \\ 0 & -2 & -5 \end{bmatr... ...ix}x_1 \\ x_2 \\ x_3 \end{bmatrix} = \begin{bmatrix}2 \\ -3 \\ -5 \end{bmatrix}$$

$-\frac{1}{2}\times$2行

$$\begin{bmatrix}1 & 2 & 3 \\ 0 & 1 & \frac{3}{2} \\ 0 & -2 & -5 \e... ... x_2 \\ x_3 \end{bmatrix} = \begin{bmatrix}2 \\ \frac{3}{2} \\ -5 \end{bmatrix}$$

1行$-2\times$2行

$$\begin{bmatrix}1 & 0 & 0 \\ 0 & 1 & \frac{3}{2} \\ 0 & -2 & -5 \e... ...x_2 \\ x_3 \end{bmatrix} = \begin{bmatrix}-1 \\ \frac{3}{2} \\ -5 \end{bmatrix}$$

3行$+2\times$2行

$$\begin{bmatrix}1 & 0 & 0 \\ 0 & 1 & \frac{3}{2} \\ 0 & 0 & -2 \en... ...x_2 \\ x_3 \end{bmatrix} = \begin{bmatrix}-1 \\ \frac{3}{2} \\ -2 \end{bmatrix}$$

$-\frac{1}{2}\times$3行

$$\begin{bmatrix}1 & 0 & 0 \\ 0 & 1 & \frac{3}{2} \\ 0 & 0 & 1 \end... ... x_2 \\ x_3 \end{bmatrix} = \begin{bmatrix}-1 \\ \frac{3}{2} \\ 1 \end{bmatrix}$$

$2行-\frac{3}{2}\times$3行

$$\begin{bmatrix}1 & 0 & 0 \\ 0 & 1 & 0 \\ 0 & 0 & 1 \end{bmatrix} ... ...rix}x_1 \\ x_2 \\ x_3 \end{bmatrix} = \begin{bmatrix}-1 \\ 0 \\ 1 \end{bmatrix}$$

これで，ガウス・ジョルダン法による対角化の作業は完了である．これか ら， $(x_1,x_2,x_3)=(-1,0,1)$と分かる．

3.3 ガウス・ジョルダン法のC言語の関数

ピボット選択は行わないで，逆行列も求めないのガウス・ジョルダン法で連立方程式を計 算するプログラムを示す．このプログラムの動作は，次の通りである．

• 仮引数「n」は，解くべき連立方程式の未知数の数である．
• 仮引数の配列「a」と「b」は，係数行列 $\boldsymbol{A}$と非同次項 $\boldsymbol{b}$である． 値は，呼び出し元からにアドレス渡しで送られる．
  • 係数行列は，配列「a[1][1]」〜「a[n][n]」に格納されている．
  • 非同次項は，配列「b[1]」〜「b[n]」に格納されている．
• 連立方程式の解 $\boldsymbol{x}$は，配列「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
著者: 山本昌志