Hi,

I've been struggling with this for awhile ( http://en.wikipedia.org/wiki/Banach%27s_matchbox_problem)

and I put together this little bit of Python code

```
def foo(n=3):
from numpy import rand
x = [1,]*n # one box of matches with elements 1
y = [0,]*n # the other box with elements 0
c=[]
while x and y:
if rand()<0.5: # 1/2 probability of picking either box
c.append(x.pop())
else:
c.append(y.pop())
if x: return len(x) # num matches in non-empty box
if y: return len(y)
```

which simulates the act of picking a match from one box or the other ($1$ or $0$) for this code. The problem is that I can't get the analytical solution for the probability of having $k$ matches left in the remaining box:

$ 1/2^\{(2 n -k )\} Binomial(2 n -k, n) $

For example, if I run 1000 cases for the case where $n=3$ as in

```
>>> y = [ foo() for i in range ( 1000)]
```

I get 1 match left in the remaining box 382 times, 2 matches left 375 times, and 3 matches left 243 times. Thus, my estimate of the probability of getting 1 match left in the remaining box is 382/1000, but the analytical solution is

$ 1/2^\{(2 * 3 -1 )\} Binomial(2 *3 -1, 3)= 0.3125$ vs. 0.382

$ 1/2^\{(2 *3 -2 )\} Binomial(2 *3 -2, 3) = 0.25$ vs. 0.375

$ 1/2^\{(2 *3 -3 )\} Binomial(2 *3 -3, 3) = 0.125$ vs. 0.243

which is nowhere near my estimates, and I have tried this with 10e6 trials without getting any closer.

Any help appreciated.