Birthday paradox

birthday_paradox.py

import numpy as np
from tqdm import tqdm
import matplotlib.pyplot as plt

res = []

N = 1000
Bs = np.ceil(np.logspace(1, 6, 24)).astype(np.int32)

for B in tqdm(Bs):
    n = 2
    while True:
        d = 0

        xs = np.sort(np.random.randint(B, size=(N, n)), axis=1)
        d = sum(map(lambda r: r[:-1][r[1:] == r[:-1]].size > 0, xs))

        if (1.0 * d / N) >= 0.5:
            break
        n += 1

    res += [(B, n)]


x, y = zip(*res)

plt.semilogx(x, y)
plt.grid()
plt.show()