terras6.rb

def f2(n)
	return n.odd? ? (n * 3 + 1) / 2 : n / 2

end

def re(n, t)

	s = (n % 2).to_s
	n1 = n
	t.times \
	{
		n = f2(n)
		s << (n % 2).to_s
		break if (n == 1)
	}
	return s, n
end

def adj(x, m, n, s2, p2, p, x3)
	s2.replace(p[0...x] + (x3 % 2).to_s)
	m = (s2[x, 1] == p[x, 1])
	n = p2.reverse.to_i(2)
	p2[x, 1] = m ? '11' : '01'
	return m, n, x3
end

def terras121(p)
	
	p2 = ['01', '11'][p[0, 1].to_i]
	n = 1
	n3 = 0
	m = p2 == '11'
	x3 = (p2.reverse.to_i(2)) >> 1

	(1...p.length).each \
	{
		|x|

		if (!m)
			ns = n.to_s(2)
			ns[0, 1] = ''
			m1 = x3 - 3**n3
			n1 = f2(m1)
		end

		n3 = p[0...x].split('').select { |z| z == '1' }.size
		
		s2 = ''
		m, n, x3 = adj(x, m, n, s2, p2, p, m ? (f2(x3) + 3**n3) : (n1 + 3**n3))
	}

	n = p2.reverse.to_i(2)
	s, n2 = re(n, p.length - 1)
	raise if (s != p)

	return n
end

def dense(w, d)
	w2 = w - 1
	
	a = (0...w2).to_a
	s = '0' * w2
	(1..(d * w - 1)).map { a.delete_at(rand(a.size)) }.each { |x| s[x, 1] = '1' }
	return ('1' + s)
	
end

def d(s)
        c = s.split('').select { |x| x == '1' }.size
	d = c.to_f / s.length
        return d
        
end

def cut(b, w, c)

	a = (1...w).to_a
	p = (1..c).map { a.delete_at(rand(a.size)) }
	
	s = b.to_s * w
	p.sort.reverse.each { |x| s[x, 0] = '-' }
	return s.split('-')
end

def entropy(w, c)
	c = [1, c].max
	c = [w - 1, c].min
	
	c2 = (c + 1) / 2
	w2 = w / 2
	l1 = c.even? ? cut(1, w2, c2) : cut(1, w2, c2 - 1)
	l0 = cut(0, w2, c2 - 1)
	l0 << '' if (c.even?)
	return [l1, l0].transpose.join
	
end

def count(m, n)
	n1 = n
	l = [n]
	while (n != 1)
		n = f2(n)
		l << n
	end
	cg = l.index{ |x| x < l[0] }
	cm1 = (0...l.size).max_by { |x| l[x] }
	return {"#{m}_{cg}" => cg, "#{m}_{cm1}" => cm1, "#{m}_{c}" => l.size}
	
end

c = 200
w = 50
f = File.open('gnuplot.cmd', 'w')
f.puts('$dat << eof')
vs = ['x', 'd_{cg}', 'd_{c}', 'e_{c}', 'd_{cm1}', 'e_{cg}', 'e_{cm1}']
f.puts(vs.join("\t"))
c.times \
{
	|x|
	de = x.to_f / (c - 1)
	
	pd = dense(w, de)
	n = terras121(pd)

	y = count('d', n)
	
	pe = entropy(w, (de * w).to_i)
	n = terras121(pe)
	
	y.merge!(count('e', n))
	
	f.puts(([x] + y.values_at(*vs)).join("\t"))
}
f.puts('eof')
f.puts('set ytics nomirror; set y2tics; plot \\')
vs[1..-1].each_with_index \
{
	|x, i|
	ax = ''
	ax = 'axes x1y2' if ((3..5).member?(i))
	f.puts("$dat using (#{i * c} + column('x')):(column('#{x}')) with line title '#{x}' #{ax},\\")
}
f.close