In [1]:
from __future__ import print_function
from collections import defaultdict
from math import factorial
import sdm as sdmlib
In [2]:
bits = 1000
sample = 1000000
radius = 451
In [3]:
_comb_cache = {}
def comb(a, b):
if a < 0:
return 0
if b == 0:
return 1
ret = _comb_cache.get((a, b), None)
if ret is None:
ret = comb(a-1, b) + comb(a-1, b-1)
_comb_cache[(a, b)] = ret
return ret
def comb(a, b):
return factorial(a) // factorial(b) // factorial(a-b)
In [38]:
def genM(n, h):
M = defaultdict(int)
for i in range(n):
x = h+i
y = i
if 0 <= x < n and 0 <= y < n:
M[(x, y)] = comb(n-h, i)
x = h-i
y = i
if 0 <= x < n and 0 <= y < n:
M[(x, y)] = comb(h, i)
for j in range(1, h+1):
mult = M[(h-j, j)]
for i in range(1, n):
x = h+i
y = i
nx = x - j
ny = y + j
if 0 <= x < n and 0 <= y < n and 0 <= nx < n and 0 <= ny < n:
#print('i={} j={} ({}, {}) ({}, {})'.format(i, j, x, y, nx, ny))
M[(nx, ny)] = mult * M[(x, y)]
return M
In [35]:
def printM(M):
for i in range(n):
for j in range(n):
y = n-1-i
x = j
if (M[(x, y)]):
print('{:3d}'.format(M[(x, y)]), end=' ')
else:
print(' ', end=' ')
print('')
In [39]:
def phi2_fn(bits, sample, radius, d):
total = 0
M = genM(bits, d)
for i in range(radius+1):
for j in range(radius+1):
total += M[(i, j)]
return 1.0 * sample * total / (2**bits)
In [7]:
_phi_fn_cache = {}
In [8]:
def phi_fn(n, H, r, d, steps=500):
key = (n, H, r, d, steps)
if key in _phi_fn_cache:
return _phi_fn_cache[key]
v = []
for _ in range(steps):
bs1 = sdmlib.Bitstring.init_random(n)
bs2 = bs1.copy()
bs2.flip_random_bits(d)
selected1 = address_space.scan_thread2(bs1, r)
selected2 = address_space.scan_thread2(bs2, r)
x = len(set(selected1) & set(selected2))
v.append(x)
mu = 1.0*sum(v)/len(v)
_phi_fn_cache[key] = mu
return mu
In [10]:
address_space = sdmlib.AddressSpace.init_random(bits, sample)
In [40]:
h = 102
a = phi_fn(bits, sample, radius, h, steps=200)
b = phi2_fn(bits, sample, radius, h)
(a, b, a-b)
Out[40]:
In [15]:
Out[15]:
In [17]:
for h in range(0, 1001):
a = phi_fn(bits, sample, radius, h, steps=20)
b = phi2_fn(bits, sample, radius, h)
print('{} {} {}'.format(a, b, abs(a-b)))
In [ ]: