Skip to content

Instantly share code, notes, and snippets.

Embed
What would you like to do?
Ruby script to calculate an inverse matrix by cofactor matrix.
#! /usr/local/bin/ruby
#***************************************************************
# Ruby script to calculate an inverse matrix by cofactor matrix.
#***************************************************************
#
class TestInverseMatrix
#MTX = [
# [2]
#]
#MTX = [
# [1, 2],
# [3, 4],
#]
#MTX = [
# [1, 4, 7],
# [8, 5, 2],
# [3, 9, 6]
#]
#MTX = [
# [1, 3, 5, 7],
# [4, 2, 8, 6],
# [9, 3, 7, 5],
# [4, 0, 6, 8]
#]
MTX = [
[1, 3, 5, 7, 9],
[4, 2, 8, 6, 0],
[9, 3, 7, 5, 1],
[4, 0, 6, 8, 2],
[3, 6, 9, 2, 5]
]
def exec
puts "Matrix A ="
MTX.each { |row| p row }
puts "det(A) = #{MTX.det}"
puts "Cofactor Matrix of A ="
MTX.cofactor_matrix.each { |row| p row }
puts "Inverse Matrix of A ="
MTX.inverse_matrix.each { |row| p row }
rescue => e
$stderr.puts "[#{e.class}] #{e.message}"
e.backtrace.each { |tr| $stderr.puts "\t#{trace}" }
exit 1
end
end
class Array
# 行列式の計算
# * 当メソッドは再帰的に呼ばれる
def det
# 以下の場合は例外スロー
# - 自身配列が Array クラスでない
# - 自身配列が空
# - 自身配列が正方行列でない
raise "Self is not a Array class!" unless self.class == Array
raise "Self array is nil!" if self.size == 0
raise "Self array is not a square matrix!" unless self.size == self[0].size
m = self
n = m.size
# 3行3列以下の行列の場合、
# => Sarrus による計算
case n
when 1; return m[0][0]
when 2; return m[0][0] * m[1][1] \
- m[0][1] * m[1][0]
when 3; return m[0][0] * m[1][1] * m[2][2] \
+ m[0][1] * m[1][2] * m[2][0] \
+ m[0][2] * m[1][0] * m[2][1] \
- m[0][0] * m[1][2] * m[2][1] \
- m[0][1] * m[1][0] * m[2][2] \
- m[0][2] * m[1][1] * m[2][0]
end
# 4行4列以上の行列の場合、
# => 第1行についての余因子展開(ラプラス展開)による計算
return n.times.map { |j| m[0][j] * cofactor(m, 0, j) }.sum
end
# 余因子行列の計算
def cofactor_matrix
# 以下の場合は例外スロー
# - 自身配列が Array クラスでない
# - 自身配列が空
# - 自身配列が正方行列でない
# - 自身配列のサイズが1x1の場合
raise "Self is not a Array class!" unless self.class == Array
raise "Self array is nil!" if self.size == 0
raise "Self array is not a square matrix!" unless self.size == self[0].size
raise "Can not calculate a cofactor matrix!" if self.size == 1
m = self
n = m.size
return n.times.map do |i|
n.times.map { |j| cofactor(m, i, j) }
end.transpose
end
# 逆行列の計算
# * 内部で余因子行列の計算も行っている
def inverse_matrix
# 以下の場合は例外スロー
# - 自身配列が Array クラスでない
# - 自身配列が空
# - 自身配列が正方行列でない
raise "Self is not a Array class!" unless self.class == Array
raise "Self array is nil!" if self.size == 0
raise "Self array is not a square matrix!" unless self.size == self[0].size
m = self
return [[1.0 / m[0][0]]] if m.size == 1
d = m.det
return m.cofactor_matrix.map do |row|
row.map { |col| col / d.to_f }
end
end
private
# 余因子(正方行列 A の (i, j) に対する)
# [方法-1] 1行ずつチェックしながら生成する方法
# [方法-2] 引数で与えた行列(配列)を一旦いわゆる「深いコピー」で
# 複製してから、 i 行と j 列を削除する方法
#
# @param m: 行列配列
# @param i: 行インデックス
# @param j: 列インデックス
# @return m2: 余因子
def cofactor(m, i, j)
# ====[ 方法-1 ]==>
#m2 = []
#n = m.size
#n.times do |k|
# next if k == i
# row = []
# n.times do |l|
# next if l == j
# row << m[k][l]
# end
# m2 << row
#end
# <===[ 方法-1 ]===
# ====[ 方法-2 ]==>
m2 = Marshal.load(Marshal.dump(m))
m2.delete_at(i)
m2 = m2.transpose
m2.delete_at(j)
m2 = m2.transpose
# <===[ 方法-2 ]===
return (-1)**(i+j) * m2.det
rescue => e
raise
end
end
TestInverseMatrix.new.exec if __FILE__ == $0
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment
You can’t perform that action at this time.