data31d.rb

def f2(n, o2, e)

	while (n.odd?)
		n = (n * 3 + 1) / 2
		break if (o2)
	end

	while (n.even?) do
		n /= 2 if (n.even?)
		break if (e)
	end
	
	return n

end

def adv(x)
	
	n1 = n = x['n']
	x['t'].times { n = f2(n, x['o2'], x['e']) }
	return d(n.to_s(2))
	
end

def ratio(x)
	
	n1 = n = x['n']
	x['t'].times { n = f2(n, x['o2'], x['e']) }
	return n.to_s(2).length.to_f / n1.to_s(2).length
	
end

def coef(x, y)

        x1 = x2 = y1 = y2 = xy = 0.0

        n = x.size
        n.times \
        {
                |i|

                x1 += x[i]
                y1 += y[i]
                xy += x[i] * y[i]
                x2 += x[i] ** 2
                y2 += y[i] ** 2
        }

        d = n * x2 - x1 ** 2
        b = (y1 * x2 - x1 * xy) / d
        m = (n * xy - x1 * y1) / d

	b = 1 if (b.nan?)
	m = 1 if (m.nan?)
        return b, m
end

def sum(l)

	t = 0.0
	l.each { |x| t += x }
	return t
end

def av(l)
	return nil if (l.empty?)
	return sum(l) / l.size
end

def corr(l, y1, yp)
	
	xav = av(l.map { |x| x[y1] })
	yav = av(l.map { |x| x[yp] })

	tx = ty = txy = e = 0.0
	l.each \
	{
		|z|
		
		x = z[y1]
		y = z[yp]
		
		txy += (x - xav) * (y - yav)
		tx += (x - xav) ** 2
		ty += (y - yav) ** 2
		e += (x - y) ** 2
	}
	r = txy / (Math.sqrt(tx) * Math.sqrt(ty))
	e /= l.size
	return r
	
end

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

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

def test(t, o, o2, e)
	w = 100
	c = 500

	l = []
	(c + 1).times \
	{
		|i|

		d =  i.to_f / c
		n = dense({'d' => d, 'w' => w, 'o' => o})
		a = {'n' => n, 'o2' => o2, 'e' => e, 't' => t}
		d2 = adv(a)
		r = ratio(a)
		l << {'d' => d, 'd2' => d2, 'wr' => r}
	}
	return l
end

def out(fn, a)
	return if (a.nil?)
	f = File.open(fn, 'a')
	f.puts(a.keys.join("\t")) if (f.size == 0)
	f.puts(a.values.join("\t"))
	f.close
end

def out2(fn, l, t)
	f = File.open(fn, 'a')
	l.each { |x| x.merge!({'t' => t}) }
	f.puts(l[0].keys.join("\t")) if (f.size == 0)
	l.each { |x| f.puts(x.values.join("\t")) }
	f.puts
	f.close
end

def fit(l, k, t, a)
	b, m = coef(l.map { |x| x['d'] }, l.map { |x| x[k]})
	l.each { |x| x['y'] = m * x['d'] + b }

	r = corr(l, k, 'y')
	
	a.merge!(k => {'t' => t, 'm' => m, 'b' => b, 'r' => r})
end

def model(t)
	tf = [false, true]

	z= (0...8).map { |x| (1..3).map { |y| tf[x.to_s(2)[-y, 1].to_i] } }

	$stderr.puts(z[6].inspect)
	
	l = test(t, *z[6])

	a = {}
	fit(l, 'd2', t, a)
	fit(l, 'wr', t, a)
	return a
end

def detect(a, d, w, f)
	c = 500
	i = 0
	c.times \
	{
		break if (w < 1.0)
		f.puts([i, d, w].join("\t"))

		d2 = d * a['d2']['m'] + a['d2']['b']
		w2 = d * a['wr']['m'] + a['wr']['b']

		d = d2
		w *= w2
		
		i += 1
	}
	f.puts
	return i == c ? 0 : i
end

def plot(a)
	
	d = 0.0
	c = 100
	f = File.open('out.txt', 'w')
	while (d <= 1.0)
		d += 1.0 / c
		
		w = d * a['wr']['m'] + a['wr']['b']
		d1 = d * a['d2']['m'] + a['d2']['b']
		f.puts([d, (d1 - d), (w - 1.0)].join("\t"))

	end
	f.close
	
end

def plot2(a)
	
	d = 0.0
	c = 10
	f = File.open('out2.txt', 'w')
	while (d <= 1.0)
		d += 1.0 / c
		
		detect(a, d, 3.0, f)

	end
	f.close
	
end

a = model(1)
$stderr.puts(a.inspect)

plot(a)
plot2(a)