Problem of the Week 1220

Expected Flip-Count

Solution

I received only one solution, from Mark Rickert, as below. The key point is that the number can be computed from the set of all 281,589 integer partitions of 52. The answer, as obtained both my Mark and me, is:

2897286937221611262532344522638228975976052384126419605213276991007651

---

190646595858594622078470596206045267396138831042087747584000000000000

This is about 15.1972.

Richard Stong proved that the asymptotic behavior of this function is linear: c(n + 1)/2, where c is given by this integral of an integral:

1 − Integrate[ x E^ (Ei(-x) - x) , {x, 0, Infinity}] / 2

Here, E is the base of natural logarithms, and Ei is the standard ExponentialIntegral function.

One can evaluate this in Mathematica as follows:

1 - (1/2)*NIntegrate[x E ^ (ExpIntegralEi[-x] - x), {x, 0, Infinity}]

For n = 52, this gives 15.1925, an error of 0.0003 from the true value above. Nice!

Those without Mathematica can get these results by entering this into wolframalpha.com:

1 - (1/2)*Integrate[x E ^ (ExpIntegralEi[-x] - x), {x, 0, Infinity}]

Mark's Solution

For 1220, I get 15.1971606111 for 52 cards and 2.60905877976 for 8 cards.

For each partition of 52, I compute how many permutations fall into that category, split the largest cycle, and compute the average.

I used the following Python program that took about 30 seconds for 52 cards:

from math import factorial

'''
Partition function courtesy of Philip Pham:
https://code.activestate.com/recipes/
        218332-generator-for-integer-partitions/
'''
def P(n):
    # base case of recursion: zero is the sum of the empty list
    if n == 0:
    yield []
    return

    # modify partitions of n-1 to form partitions of n
    for p in P(n-1):        
    p.append(1)
    yield p
    p.pop()
    if p and (len(p) < 2 or p[-2] > p[-1]):
      p[-1] += 1
        yield p

N = 52
total = 0.0
for p in P(N):
    sub = factorial(N)
    for i in p:
    #~ sub /= factorial(i)
    #~ sub *= factorial(i-1)
    sub /= i
    for i in xrange(1, N+1):
    c = p.count(i)
    if c > 1:
      sub /= factorial(c)
    for i in xrange(len(p)):
    if i == 0 and p[i] > 1:
      total += ((p[i]/2)**2) / float(N) * sub
        total += (((p[i]+1)/2)**2) / float(N) * sub
	else:
	  total += (p[i]**2) / float(N) * sub
    #~ print p, sub, total
print total / factorial(N)

Here is my Mathematica code that computes the same thing. This code is as short as I can get it, and can just about be called a one-liner. As with Mark's method, it starts with the set of all integer partitions of 52, and uses some delicate work to quickly get the number of permutations having cycle-lengths equal to each of the partitions. If you want Mathematica code that is easier to understand, just ask me. The part with the Quotient is where the largest entry in the partition gets cut in half. This takes 4.5 seconds for the n = 52 case.

n = 52;
(1/n) Sum[Total[Join[Rest[p], Quotient[p[[1]] + {0, 1}, 2]]^2] /
          ((Times @@ p ) * (Times @@ ((Length /@ Split[p])!))),
       {p, IntegerPartitions[n]}]
       
2897286937221611262532344522638228975976052384126419605213276991007651/
   190646595858594622078470596206045267396138831042087747584000000000000

As regards the harder problem I set, to find a strategy that beats the cycle-shortening strategy: at the time I posed it, I thought I had a better strategy. Then I found an error. Working with Larry Carter and Mark Rickert, we have indeed found a strategy that appears as if it will always have a smaller expected flip-count than the basic strategy, except for 2 and 3, in which case it ties with the basic strategy. Details will have to wait until the dust settles.

[Back to Problem 1220]

© Copyright 2016 Stan Wagon. Reproduced with permission.