construct45.rb

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

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 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)

	return n
end

def seq(n)

	n1 = n
	l = []
	begin
		l << n
		n = f2(n)
	end while (n != 1)
	
	return l
end

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

def dp(x)
	d(x.to_s(2)) >= 0.5
end

def avg(c)
	return [c['t'], c['c'].to_f / c['t']]
end

def inc(c, x, l, c1)
	c['c'] += x
	c['t'] += 1
	l << [avg(c), c1]
end

	
def seed(d, c)
	
	n = terras121(dense(50, d))

	l = seq(n)
	
	
	l.each_with_index \
	{
		|x, j|
		puts([j, d(n.to_s(2)), c].join("\t"))
		n *= 3 if (x.odd?)
	}
	puts
end

def seed2(d, c)
	
	n = terras121(dense(50, d))

	j = 0
	da = 0.0
	50.times \
	{
		|j|
		da += d(n.to_s(2))
		puts([j, da / (j + 1), c].join("\t"))
		n *= 3
	}
	puts
end

def seed3(d, c, l)
	
	n = terras121(dense(50, d))

	j = t = 0
	da = 0.0
	100.times \
	{
		|j|
		d = d(n.to_s(2))
		t += 1 if (d >= 0.5)
		l << [j, (0.5 - t.to_f / (j + 1)).abs, c]
		n *= 3
	}
	l << []
end

def graph(f, l, t)
	f.puts('$dat << eof')
	l.each { |x| f.puts(x.join("\t")) }
	f.puts('eof')
	f.puts("set title '#{t}'; plot $dat using 1:2:3 with point linecolor palette pt 5 ps 0.75; unset title; pause -1;")
	
end

def range(r, f)
	l = []
	r.each { |x| seed3(x.to_f / 99, (50 - x).abs, l) }
	graph(f, l, "r = #{r.to_s}")
end

def data()
	
	f = File.open('gnuplot.cmd', 'w')
	r1 = (0..99)
	range(r1, f)
	f.close
	
end

data()