quolc

# efficient with-replacement-sampling algorithm
# cf. http://blog.unnono.net/2010/07/blog-post.html

probs = [0.2, 0.3, 0.4, 0.1]
n = len(probs)

# pre-computation of a[] and p[]
a = [-1]*n
p = [0]*n
fixed = 0