I defined a very simple function for the Birthday Problem:
from math import comb, factorial
import numpy as np
def birthday(k):
return 1-((factorial(k)*comb(365,k))/(365**k))
The function works, since k=23
gives:
[in] birthday(23)
[out] 0.5072972343239854
I want to iterate over a list with k
going from 1 to 50 and save the probability outcomes to a list prob
, like this:
klist=np.arange(1,51)
prob=[]
for k in klist:
prob.append(birthday(k))
For k=1 up to k=7, there is no problem, but from k=8
onwards, suddenly the output doesn't make sense anymore with very large negative values. What am I doing wrong here? (showing output for k=1 to 10)
[ 1 2 3 4 5 6 7 8 9 10 ]
[0.0, 0.002739726027397249, 0.008204165884781345, 0.016355912466550326, 0.02713557369979358, 0.040462483649111536, 0.056235703095975365, -203.08817475498518, -20769.916905383445, -11786425.811859423]
However, just running:
[in] birthday(10)
Gives the correct:
[out] 0.11694817771107768
CodePudding user response:
The problem is somehow in using this:
klist=np.arange(1,51)
That makes k an numpy.int32
for every call to birthday
. And if you try this:
print(birthday(np.int32(10)))
You'll find that has the same problem.
This works though:
from math import comb, factorial
import numpy as np
def birthday(k):
return 1-((factorial(k)*comb(365,k))/(365**k))
prob=[]
for k in range(51):
prob.append(birthday(k))
A 32-bit integer appears to be not large enough to hold the values birthday()
computes and it would appear it overflows.