コンピュータは有限桁の数値しか扱う事はできないので、桁数の多い場合や無限小数の場合は四捨五入され切り捨てられます。なので実際の数値とは多少の誤差が生じますが、これを丸め誤差といいます。
なので、コンピュータを用いて連立方程式等の計算すると多少の誤差が生じますが、誤差が最小になるように高精度で計算をすれば限りなく正しい解を求めることができます。
今回は、連立一次方程式を解くための「ヤコビ法」と「ガウス・ザイデル法」についてプログラムと共に紹介します。
ヤコビ法，ガウス・ザイデル法
ヤコビ法 -Jacobi Method-
ヤコビ法とは、n元の連立一次方程式の連立一次方程式を反復法で解く手法の１つです。
を反復法で解く手法の１つです。元の方程式である[Ax = b]を[r = b - Ax]と変形し、rの値が小さくなるようにxの値を反復的に変えていき、r = 0を満たすxを見つけることができれば、方程式が解けたことになります。
例えば以下のような3元の連立一次方程式を解くことを考える。
回目の反復で得られたの値をと書く。 初期値は(0,0,0)とする。
は(0,0,0)とする
という反復を繰り返していく。また、収束判定として、適当な正数εに対して、
また、
が成り立てば終了する。収束判定の値であるεを厳しくすれば精度が高まります。
次のような簡単な3元連立方程式を解いてみる。
3a + 2b + c = 10
a + 4b + c = 12
2a + 2b + 5c = 21
以下にPythonのプログラムです。
class Jacobi():
def calculation(self):
coefficient_matrix = [[3,2,1],[1,4,1],[2,2,5]]
y = [10,12,21]
convergence_value = 10 ** -16
length = len(y)
x = [0] * length
x2 = [0] * length
counter = 1;
while True:
for i in range(length):
x2[i] = (y[i] - coefficient_matrix[i][(i+1)%length]*x[(i+1)%length] - coefficient_matrix[i][(i+2)%length]*x[(i+2)%length]) / coefficient_matrix[i][i]
judge_counter = 0;
for i in range(length):
if abs(x2[i] - x[i]) < convergence_value:
judge_counter = judge_counter + 1
if length == judge_counter:
break
print("------" + str(counter) + "回目" + "------")
counter = counter + 1
for i in range(length):
print(x2[i])
for i in range(length):
x[i] = x2[i]
if __name__ == "__main__":
instance = Jacobi()
instance.calculation()実行結果：
ガウス・ザイデル法 -Gauss Seidel-
ガウス・ザイデル法も同様にn元の連立一次方程式を解く反復法ので解く手法の１つである。
以下のように、k回反復した後に、k+1回目の解を求める際にすでに求めた解を次の方程式に利用する。
このようにすることで、ガウス・ザイデル法はヤコビ法よりも収束の速度が速くなります。またガウス・ザイデル法では反復計算で書き換えた変数を使うだけであるが、ヤコビ法では「1つ前の値」を記憶する必要があるので、プログラム上でもガウス・ザイデル法の方が簡単であることが言える。
同様の3元連立一次方程式を解きます。以下がPythonのプログラムです。
class Gauss_Seidel():
def calculation(self):
coefficient_matrix = [[3,2,1],[1,4,1],[2,2,5]]
y = [10,12,21]
convergence_value = 10 ** -16
length = len(y)
x = [0] * length
x_before = [0] * length
counter = 1;
while True:
for i in range(length):
x_before[i] = x[i]
x[i] = (y[i] - coefficient_matrix[i][(i+1)%length]*x[(i+1)%length] - coefficient_matrix[i][(i+2)%length]*x[(i+2)%length]) / coefficient_matrix[i][i]
judge_counter = 0;
for i in range(length):
if abs(x_after[i] - x[i]) < convergence_value:
judge_counter = judge_counter + 1
if length == judge_counter:
break
print("------" + str(counter) + "回目" + "------")
counter = counter + 1
for i in range(length):
print(x[i])
if __name__ == "__main__":
instance = Gauss_Seidel()
instance.calculation()実行結果：
ヤコビ法より、計算回数が非常に速いことがわかります。
ヤコビ法とガウス・ザイデル法の特徴
- ガウス・ザイデル法は、ヤコビ法よりも収束する速度が速く、計算量が少ない。
- プログラムもガウス・ザイデル法の方が簡潔である。（反復計算で書き換えた変数をそのまま用いるため）
しかし、ガウス・ザイデル法は書き換えた変数を用いるため逐次計算となる。そのため並列計算には向かない。また、ヤコビ法は各式で独立に計算ができるため、並列計算に向いている。
以上から、一概にガウス・ザイデル法が良いとは言えません。また、ヤコビ法やガウス・ザイデル法を緩和係数を用いて組み合わせた「SOR法」もありますが、今回は割愛します。
