情報処理演習 第4回 行列の直接解法 「Ruby 言語による Gauss-Joran の解法」 (講義ノート)

Ruby 言語による Gauss-Jordan の解法例は、例えば次のようになる。

#!/usr/koeki/bin/ruby 
require 'mathn'

def tr_mat(a)
  tr = Array.new
  for i in 0 .. a[0].length - 1
    col = Array.new
    for j in 0 .. a.length - 1
      col << a[j][i]
    end
    tr << col
  end
  tr
end

def id_mat(a)
  id = Array.new
  for i in 0 .. a.length - 1
    col = Array.new
    for j in 0 .. a[0].length - 1
      if j == i
        col << 1
      else
        col << 0
      end
    end
    id << col
  end
  id
end

def mul(p,q)
  multi = Array.new
  mul_p = tr_mat(p)  
  for i in 0 .. mul_p.length - 1
    col = Array.new
    for j in 0 .. q[0].length - 1
      sum = 0
      for k in 0 .. mul_p[0].length - 1
        sum += mul_p[i][k] * q[k][j]
      end
      sum
      col << sum
    end
    if col != Array.new(mul_p.length,0)
      multi << col
    end
  end
  multi
end

def p_mat(a,i,c)
  p_m = id_mat(a)
  p_m[i][i] = c
  mul(p_m, a)
end

def q_mat(a,i,j,c)
  q_m = id_mat(a)
  q_m[i][j] = c
  mul(q_m, a)
end

def r_mat(a,i,j)
  r_m = id_mat(a)
  k = r_m[i]
  l = r_m[j]
  r_m[j] = k
  r_m[i] = l
  r_m
  mul(r_m, a)
end

def pivot(a,i)
  c = Array.new
  for j in 0 .. a[0].length - 1
    c << a[i][j] / a[i][i]
  end
  a[i] = c
  a
end

mat = [[1,2,3,1],[4,5,6,2],[7,8,0,3]]

for i in 0 .. mat.length - 1
  mat = pivot(mat,i)
  for j in 0 .. mat.length - 1
    if j != i
      c = -1 * mat[j][i]
      mat = q_mat(mat,i,j,c)
    end
  end
end

ans = Array.new
for i in 0 .. mat.length - 1
  ans << mat[i].pop
end

p ans

tr_mat は転置行列を作る。 id_mat は単位行列を作る。 mul は行列の積を計算する。 上記の行基本変形に対応する行列をそれぞれ、p_mat, q_mat, r_mat とした。