#######################################################
#                                                     #
#  This is test script for 'gaussian-elimination.rb'  #
#                                                     #
#######################################################
require 'test/unit'
require 'algebra'

class TestGaussianElimination < Test::Unit::TestCase
  def setup
    @G = MatrixAlgebra(Rational, 2, 2)
    @R = MPolynomial(Rational)
    @F2 = ResidueClassRing(Integer, 2)
    @M = SquareMatrix(@F2, 8)
    @MR = SquareMatrix(Rational, 5)
    @MI = SquareMatrix(Integer, 5)
    @MR56 = MatrixAlgebra(Rational, 5, 7)
    @MI56 = MatrixAlgebra(Integer, 5, 7)
  end

  def test_left_elimination
    a0 = [
      [3, 2, 0, -1, 7, 2, 3],
      [-1, -5, 0, 2, 5, 1, 2],
      [2, -3, 0, 1, 12, 3, 5],
      [9, 2, 0, -8, 3, 4, 2],
      [7, 5, 0, -9, -9, 1, -3]
    ]
    # puts
    a = @MR56[*a0]
    # a.display
    # puts "=>"
    a.left_eliminate!
    # a.display
    b0 = [
      [1,   0,   0, 0, Rational(319) / 85, Rational(82) / 85, Rational(134) / 85],
      [0,   1,   0, 0, -Rational(4) / 17, -Rational(3) / 17, -Rational(2) / 17],
      [0,   0,   0, 1, Rational(322) / 85, Rational(46) / 85, Rational(127) / 85],
      [0,   0,   0, 0,   0,   0,   0],
      [0,   0,   0, 0,   0,   0,   0]
    ]
    assert_equal(@MR56[*b0], a)

    # puts
    #
    a = @MI56[*a0]
    # a.display
    # puts "=>"
    a.left_eliminate_euclidian!
    # a.display
    b0 = [
      [1,   5,   0, -2, -5, -1, -2],
      [0,   1,   0, -20, -76, -11, -30],
      [0,   0,   0, -85, -322, -46, -127],
      [0,   0,   0,   0,   0,   0,   0],
      [0,   0,   0,   0,   0,   0,   0]
    ]
    assert_equal(@MI56[*b0], a)
  end

  def test_determinant_by_elimination
    a0 = [
      [3, 2, -1, 7, 2],
      [-1, -5, 2, 5, 1],
      [2, 8, 2, -3, 1],
      [9, 2, -8, 3, 4],
      [3, -1, 3, 4, 5]
    ]
    # puts
    d = @MR[*a0].determinant_by_elimination
    # p [d, d.class]
    assert_equal(1383, d)
    d = @MI[*a0].determinant_by_elimination_euclidian
    # p [d, d.class]
    assert_equal(1383, d)
  end

  def test_kernel_basis_2x2
    # puts
    a = @G[[0, 0], [1, -1]]
    k = a.kernel_basis
    # p Hash[a, k]
    assert_equal([Vector(Rational, 2)[1, 1]], k, (Hash[a, k]).to_s)
  end

  def test_kernel_basis_8x8_Z2
    # puts "Kernel of"
    m = @M[
      '10001101'.scan(/./).collect(&:to_i),
      '00000011'.scan(/./).collect(&:to_i),
      '01001011'.scan(/./).collect(&:to_i),
      '00001111'.scan(/./).collect(&:to_i),
      '00101011'.scan(/./).collect(&:to_i),
      '00001001'.scan(/./).collect(&:to_i),
      '00011011'.scan(/./).collect(&:to_i),
      '00000111'.scan(/./).collect(&:to_i)
    ]
    cm = m - m.unity
    # cm.display
    # puts "is"
    k = cm.kernel_basis
    # k.each do |b|
    #   puts "Vector(F2, 8)#{b}"
    # end
    assert_equal([
                   Vector(@F2, 8)[1, 0, 0, 0, 0, 0, 0, 0],
                   Vector(@F2, 8)[0, 1, 1, 1, 1, 1, 1, 0],
                   Vector(@F2, 8)[0, 1, 1, 0, 1, 0, 0, 1]
                 ], k)
  end
end

#  R = MPolynomial(Integer)
# a, b, c, d, e, f, g, h, i, j, k, l, x, y, z = R.vars("abcdefghijklxyz")
# M3 = MatrixAlgebra(R, 3, 4)
# m = M3[
#  [a, b, c, d],
#  [e, f, g, h],
#  [i, j, k, l]
# ]
# m.display
# m.swap_r(0, 2).display
#  m.set_column(1, [1,2,3]).display
# m.swap_c(0, 2).display
# m.multiply_r(2, x).display
# m.multiply_c(0, x).display
# m.mix_r(0, 2, x).display
# m.mix_c(0, 2, x).display
