A Specific Combination

Question

This code is going to be included in my Minesweeper Probabilities project. (so far it only exists on the develop-branch)

The purpose is to return a specific combination, let's say for example that you have 7 elements and want to return 3 of them, simple combinatorics tells us that there are \$\binom{7}{3} = \frac{7*6*5}{3*2*1} = 35\$ combinations to do this, so let's say that we want to return an arbitrary combination, let's say combination number \$20\$.

So we have to decide the order for the combinations. Let's take the elements [0 1 2 3 4 5 6], the order that felt the most natural to me is:

012 013 014 015 016
023 024 025 026
034 035 036
045 046
056

123 124 125 126
134 135 136
145 146
156

234 235...

So when we have written these, we can see that 134 is the 20th combination.

---

But without actually writing them all down, how can we come to the conclusion that 134 is the 20th combination?

There are \$\binom{6}{2} = 15\$ combinations where \$0\$ is the first number. As \$20 > 15\$, we know that 0 is not the first number. So we reduce 20 by 15, increase the nextNumber variable and continue.

Then we have nextNumber == 2 and combination == 5, as there are \$\binom{5}{2}\$ combinations where \$1\$ is the first number, we know that 1 is the first number, reduce the remainingSize and continue on our way...

Then there are 4 combinations were 2 is the next number [123, 125, 126, 127], but as 5 > 4, 2 is not the next number. So again we reduce the combination we're looking for (5) by the number of combinations for 2 being the next number (4), so we get combination = 1 and also increase nextNumber to 3.

Now, as combination == 1 the remaining numbers are 3 and 4.

I hope the code is pretty self-explanatory.

---

Code

As I need to support very big combinations, such as \$\binom{256}{51} = 2.034*10^{54}\$, I assumed that double will not have enough precision, so I went with BigInteger. I hope this doesn't cause too much impact on performance though.

Can also be found on github.

public static BigInteger nCrBigInt(int n, int r) {
    if (r > n || r < 0) {
        return BigInteger.ZERO;
    }
    if (r == 0 || r == n) {
        return BigInteger.ONE;
    }
    if (r > n / 2) {
        // As Pascal's triangle is horizontally symmetric, use that property to reduce the for-loop below
        r = n - r;
    }

    BigInteger value = BigInteger.ONE;

    for (int i = 0; i < r; i++) {
        value = value.multiply(BigInteger.valueOf(n - i)).divide(BigInteger.valueOf(i + 1));
    }

    return value;
}

public static int[] specificCombination(final int elements, final int size, final BigInteger combinationNumber) {
    if (combinationNumber.signum() != 1) {
        throw new IllegalArgumentException("Combination must be positive");
    }
    if (elements < 0 || size < 0) {
        throw new IllegalArgumentException("Elements and size cannot be negative");
    }

    int[] result = new int[size];

    int resultIndex = 0;
    int nextNumber = 0;
    BigInteger combination = combinationNumber;
    int remainingSize = size;
    int remainingElements = elements;

    while (remainingSize > 0) {
        BigInteger ncr = nCrBigInt(remainingElements - 1, remainingSize - 1);
        if (ncr.signum() == 0) {
            throw new IllegalArgumentException("Combination out of range: " + combinationNumber + " with " + elements + " elements and size " + size);
        }
        if (combination.compareTo(ncr) <= 0) {
            result[resultIndex] = nextNumber;
            remainingSize--;
            resultIndex++;
        }
        else {
            combination = combination.subtract(ncr);
        }
        remainingElements--;
        nextNumber++;
    }

    return result;
}

Test

@Test
public void specificCombinationVeryBig() {
    int[] result = Combinatorics.specificCombination(256, 51, BigInteger.valueOf(Long.MAX_VALUE - 42));
    assertArrayEquals(new int[]{ 0, 1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25,26,27, 28, 29, 30, 31, 32, 33,
            34, 35, 36, 37, 38, 52, 73, 94, 99, 132, 163, 169, 179, 190, 214, 227, 230 }, result);
}

@Test
public void specificCombinations() {
    assertArrayEquals(new int[]{ 0, 1, 2, 3 }, Combinatorics.specificCombination(8, 4, BigInteger.ONE));
    assertArrayEquals(new int[]{ 0, 1, 2, 4 }, Combinatorics.specificCombination(8, 4, BigInteger.valueOf(2)));
    assertArrayEquals(new int[]{ 0, 2, 4, 5 }, Combinatorics.specificCombination(8, 4, BigInteger.valueOf(20)));

    assertArrayEquals(new int[]{ 1, 4, 6 }, Combinatorics.specificCombination(7, 3, BigInteger.valueOf(24)));
    assertArrayEquals(new int[]{ 1, 5, 6 }, Combinatorics.specificCombination(7, 3, BigInteger.valueOf(25)));
    assertArrayEquals(new int[]{ }, Combinatorics.specificCombination(7, 0, BigInteger.ONE));
}

@Test(expected = IllegalArgumentException.class)
public void specificCombinationOutOfRange() {
    Combinatorics.specificCombination(7, 3, BigInteger.valueOf(36));
}

@Test(expected = IllegalArgumentException.class)
public void specificCombinationZeroOrLess() {
    Combinatorics.specificCombination(7, 3, BigInteger.ZERO);
}

@Test(expected = IllegalArgumentException.class)
public void specificCombinationWithNegativeSize() {
    Combinatorics.specificCombination(7, -1, BigInteger.ONE);
}

@Test(expected = IllegalArgumentException.class)
public void specificCombinationWithTooBigSize() {
    Combinatorics.specificCombination(7, 8, BigInteger.ONE);
}

@Test(expected = IllegalArgumentException.class)
public void specificCombinationWithNegativeElements() {
    Combinatorics.specificCombination(-1, 3, BigInteger.ONE);
}

@Test
public void nCrBigInt() {
    assertEquals(BigInteger.valueOf(28), Combinatorics.nCrBigInt(8, 2));
    assertEquals(BigInteger.valueOf(28), Combinatorics.nCrBigInt(8, 6));
    assertEquals(BigInteger.valueOf(70), Combinatorics.nCrBigInt(8, 4));
    assertEquals(BigInteger.valueOf(56), Combinatorics.nCrBigInt(8, 3));
    assertEquals(BigInteger.valueOf(35), Combinatorics.nCrBigInt(7, 3));
    assertEquals(BigInteger.ZERO, Combinatorics.nCrBigInt(1, -1));
    assertEquals(BigInteger.ZERO, Combinatorics.nCrBigInt(0, 1));

    for (int i = 0; i < 100; i++) {
        assertEquals(BigInteger.ONE, Combinatorics.nCrBigInt(i, 0));
        assertEquals(BigInteger.ONE, Combinatorics.nCrBigInt(i, i));
    }
}

Primary Concerns

• Is there a better way? Specifically with regards to performance.
• I have added some test cases, are there any edge-cases I have missed? (I don't think there are)

Answer 1

A significant part of the algorithm is to compute the binomial coefficients

while (remainingSize > 0) {
    BigInteger ncr = nCrBigInt(remainingElements - 1, remainingSize - 1);
    // ...
}

However, in each loop iteration, either

• remainingElements is decremented, or
• both remainingElements and remainingSize are decremented,

which means that you can take advantage of the following identities for binomial coefficients:

 C(n-1, k-1) = C(n, k) * k / n
 C(n-1, k)   = C(n, k) * (n - k) / n

and compute the next value of ncr from the previous value, instead of calling nCrBigInt() again. That is one multiplication and division instead of min(k, n-k).

Answer 2

I like the concept, and the code.

Given that your primary concern is about performance, I recommend three things:

1. when the solution fits inside the scope of long, use it instead of BigInteger. Basically, implement the solution once using BigInteger, and again using long. As you reduce the problem in your while loop, switch to the long version.
2. Using your while loop instead of recursion is probably faster, but, as it happens, recursion can have advantages, especially as methods are heavily used, and get compiled, and in-lined. Your algorithm is suitable for recursion, and you can create your result array, and pass it as an argument to the recursive method. I would try to use the recursion, and benchmark to compare.
3. You can cache (or precompute) the results of the nCrBigInt function. A cache with a stepped-size array, of size vs. get, would be easy to populate....

For example, you could have a pair of functions and a data set:

private static BigInteger[][] nCr = new BigInteger[8192][]; // overspecify the size, can be generous....

public static BigInteger nCrBigInt(int n, int r) {
    if (r > n || r < 0) {
        return BigInteger.ZERO;
    }
    if (r == 0 || r == n) {
        return BigInteger.ONE;
    }
    if (r > n / 2) {
        // As Pascal's triangle is horizontally symmetric, use that property to reduce the for-loop below
        r = n - r;
    }
    if (nCr.length <= n) {
        nCr = Arrays.copyOf(nCr, n + 128); // extend the size
    }
    if (nCr[n] == null) {
        nCr[n] = new BigInteger[n / 2 + 1];
    }
    if (nCr[n][r] == null) {
        nCr[n][r] = nCrBigIntCalc(n, r);
    }
    return nCr[n][r];
}

private static BigInteger nCrBigIntCalc(int n, int r) {

    BigInteger value = BigInteger.ONE;

    for (int i = 0; i < r; i++) {
        value = value.multiply(BigInteger.valueOf(n - i)).divide(BigInteger.valueOf(i + 1));
    }

    return value;
}        

If you need this to be concurrent, you can adapt it. The bottom line is that this will, for the most part, remove a log(n) or so order of complexity.

A pre-computed result would be relatively easy to compute, and store as a constant value too.