Ruby - 連立方程式解法（ガウスの消去法（ピボット選択））！

かつて、連立方程式を「ガウスの消去法」で解くアルゴリズムを Ruby で実装したことを紹介しました。

• Ruby - 連立方程式解法（ガウスの消去法）！

しかし、計算途中で対角成分がゼロになるケースでは計算ができませんでした。

今回はその問題を解決すべく、「ガウスの消去法（ピボット選択）」で解くアルゴリズムを実装してみました。（実際、連立方程式を使うような自作ツールではこのピボット選択のガウスの消去法を使用しているが、ブログでは紹介していなかった）

0. 前提条件

• Debian GNU/Linux 11.3 での作業を想定。
• Ruby 3.1.1 での実行を想定。
• 連立方程式の解法（ガウスの消去法（ピボット選択））についての説明は割愛。（Web 上等で容易に確認可能）

1. Ruby スクリプト作成

• Shebang ストリング（1行目）では、フルパスでコマンド指定している。（当方の慣習）

File: gauss_elimination_pivot.rb

#! /usr/local/bin/ruby
#*********************************************
# 連立方程式の解法 ( ガウスの消去法（ピボット選択） )
#*********************************************
#
class GaussElimination
  def initialize
    # 係数
    @a = [
      [1.0, 2.0, 7.0, 6.0,  6.0],
      [2.0, 4.0, 4.0, 2.0,  2.0],
      [1.0, 8.0, 5.0, 2.0, 12.0],
      [2.0, 4.0, 3.0, 3.0,  5.0]
    ]
    # 次元の数
    @n = @a.size
  end

  # 計算・結果出力
  def exec
    # 元の連立方程式をコンソール出力
    display_equations
    # 前進消去
    (@n - 1).times do |k|
      raise unless pivot(k)
      (k + 1).upto(@n - 1) do |i|
        d = @a[i][k] / @a[k][k].to_f
        (k + 1).upto(@n) do |j|
          @a[i][j] -= @a[k][j] * d
        end
      end
    end
    # 後退代入
    (@n - 1).downto(0) do |i|
      d = @a[i][@n]
      (i + 1).upto(@n - 1) do |j|
        d -= @a[i][j] * @a[j][@n]
      end
      @a[i][@n] = d / @a[i][i].to_f
    end
    # 結果出力
    display_answers
  rescue => e
    raise
  end

private

  # ピボット選択
  def pivot(k)
    pv = k
    v_max = @a[k][k].abs
    (k + 1).upto(@n - 1) do |i|
      if @a[i][k].abs > v_max
        pv = i
        v_max = @a[i][k]
      end
    end
    if k != pv
      dummy   = @a[k]
      @a[k]  = @a[pv]
      @a[pv] = dummy
    end
    return true
  rescue => e
    return false
  end

  # 元の連立方程式をコンソール出力
  def display_equations
    @n.times do |i|
      @n.times { |j| printf("%+dx%d ", @a[i][j], j + 1) }
      puts "= %+d" % @a[i][@n]
    end
  rescue => e
    raise
  end

  # 結果出力
  def display_answers
    0.upto(@n - 1) { |k| puts "x%d = %f" % [k + 1, @a[k][@n]] }
  rescue => e
    raise
  end
end

if __FILE__ == $0
  begin
    obj = GaussElimination.new
    obj.exec
  rescue => e
    $stderr.puts "[#{e.class}] #{e.message}"
    e.backtrace.each{ |tr| $stderr.puts "\t#{tr}" }
  end
end

• Gist - Ruby script to solve simultaneous equations with Gauss elimination method(pivot).

2. 実行

まず、実行権限を付与。

$ chmod +x gauss_elimination_pivot.rb

そして、実行

$ ./gauss_elimination_pivot.rb
+1x1 +2x2 +7x3 +6x4 = +6
+2x1 +4x2 +4x3 +2x4 = +2
+1x1 +8x2 +5x3 +2x4 = +12
+2x1 +4x2 +3x3 +3x4 = +5
x1 = -3.000000
x2 = 2.000000
x3 = -1.000000
x4 = 2.000000

ちなみに、今回の例は、ピボット選択しない場合には計算ができないケース。
試しに 32行目(raise ...) をコメントアウトすると、以下のように計算できなくなる。

$ ./gauss_elimination_pivot.rb
+1x1 +2x2 +7x3 +6x4 = +6
+2x1 +4x2 +4x3 +2x4 = +2
+1x1 +8x2 +5x3 +2x4 = +12
+2x1 +4x2 +3x3 +3x4 = +5
x1 = NaN
x2 = NaN
x3 = NaN
x4 = NaN

---

以上。