How to calculate the probability of getting the sum X using N six-sided dice

Viewed 952

The Challenge: For example, what is the probability of getting the sum of 15 when using 3 six-sided dice. This can be for example by getting 5-5-5 or 6-6-3 or 3-6-6 or many more options.

A brute force solution for 2 dice - with complexity of 6^2:

Assuming we had only 2 six-sided dice, we can write a very basic code like that:

public static void main(String[] args) {
   System.out.println(whatAreTheOdds(7));
}

public static double whatAreTheOdds(int wantedSum){
    if (wantedSum < 2 || wantedSum > 12){
        return 0;
    }

    int wantedFound = 0;
    int totalOptions = 36;

    for (int i = 1; i <= 6; i++) {
        for (int j = 1; j <= 6; j++) {
            int sum = i+j;
            if (sum == wantedSum){
                System.out.println("match: " + i  + " " + j );
                wantedFound +=1;
            }
        }
    }

    System.out.println("combinations count:" + wantedFound);
    return (double)wantedFound / totalOptions;
}

And the output for 7 will be:

match: 1 6

match: 2 5

match: 3 4

match: 4 3

match: 5 2

match: 6 1

combination count:6

0.16666666666666666

The question is how to generalize the algorithm to support N dice:

public static double whatAreTheOdds(int wantedSum, int numberOfDices)

Because we can't dynamically create nested for loops, we must come with a different approach.

I thought of something like that:

 public static double whatAreTheOdds(int sum, int numberOfDices){

    int sum;
    for (int i = 0; i < numberOfDices; i++) {
        for (int j = 1; j <= 6; j++) {

        }
    }
}

but failed to come up with the right algorithm.

Another challenge here is - is there a way to do it efficiently, and not in a complexity of 6^N?

4 Answers

Here is a recursive solution with memoization to count the combinations.

import java.util.Arrays;
import java.lang.Math;

class Dices {
    public static final int DICE_FACES = 6;

    public static void main(String[] args) {
        System.out.println(whatAreTheOdds(40, 10));
    }

    public static double whatAreTheOdds(int sum, int dices) {
        if (dices < 1 || sum < dices || sum > DICE_FACES * dices) return 0;

        long[][] mem = new long[dices][sum];
        for (long[] mi : mem) {
            Arrays.fill(mi, 0L);
        }
        long n = whatAreTheOddsRec(sum, dices, mem);
        return n / Math.pow(DICE_FACES, dices);
    }

    private static long whatAreTheOddsRec(int sum, int dices, long[][] mem) {
        if (dices <= 1) {
            return 1;
        }
        long n = 0;
        int dicesRem = dices - 1;
        int minFace = Math.max(sum - DICE_FACES * dicesRem, 1);
        int maxFace = Math.min(sum - dicesRem, DICE_FACES);
        for (int i = minFace; i <= maxFace; i++) {
            int sumRem = sum - i;
            long ni = mem[dicesRem][sumRem];
            if (ni <= 0) {
                ni = whatAreTheOddsRec(sumRem, dicesRem, mem);
                mem[dicesRem][sumRem] = ni;
            }
            n += ni;
        }
        return n;
    }
}

Output:

0.048464367913724195

EDIT: For the record, the complexity of this algorithm is still O(6^n), this answer just aims to give a possible implementation for the general case that is better than the simplest implementation, using memoization and search space prunning (exploring only feasible solutions).

As Alex's answer notes, there is a combinatorial formula for this:

formula

In this formula, p is the sum of the numbers rolled (X in your question), n is the number of dice, and s is the number of sides each dice has (6 in your question). Whether the binomial coefficients are evaluated using loops, or precomputed using Pascal's triangle, either way the time complexity is O(n2) if we take s = 6 to be a constant and X - n to be O(n).


Here is an alternative algorithm, which computes all of the probabilities at once. The idea is to use discrete convolution to compute the distribution of the sum of two random variables given their distributions. By using a divide and conquer approach as in the exponentiation by squaring algorithm, we only have to do O(log n) convolutions.

The pseudocode is below; sum_distribution(v, n) returns an array where the value at index X - n is the number of combinations where the sum of n dice rolls is X.

// for exact results using integers, let v = [1, 1, 1, 1, 1, 1]
// and divide the result through by 6^n afterwards
let v = [1/6.0, 1/6.0, 1/6.0, 1/6.0, 1/6.0, 1/6.0]

sum_distribution(distribution, n)
    if n == 0
        return [1]
    else if n == 1
        return v
    else
        let r = convolve(distribution, distribution)
        // the division here rounds down
        let d = sum_distribution(r, n / 2)
        if n is even
            return d
        else
            return convolve(d, v)

Convolution cannot be done in linear time, so the running time is dominated by the last convolution on two arrays of length 3n, since the other convolutions are on sufficiently shorter arrays.

This means if you use a simple convolution algorithm, it should take O(n2) time to compute all of the probabilities, and if you use a fast Fourier transform then it should take O(n log n) time.

You might want to take a look at Wolfram article for a completely different approach, which calculates the desired probability with a single loop.

The idea is to have an array storing the current "state" of each dice, starting will every dice at one, and count upwards. For example, with three dice you would generate the combinations:

111
112
...
116
121
122
...
126
...
665
666

Once you have a state, you can easily find if the sum is the one you are looking for.

I leave the details to you, as it seems a useful learning exercise :)

Related