mt19937-period.rb

P = 19937
N = 624
M = 397
W = 32
R = 31
MATRIX_A = 0x9908b0df
UMASK = 0x80000000
LMASK = 0x7fffffff
MAG01 = [0, MATRIX_A]

class State
  def initalize(vector)
    @vector = vector
    @index = 0 
  end
  
  def next
    @vector[@index] = at(M) ^ mult_a((at(0) & UMASK) | (at(1) & LMASK))
    @index = (@index + 1) % N
  end
  
  def b()
    at(1) & 1
  end
  
  def at(i)
    @vector[(@index + i) % N]
  end
  
  def mult_a(v)
    (v >> 1) ^ MAG01[v & 1]
  end
end

def phi_inverse(bits)
  x = [0] + bits.to_a
  P.downto(N) do |i|
    y = (x[i] ^ x[i-N+M] ^ ((x[i-N] & 1) * MATRIX_A)) << 1
    x[i-N+1] = (x[i-N+1] & UMASK) | (y & LMASK) | (x[i-N+1] & 1) 
    x[i-N] = (y & UMASK) | (x[i-N] & LMASK)
  end
  new State(y[0, P])
end

def generate(state)
  bits = []
  P.times do
    bits << state.b()
    state.next()
    state.next()
  end
  bits
end

def test()
  x = initial = Array.new(P) {|i| i == 1 ? 1 : 0 }
  P.times do
    state = phi_inverse(x)
    x = generate(state)
  end
  x == initial
end