2

Disclaimer: I'm a computer programmer more than a mathematician, so reading text like that of the answer to this question is a little (read: a lot) over my head.

I've written an algorithm (brute-force, O(nk) running time, where n is the number of dice, and k is the number of sides) that calculates the probability of rolling each possible value of a specified number of dice, each of may have a different number of faces. For example, rolling d8 + 2d20 (1 8-sided die and 2 20-sided dice) has a minimum roll of 3 and a maximum roll of 48, each of which have a probability of 0.0003125. The values 25 and 26 are most likely to show up, each with a probability of 0.045.

I'm wondering if there's a way, without doing the full calculation, to determine--or make a conservative estimate of--this maximum value of 0.045 in polynomial time. I know that this is also the probability of rolling the middle value, or ceiling((max - min) / 2) + min.

I've written a probability distribution graph using jQuery Flot. Since this is such a long running calculation for large numbers of dice with large numbers of sides (e.g. 4d100), I calculate the distribution in chunks, and update the graph periodically. I'd like to have a maximum probability calculated ahead of time so the axes don't change as I make updates to the data.

Calculating the total number of permutations is obviously easy--multiple each die's number of sides together--so calculating the probability of any one permutation is 1 / totalPermutations. That's about as far as I got. Help?

Community
  • 1
dx_over_dt
  • 143
  • Not polynomial in the size of the input, but you could divide and conquer. When you have distributions A and B, the result for the sum is a straightforward convolution, taking time (length A)(length B). Since in practice the lengths are less than nk, one has an algorithm of O(nnkkn) runtime and O(nk) space. –  Jun 19 '15 at 18:11
  • That sounds like that would get you the entire distribution? –  Jun 19 '15 at 18:12
  • Yes, and then run through to find the max. The point is to not calculate distinct throws, but just consider distinct sums. If you have regular dice, you can bin the really low probability results together to save time and space. You can also use theory (multinomial distribution? I'm not sure) to give you a formula. I don't know what your dice probabilities look like, so I am unsure how much theory can help. –  Jun 19 '15 at 18:19
  • If all your dice were two sided and uniform, the answer for n dice would be the central binomial coefficient, which has probability about 1/sqrt(n). If all your cases are regular like your example, you will have a similar max of O(sqrt(n)), and a multinomial distribution will help. I'm hazy on the details though. In any case, you don't need n^k time for this unless k is 4 or less. –  Jun 19 '15 at 18:30
  • This isn't research-level mathematics without some indication that straightforward approaches are not adequate. You can use dynamic programming: Let a[i,j] be the probability that the first j dice have a sum of i. The simple approach to calculating this takes $O(n^2k^2)$ because you calculate nk values n times, and each sum is over k terms. You can speed that up by using the fast Fourier transform, or by using the binomial theorem if dice are repeated. You could estimate the maximum values using a normal approximation instead. – Douglas Zare Jun 19 '15 at 18:44
  • The probability generating function for the sum of $n$ dice of $d_1, \ldots, d_n$ sides respectively is $\prod_{j=1}^n \sum_{i=1}^{d_j} x^i/d_j $. The probability of a sum of $s$ is the coefficient of $x^s$. –  Jun 19 '15 at 18:54
  • @RobertIsrael I could follow (or at least take your word for) all of that except for x... What is x in each iteration? – dx_over_dt Jun 19 '15 at 19:20
  • There's no iteration. The probability generating function is a polynomial in the indeterminate $x$. The coefficients of this polynomial give the probabilities. Maybe I'd better expand this into a solution. – Robert Israel Jun 19 '15 at 21:31

1 Answers1

2

The probability generating function for the sum of $n$ dice of $d_1, \ldots, d_n$ sides respectively is

$$ G(x) = \prod_{j=1}^n \sum_{i=1}^{d_j} \dfrac{x^i}{d_j}$$

The probability of a sum of $s$ is the coefficient of $x^s$ in this polynomial.

You can also write $$\sum_{i=1}^d x^i = x \dfrac{x^{d} - 1}{x - 1}$$

For example, for four dice, each with $100$ sides,

$$ G(x) = 100^{-4} x^4 \left(\dfrac{x^{100} - 1}{x - 1}\right)^{4}$$

Consider the numerator and denominator $$\eqalign{N(x) &= 100^{-4} x^4 (x^{100} - 1)^4 = 100^{-4} (x^{404} - 4 x^{304} + 6 x^{204} - 4 x^{104} + x^4)\cr D(x) &= (x - 1)^4 = x^4 - 4 x^3 + 6 x^2 - 4 x + 1\cr}$$ Now $G(x) = \sum_{s=4}^{400} p_s x^s$ where $p_s$ is the probability of getting a sum of $s$. Since $G(x) D(x) = N(x)$, the coefficient of $x^s$ in $N(x)$ is $$ r_s = p_{s} - 4 p_{s-1} + 6 p_{s-2} - 4 p_{s-3} + p_{s-4}$$ which gives you a recurrence relation you can use to calculate all the $p_s$, starting with $p_0 = \ldots = p_3 = 0$, where $r_s = 0$ except for $$r_4 = r_{404} = 10^{-8}, r_{104} = r_{304} = -4 \times 10^{-8}, r_{204} = 6 \times 10^{-8}$$

Robert Israel
  • 448,999