Javaプログラミング

連立方程式

• 前進消去
着目する行から非対角項を対角項で割ったものを掛け、差し引くことで非対角項より下の行列の領域、下三角領域を全て順次０にし、上三角行列にする。
\[ \left[ \begin{array}{cccc|c} a_{11} &a_{12} & \cdots &a_{1n} &b_{1}\\ a_{21} &a_{22} & \cdots &a_{2n} &b_{2}\\ \vdots & \vdots & \ddots & \vdots & \vdots \\ a_{n1} &a_{n2} & \cdots &a_{nn} &b_{n}\\ \end{array} \right] \Rightarrow \left[ \begin{array}{cccc|c} a_{11}' &a_{12}' & \cdots &a_{1n}' &b_{1}'\\ 0 &a_{22}' & \cdots &a_{2n'} &b_{2}'\\ \vdots & \vdots & \ddots & \vdots & \vdots \\ 0 &0 & \cdots &a_{nn}' &b_{n}'\\ \end{array} \right] \] 対角項\(a_{kk}\) \((k=1\dots n)\)に着目する
\[ \left[ \begin{array}{ccccc|c} \ddots &　\vdots & \vdots& \cdots &\vdots &\vdots\\ \cdots &a_{kk} &a_{kj} & \cdots &a_{kn} &b_{k}\\ \cdots &a_{ik} &a_{ij} & \cdots &a_{in} &b_{i}\\ \vdots &　\vdots & \vdots& \ddots &\vdots &\vdots\\ \cdots & a_{nk} &\cdots&\cdots &a_{nn} &b_{n}\\ \end{array} \right] \] まずは行基本変形により対角項\(a_{kk}\)のより下の列要素を全て０にすることを考える。 そのためにはまず \[\alpha=\frac{a_{ik}}{a_{kk}}~~~~~(i=k+1\dots n)\] とし、 \[a_{i}{k}=a_{ik}-\alpha a_{kk}=0\] の操作をして\(a_{i}{k}\)の値を０にする。行基本変形の操作なので同じ行の値にも影響を与える。 すなわち \[a_{ij}'=a_{ij}-\alpha a_{kj}~~~~~(j=k+1\dots n)\] \[b_{i}'=b_{i}-\alpha b_{k}\] となる。　
これを最下行\((i=n,k=n-1)\)までくりかえすと以下の上三角行列が得られる。 以上が前進消去法である。 \[ \left[ \begin{array}{cccc|c} a_{11}' &a_{12}' & \cdots &a_{1n}' &b_{1}'\\ 0 &a_{22}' & \cdots &a_{2n}' &b_{2}'\\ \vdots & \vdots & \ddots & \vdots & \vdots \\ 0 &0 & \cdots &a_{nn}' &b_{n}'\\ \end{array} \right] \] この方法で以下の連立方程式を解いてみる。 \[ \left\{ \begin{array}{1} 6x_{1}+5x_{2}+4x_{3}=8\\ 12x_{1}+13x_{2}+10x_{3}=16\\ 18x_{1}+21x_{2}+17x_{3}=27\\ \end{array} \right. \] 拡大係数行列に変換して \[ \left[ \begin{array}{ccc|c} a_{11} &a_{12} &a_{13} &b_{1}\\ a_{21} &a_{22} &a_{23} &b_{2}\\ a_{31} &a_{32} &a_{33} &b_{3}\\ \end{array} \right] ＝ \left[ \begin{array}{ccc|c} 6 & 5 & 4 & 8\\ 12 &13 &10 &16\\ 18 &21 &17 &27\\ \end{array} \right] \]
まず\(k=1\)のときの対角要素\(a_{11}\)である6に着目する、これより下の列要素を全て0にしたいので \(i=2\)行目の要素12より\(α=\frac{12}{6}=2\)を得る。 \[a_{21}'=a_{21}-\alpha a_{11}=12-6\cdot2=0\]となり、残りの行要素 \[a_{22}'=a_{22}-\alpha a_{12}=13-10=3\] \[a_{23}'=a_{23}-\alpha a_{13}=10-8=2\] \[b_{2}'=b_{2}-\alpha b_{1}=16−16=0\] となり同様の計算を\(i=3,k=2\)行目の作業を行うと
\[ \left[ \begin{array}{ccc|c} a_{11} &a_{12} &a_{13} &b_{1}\\ a_{21} &a_{22} &a_{23} &b_{2}\\ a_{31} &a_{32} &a_{33} &b_{3}\\ \end{array} \right] ＝ \left[ \begin{array}{ccc|c} 6 & 5 & 4 & 8\\ 0 & 3 & 2 & 0\\ 0 & 6 & 5 & 3\\ \end{array} \right] \] が得られる。さらに\(k=1\)のときの対角要素\(a_{11}\)に着目し、同様の計算をおこなうと \[ \left[ \begin{array}{ccc|c} a_{11} &a_{12} &a_{13} &b_{1}\\ a_{21} &a_{22} &a_{23} &b_{2}\\ a_{31} &a_{32} &a_{33} &b_{3}\\ \end{array} \right] ＝ \left[ \begin{array}{ccc|c} 6 & 5 & 4 & 8\\ 0 & 3 & 2 & 0\\ 0 & 0 & 1 & 3\\ \end{array} \right] \] が得られる。これで上三角行列が得られたので前進消去は完了した。
• 後退代入
前進消去で得られた、上三角行列の下から順次代入して解を求める。 上三角行列になった場合、最下行の一つの要素\(a_{nn}'\)に対する\(b_{n}'\)は解になっているはずである。（あるいは両辺が定数倍された状態。） \[ \left[ \begin{array}{ccccc|c} a_{11}' &a_{12}' & \cdots &a_{1(n-1)}' &a_{1n}' &b_{1}'\\ 0 &a_{22}' & \cdots &a_{2(n-1)}' &a_{2n}' &b_{2}'\\ \vdots & \vdots & \ddots & \vdots & \vdots & \vdots \\ 0 &0 & \cdots &a_{(n-1)(n-1)}' &a_{(n-1)n}' &b_{n-1}'\\ 0 &0 & \cdots &0 &a_{nn}' &b_{n}'\\ \end{array} \right] \] よって。 \[x_{n}=\frac{b_{n}'}{a_{nn}'}\] として解 \(x_{n}\) を得る。 これを上の行\(a_{(n-1)n}'\)に代入操作をする。 \[ x_{n-1}=\frac{b_{n-1}'-a_{(n-1)n}'x_{n}}{a_{(n-1)(n-1)}'} \] これを上の行に向かって最初の行になるまで繰り返す。
一般化すると \[ x_k=\frac{b'_k-\sum_{j=k+1}^{n} a'_{kj}x_j}{a'_{kk}} ~~~~~ (k=n-1, \cdots 1) \] となる。 以上が後退代入法である。 この方法を使って先ほどの続きを計算すると \[ \left[ \begin{array}{1} a_{11}' & a_{12}' & a_{13}'\\ 0 & a_{22}' & a_{23}'\\ 0 & 0 & a_{33}'\ \end{array} \right] \left[ \begin{array}{1} x_1\\ x_2\\ z_3\\ \end{array} \right] = \left[ \begin{array}{1} b_{1}'\\ b_{2}'\\ b_{3}'\\ \end{array} \right] \] \[ \left[ \begin{array}{1} 6 & 5 & 4\\ 0 & 3 & 2\\ 0 & 0 & 1\\ \end{array} \right] \left[ \begin{array}{1} x_1\\ x_2\\ z_3\\ \end{array} \right] = \left[ \begin{array}{1} 8\\ 0\\ 3\\ \end{array} \right] \] \[ x_3=3 \] \[ x_2=\frac{b_2'-a_{23}'x_3}{a_{22}'} = \frac{0 - 2 \times 3}{3} = -2 \] \[ x_1=\frac{b_1' - (a_{12} x_2 + a_{13}x_3)}{a_{11}'} = \frac{8 - 5 \times (-2) -4 \times 3}{6} = 1 \] の様に解が得られる \[ \therefore x=1,y=-2,z=3 \]
• 解法一般化
これまでの解法を再度まとめると以下の様になる。
• 前進消去
\[ \alpha = \frac{a_{ik}}{a_{kk}} ~~~~~ (i=k+1, \cdots ,n) \] \[ a_{ij}=a_{ij}-\alpha \times a_{kj} ~~~~~ (j=k+1, \cdots , n) \] \[ b_j=b_i-\alpha \times b_j ~~~~~ (i=j+1, \cdots , n) \]
• 後退代入
\[ x_n=\frac{b'_n}{a'_{nn}}, ~~~~~ x_k=\frac{b'_k-\sum^{n}_{j=k+1}a'_{kj}x_j}{a'_kk} ~~~~~ (k=n-1, \cdots 1) \]
• ピボット選択
プログラムの実装においておおよその計算結果は同じだが、係数行列の対角項の絶対値が大きい時、計算誤差が小さくなることが知られている。 そのため対角項の絶対値が大きくなる様に行の並び替えを行う。この操作をピボット選択という。

//前進消去
for(int i=0;i<n;i++) {
	for(int j=i+1;j<n;j++) {
		double alpha=a[j][i]/a[i][i];
		for(int k=0;k<m;k++) {
			a[j][k]-=alpha*a[i][k];
		}
	}
}

//後退代入
x[n-1]=(a[n-1][n])/a[n-1][n-1];

for(int i=n-2;i>=0;i--) {
	x[i]=a[i][n];
	for(int j=n-1;j>=i+1;j--) {
		x[i]-=(a[i][j]*x[j]);
	}
	x[i]/=a[i][i];
}