# How to calculate probability distribution of rolling n dice? (with a twist!)

I'm trying to work on a project in code that asks this question. I'm currently able to get the right answer by breaking it up, but I'm wondering if there's a mathematical formula I can address this problem to as I'm actually facing computational resource limitations.

For example, the amount of combinations for 50 dice is 50^6, which is 15,625,000,000 which is something my computer can't handle.

I'm trying to find the probability distribution of getting a certain value of n dice where no single value of any dice is 1.

For example, with a single dice, the probability of getting a 1,2,3,4,5, or 6 is 1/6.

For two dice, the probability of getting a total value of 4 or 12 is 1/36 (I ignore the case of 2 and 3 since one of the dice has to have a value of 1). The probability of a 5 or 11 WITHOUT a dice having a value of one is 1/18, 6 or 10 WITHOUT a dice having a value of 1 is 1/12, etc. etc.

Is anyone familiar with generating the probability distribution that all the possible values that can be rolled with n dice such that no dice contains a 1?

• You have a poor computer indeed if it cannot handle an 11 digit number! Since most computers are binary, you might find it convenient to represent $50=10^2/2$, for then $$50^6=10^{12}/2^6=0.015625\times 10^{12}$$ shows that you only need to deal with a six-digit number separately from the exponent. Regardless, there are many solutions, including those given at stats.stackexchange.com/questions/3614 and stats.stackexchange.com/questions/24385.
– whuber
Mar 11 '18 at 14:38

Solve the problem for fifty 5-sided dices (and add 50 to the dice rolls).

You could do the first part (the standard dice roll problem) by

• approximation with a normal distribution ( https://math.stackexchange.com/questions/406192/probability-distribution-of-rolling-multiple-dice )
• by the formula for the exact solution from Whuber mentioned in the comments.

• or by computation doing a convolution 49 times

rolls <- rep(0,5*50)

rolls[1:5] <- 1

for (i in 1:49) {
rolls <- c(0,rolls[1:245],0,0,0,0) +
c(0,0,rolls[1:245],0,0,0) +
c(0,0,0,rolls[1:245],0,0) +
c(0,0,0,0,rolls[1:245],0) +
c(0,0,0,0,0,rolls[1:245])
}


The probability for 5-sided dice roll being x is rolls[x]/5^50. The probability for a 6-sided dice roll being x (and none of the dices is a one) is rolls[x-50]/5^60. Or rolls[x-50]/5^50 if you wish to condition on none of the dices being 1.

The above code can be improved by using 128 bit integers (if you find that necessary)

• To what, exactly will you add 50? Certainly not to the probability. Given that the chance of the event in question is $(5/6)^{50}\approx 0.00011,$ there doesn't seem to be anything in the output of your code that directly corresponds to it. A little more explanation will help readers understand this solution.
– whuber
Mar 11 '18 at 14:33
• I intended to add the 50 to the variable. The probability for X=x is rolls[x-50]/6^50. But I actually didn't wanted to give away too much, create a completely finished and polished solution, before it was clear that this is not a self-study question. Mar 11 '18 at 15:00