That looks reasonable. The operation you are implementing is known as 'convolution' and is equivalent to multiplying polynomials. It would be a little more general if you had the input 'die' be a sequence of the count for each outcome, so d6 would be [1]*6 (or [0]+[1]*6 if you prefer). That would allow allow you to represent unfair dice and also to add not just a die to a distribution, but to add any two distributions, so you can play tricks like computing 16d6 as (d6)*2*2*2*2. (The code above is a convolution that restricts the second distribution 'die' to have only 0 and 1 coefficients.) The general convolution can be implemented much like what you have, except that you need another multiplication (to account for the fact that the coefficient is not always 0 or 1). My not particularly efficient implementation:
def vsum(seq1, seq2): return [ a + b for a, b in zip(seq1, seq2) ] def vmul(s, seq): return [ s * a for a in seq ] # Convolve 2 sequences # equivalent to adding 2 probabililty distributions def conv(seq1, seq2): n = (len(seq1) + len(seq2) -1) ans = [ 0 ] * n for i, v in enumerate(seq2): vec = [ 0 ] * i + vmul(v, seq1) + [ 0 ] * (n - i - len(seq1)) ans = vsum(ans, vec) return ans # Convolve a sequence n times with itself # equivalent to multiplying distribution by n def nconv(n, seq): ans = seq for i in range(n-1): ans = conv(ans, seq) return ans print nconv(3, [ 1 ] * 6) print nconv(3, [ 1.0/6 ] * 6) print nconv(2, [ .5, .3, .2 ]) [1, 3, 6, 10, 15, 21, 25, 27, 27, 25, 21, 15, 10, 6, 3, 1] [0.0046296296296296294, 0.013888888888888888, 0.027777777777777776, 0.046296296296296294, 0.069444444444444448, 0.097222222222222238, 0.11574074074074074, 0.125, 0.125, 0.11574074074074073, 0.097222222222222224, 0.069444444444444448, 0.046296296296296294, 0.027777777777777776, 0.013888888888888888, 0.0046296296296296294] [0.25, 0.29999999999999999, 0.29000000000000004, 0.12, 0.040000000000000008] -- http://mail.python.org/mailman/listinfo/python-list