Gauss-Seidel.py の技術ドキュメント

プログラムの動作

本プログラム Gauss-Seidel.py は、ガウス・ザイデル法を用いて特定の形式の連立一次方程式を反復的に解くことを目的としています。

主な機能は以下の通りです。

連立一次方程式の数値解法: 係数行列 $H$ と定数ベクトル $\mathbf{y}$ が与えられたとき、$H \mathbf{x} = \mathbf{y}$ を満たす解ベクトル $\mathbf{x}$ を数値的に求めます。
反復計算の進捗表示: 各反復ステップでの現在の解ベクトル $\mathbf{x}$ の値を標準出力に表示します。
収束判定: 解ベクトルが連続する反復間で指定された許容誤差 eps 以下になった場合、計算を停止し収束したことを報告します。
最大反復回数: 指定された最大反復回数 nmaxiter に達しても収束しない場合、計算を停止し、収束しなかったことを報告します。

これにより、比較的単純な構造を持つ連立一次方程式を効率的に解くことができます。

原理

本プログラムは、線形方程式系 $H \mathbf{x} = \mathbf{y}$ を解くための反復法であるガウス・ザイデル法を実装しています。ここで、$H$ は $n \times n$ の係数行列、$x$ は $n$ 次の未知数ベクトル、$y$ は $n$ 次の定数ベクトルです。

ガウス・ザイデル法は、初期推測値 $\mathbf{x}^{(0)}$ から出発し、以下の更新式を用いて解ベクトル $\mathbf{x}$ を反復的に改善していきます。

各反復 $k+1$ における $x_i$ の更新式は、通常次のように定義されます。 $$x_i^{(k+1)} = \frac{1}{H_{ii}} \left( y_i - \sum_{j < i} H_{ij} x_j^{(k+1)} - \sum_{j > i} H_{ij} x_j^{(k)} \right)$$

しかし、Gauss-Seidel.py の実装では、各行 $i$ において、非対角項として $H_{i, (i+1)\%n}$ と $H_{i, (i+2)\%n}$ のみが考慮される特別な形式の連立方程式を解いています。ここで %n はモジュロ演算を示し、インデックスが循環することを意味します。したがって、本プログラムで解かれる連立一次方程式の各行は、以下の形式を持つものと仮定されます。

\[H_{ii} x_i + H_{i, (i+1)\%n} x_{(i+1)\%n} + H_{i, (i+2)\%n} x_{(i+2)\%n} = y_i\]

この特殊な形式に基づき、各反復 $k+1$ における $x_i$ の更新式は、プログラム中で次のように実装されています。

\[x_i^{(k+1)} = \frac{1}{H_{ii}} \left( y_i - H_{i, (i+1)\%n} x_{(i+1)\%n}^{*} - H_{i, (i+2)\%n} x_{(i+2)\%n}^{*} \right)\]

ここで、$x_j^*$ は、現在の反復ループ内で既に更新された成分 (すなわち $j < i$ の場合は $x_j^{(k+1)}$) または前回の反復における成分 (すなわち $j \ge i$ の場合は $x_j^{(k)}$) を参照します。プログラム内では、x 配列がインプレースで更新されるため、x[i1] や x[i2] は、参照されるインデックス i1 や i2 が現在の i よりも小さい場合は最新の値が、大きい場合は前回の反復の値が自動的に利用されることになります。

収束判定条件: アルゴリズムは、連続する反復における解ベクトルの変化が、事前に設定された許容誤差 eps より小さくなったときに収束したと判断します。具体的には、すべての成分 $i$ について以下の条件が満たされた場合です。

\[|x_i^{(k+1)} - x_i^{(k)}| < \epsilon\]

もし最大反復回数 nmaxiter に達してもこの条件が満たされない場合、プログラムは収束しなかったと判断します。

必要な非標準ライブラリとインストール方法

このプログラムはPythonの標準ライブラリのみを使用しており、numpy などの非標準ライブラリは必要ありません。そのため、特別なインストール手順は不要です。

必要な入力ファイル

本プログラムは、外部の入力ファイルを必要としません。必要なデータ（係数行列 Hij、定数ベクトル yi、許容誤差 eps、最大反復回数 nmaxiter）は、プログラムの if __name__ == "__main__": ブロック内に直接ハードコードされています。

Hij: 係数行列 $H$ を表すリストのリスト（例: [[3, 2, 1], [1, 4, 1], [2, 2, 5]]）。
yi: 定数ベクトル $\mathbf{y}$ を表すリスト（例: [10, 12, 21]）。
eps: 収束判定のための許容誤差を表す浮動小数点数（例: 1.0e-16）。
nmaxiter: 最大反復回数を表す整数（例: 300）。

これらの値は、プログラムコードを直接編集することで変更可能です。

生成される出力ファイル

Gauss-Seidel.py は、実行時にいかなるファイルも生成しません。すべての結果は標準出力（コンソール）に表示されます。

出力される内容は以下の通りです。

各反復ごとに、その反復回数と、その時点での解ベクトル x の値。
計算が収束した場合、「Converged」というメッセージ。
計算が最大反復回数に達しても収束しなかった場合、「Not converged」というメッセージ。

コマンドラインでの使用例 (Usage)

Gauss-Seidel.py は引数を取らないため、基本的な実行コマンドは以下のようになります。

python Gauss-Seidel.py

コマンドラインでの具体的な使用例

以下のコマンドを実行します。

python Gauss-Seidel.py

このプログラムにハードコードされているデータは以下の通りです。

係数行列 $H$: [[3, 2, 1], [1, 4, 1], [2, 2, 5]]
定数ベクトル $\mathbf{y}$: [10, 12, 21]
許容誤差 eps: 1.0e-16
最大反復回数 nmaxiter: 300

この入力を元に、プログラムは以下のような出力を標準出力に表示します（出力の一部を抜粋）。

1-th iteration
  x= [3.3333333333333335, 1.3333333333333335, 2.0]
2-th iteration
  x= [2.888888888888889, 1.8333333333333335, 2.5]
3-th iteration
  x= [2.5, 2.0833333333333335, 2.8]
...
13-th iteration
  x= [2.0000000000000004, 2.0000000000000004, 3.0]
Converged

実行結果の説明:

上記の例では、プログラムは13回目の反復で収束条件を満たし、計算を停止しています。最終的な解ベクトル x は [2.0, 2.0, 3.0] に非常に近い値を示していることが分かります。各反復での x の値が逐次表示され、最終的に「Converged」メッセージが出力されることで、指定された許容誤差内で解が求められたことが確認できます。