Project Euler Problem 740

Secret Santa is a process that allows n people to give each other presents, so that each person gives a single present a

Project Euler Problem 740

Solution

Answer: 0.0189581208

Let the state before a turn be described by:

  • $c_0$: remaining people with 0 copies of their own name still in the hat,
  • $c_1$: remaining people with 1 copy,
  • $c_2$: remaining people with 2 copies.

Initially the state is $(0,0,n)$.

If only one person remains, the process fails unless that person has no copies of their own name left:

$$q = \begin{cases} 0 & (c_0,c_1,c_2)=(1,0,0)\ 1 & \text{otherwise} \end{cases}$$

We use memoized recursion over states.

For a chosen current person of type $t\in{0,1,2}$:

  • there are $2m-t$ drawable slips,
  • where $m=c_0+c_1+c_2$.

We enumerate all ways the two drawn slips can affect future players’ remaining self-slips and recurse.

The following Python program computes the probability exactly:

from functools import lru_cache
from math import comb

@lru_cache(None)
def fail(c0, c1, c2):
    m = c0 + c1 + c2

    # Last remaining person
    if m == 1:
        return 0 if c0 == 1 else 1

    ans = 0.0

    for t, ct in [(0, c0), (1, c1), (2, c2)]:
        if ct == 0:
            continue

        pcur = ct / m

        d0, d1, d2 = c0, c1, c2
        if t == 0:
            d0 -= 1
        elif t == 1:
            d1 -= 1
        else:
            d2 -= 1

        D = 2 * m - t
        R = d1 + 2 * d2
        I = D - R

        den = comb(D, 2)

        sub = 0.0

        # draw two irrelevant slips
        sub += comb(I, 2) / den * fail(d0, d1, d2)

        # draw one relevant slip
        if d2:
            sub += (2 * d2 * I) / den * fail(d0, d1 + 1, d2 - 1)

        if d1:
            sub += (d1 * I) / den * fail(d0 + 1, d1 - 1, d2)

        # draw two relevant slips
        if d2:
            sub += d2 / den * fail(d0 + 1, d1, d2 - 1)

        if d2 >= 2:
            sub += 4 * comb(d2, 2) / den * fail(d0, d1 + 2, d2 - 2)

        if d2 and d1:
            sub += 2 * d2 * d1 / den * fail(d0 + 1, d1, d2 - 1)

        if d1 >= 2:
            sub += comb(d1, 2) / den * fail(d0 + 2, d1 - 2, d2)

        ans += pcur * sub

    return ans

print(fail(0, 0, 100))

This reproduces the given values:

  • $q(3)=0.3611111111$
  • $q(5)=0.2476095994$

and computes:

$$q(100)=0.018958120754006404\ldots$$

Rounded to 10 decimal places:

Answer: 0.0189581208