Using Python to simulate Kruskal count and estimating the success probability

Introduction

If you have no idea what Kruskal count is, please refer to: A simple to learn magic with amazing effect, The Kruskal Count. I will be using Python to simulate the trick. By doing the trick many times, and counting the number of success over all the trials, I can estimate the success probabiliy.
 

The Procedure

To calculate the success probability, we will first generate a random permutation of 52 cards. For each of the 1 to 10 secret number, we will do a Kruskal count, and see which card is the last selected card. For instance:

How do we calculate the success probability in this case?

The Code

kruskal.py

import random

def shuffle(x):
    z = []
    y = x[:]
    for i in reversed(range(len(x))):
        z.append(y.pop(random.randint(0, i)))
    return z

def incremental(x):
    r = x % 13
    if 12 >= r >= 10:
        return 5
    else:
        return r + 1

def percent(x):
    a = {}
    for j in x:
        if j in a:
            a[j]+=1
        else:
            a[j]= 1
    p = 0
    for y in a.values():
        p += y**2
    return p / float(len(x)**2)

def mean(x):
        total = 0
        for i in x:
                total += i
        return total / float(len(x))

def main(n):
    pl = []
    unsuffled_decks = range(52)
    r10 = range(10)
    results = [-1 for x in r10]
    for k in range(n):
        decks = shuffle(unsuffled_decks)
        for i in r10:
            c = i
            while c < 52:
                s = decks[c]
                c += incremental(s)
            results[i] = s
        if (k % 10000 == 0):
            print 'pass: %i' % k
        pl.append(percent(results))
    return pl

print mean(main(1000000))

I will explain what each of the function is doing.

shuffle(x):

def shuffle(x):
    z = []
    y = x[:]
    for i in reversed(range(len(x))):
        z.append(y.pop(random.randint(0, i)))
    return z

The argument x is a list of different objects. And the function shuffle(x) will return a new list, which is a random permutation of the list x. Please notice that though in random package, there is a function called random.shuffle. We won't be using it since according to the documentation, even for small size x, most permutations can never be generated.

Basically, we first copy the list x to list y, so that we can change the list y without affecting the original list. Then we will randomly pick an element from y, and remove it using y.pop function, and put it in a new list z. We continue until all the element in y is removed, and added to z. z is a random permutation of the original list x.

incremental(x):

def incremental(x):
    r = x % 13
    if 12 >= r >= 10:
        return 5
    else:
        return r + 1

To explain this function, I need to explain how we are going to represent the deck in the program. We will be using [0,1,2,..,51] to represent the whole 52 cards-deck, and each card is represented by a number ranged from 0 to 51. For instance:

• ♠ 1 => 0
• ♠ 2 => 1
  
• ....
• ♠ K => 12
  
• ♥ 1 => 13
  
• ....

The suit is not important in this case, since the we only need the face value of the card. As mentioned in A simple to learn magic with amazing effect, The Kruskal Count, all picture cards have value of 5, and other have their face value. The incremental(x) function is used to calculate the value of each card. We first do a modulo 13 of the card value x. If it's value is from 0 to 9, it's face value is the value plus 1. If the modulo value is 10 to 12, it means it is a picture cards, and should have a value of 5.

main(n):

def main(n):
    pl = []
    unsuffled_decks = range(52)
    r10 = range(10)
    results = [-1 for x in r10]
    for k in range(n):
        decks = shuffle(unsuffled_decks)
        for i in r10:
            c = i
            while c < 52:
                s = decks[c]
                c += incremental(s)
            results[i] = s
        if (k % 10000 == 0):
            print 'pass: %i' % k
        pl.append(percent(results))
    return pl

In this function, it will run the Kruskal Count and compute the success probability for n times, and then return a list of the success probability.

The variable pl is used to store the probability list. And the variable results is a list of 10 items, and they stored the ending card for the choice of 1 to 10 secret numbers.

The main loop:

        decks = shuffle(unsuffled_decks)
        for i in r10:
            c = i
            while c < 52:
                s = decks[c]
                c += incremental(s)
            results[i] = s

First of all, we shuffle the unsuffled_decks using the shuffe(x) function. Then for each of the secret number from 1 to 10, we set the current index c to the secret number first. Then while the current index c doesn't excess the number of card (since we counted from zero, so the condition is less than 52),  we put the current card decks[c] to the variable s. Then c is increased by the value of the card computed by the incremental(s) function. Once the current index excess 52, the loop ended, and we put the final card s into the results list.

After we get the final cards for all secret numbers, I will compute the success probability using percent(results), and append the probability to the pl list.

Since we are computing the probability by simulation, we need to repeat the experiment many times, in order to obtain a reliable result.

Result

You can copy and paste the code to a file named kruskal.py and run it. Here is the output:

pass: 0
pass: 10000
pass: 20000
pass: 30000
pass: 40000
pass: 50000
...
pass: 970000
pass: 980000
pass: 990000
0.84181804

So the success probability is found to be 0.84181804. So, if you as magician chooses a secret number randomly and do the Kruskal Count as the audience, you will have 84.18% of ending up with the same ending card as the audience.

If you wanted to calculated the result instead of doing simulation, you can refer to another post of mine which detailed the formulation of the problem, and how to solve it using programming technique.