From: Shin-ichiro HARA <sinara@...>
Date: 2005-11-25T00:24:54+09:00
Subject: [ruby-dev:27817] Re: Matrix class is broken without mathn

原です。

In [ruby-dev:27790]
>けいじゅ＠いしつかです.

>>そういえば整数行列の逆行列はどうしよう。
>
>整数行列の逆行列が整数行列になることがあるからですか?
>
>そこまでは考えなくてよいんじゃないでしょうか? 

そうですよね。特殊なケースまでmatrix.rbで面倒を見る必要はな
いと私も思います。

、、、とは言え、もし実装したらどうなるか興味が沸いたので、やっ
てみました。

class Matrix
  def inverse_e
    size = row_size
    a = to_a
    ar = Matrix.I(size).to_a
    d = 0
    while d < size
      if left_sweep_e(d, a, ar)
        Matrix.Raise ErrNotRegular if a[d][d].zero?
        d += 1
      end
    end

    Matrix[*ar].inverse_from(Matrix[*a])
  end

  def left_sweep_e(d, a, ar)
    r_size, c_size = a.size, ar.size
    clean = true
    if i = (d...r_size).find{|x| !a[x][d].zero?}
      i = (i...r_size).find_all{|x| !a[x][d].zero?}.
          min{|x, y| a[x][d].abs <=> a[y][d].abs}
      if i > d
        a[d], a[i] = a[i], a[d]
        ar[d], ar[i] = ar[i], ar[d]
      end
      (d+1...r_size).each do |i|
        next if a[i][d].zero?
        clean = false
        q = a[i][d] / a[d][d]
        (d...c_size).each do |j|
          a[i][j] -= a[d][j]*q
        end
        (0...r_size).each do |j|
          ar[i][j] -= ar[d][j]*q
        end
      end
    end
    clean
  end
end

ここで使っているMatrix#inverse_fromは「quo版」ではなくて、
オリジナルの「/版」です。

他にも単因子を求めるアルゴリズムとか作ってみたのですが、こう
なるとmatrix.rbじゃなくてmatrix-i.rbとか別ファイルにしたほう
がいいかもしれません。