Number of combinations (N choose R) in C++

Question

Here I try to write a program in C++ to find NCR. But I've got a problem in the result. It is not correct. Can you help me find what the mistake is in the program?

#include <iostream>
using namespace std;
int fact(int n){
    if(n==0) return 1;
    if (n>0) return n*fact(n-1);
};

int NCR(int n,int r){
    if(n==r) return 1;
    if (r==0&&n!=0) return 1;
    else return (n*fact(n-1))/fact(n-1)*fact(n-r);
};

int main(){
    int n;  //cout<<"Enter A Digit for n";
    cin>>n;
    int r;
         //cout<<"Enter A Digit for r";
    cin>>r;
    int result=NCR(n,r);
    cout<<result;
    return 0;
}

Answer 1

Your formula is totally wrong, it's supposed to be fact(n)/fact(r)/fact(n-r), but that is in turn a very inefficient way to compute it.

See Fast computation of multi-category number of combinations and especially my comments on that question. (Oh, and please reopen that question also so I can answer it properly)

The single-split case is actually very easy to handle:

unsigned nChoosek( unsigned n, unsigned k )
{
    if (k > n) return 0;
    if (k * 2 > n) k = n-k;
    if (k == 0) return 1;

    int result = n;
    for( int i = 2; i <= k; ++i ) {
        result *= (n-i+1);
        result /= i;
    }
    return result;
}

Demo: http://ideone.com/aDJXNO

If the result doesn't fit, you can calculate the sum of logarithms and get the number of combinations inexactly as a double. Or use an arbitrary-precision integer library.

---

I'm putting my solution to the other, closely related question here, because ideone.com has been losing code snippets lately, and the other question is still closed to new answers.

#include <utility>
#include <vector>

std::vector< std::pair<int, int> > factor_table;
void fill_sieve( int n )
{
    factor_table.resize(n+1);
    for( int i = 1; i <= n; ++i )
        factor_table[i] = std::pair<int, int>(i, 1);
    for( int j = 2, j2 = 4; j2 <= n; (j2 += j), (j2 += ++j) ) {
        if (factor_table[j].second == 1) {
            int i = j;
            int ij = j2;
            while (ij <= n) {
                factor_table[ij] = std::pair<int, int>(j, i);
                ++i;
                ij += j;
            }
        }
    }
}

std::vector<unsigned> powers;

template<int dir>
void factor( int num )
{
    while (num != 1) {
        powers[factor_table[num].first] += dir;
        num = factor_table[num].second;
    }
}

template<unsigned N>
void calc_combinations(unsigned (&bin_sizes)[N])
{
    using std::swap;

    powers.resize(0);
    if (N < 2) return;

    unsigned& largest = bin_sizes[0];
    size_t sum = largest;
    for( int bin = 1; bin < N; ++bin ) {
        unsigned& this_bin = bin_sizes[bin];
        sum += this_bin;
        if (this_bin > largest) swap(this_bin, largest);
    }
    fill_sieve(sum);

    powers.resize(sum+1);
    for( unsigned i = largest + 1; i <= sum; ++i ) factor<+1>(i);
    for( unsigned bin = 1; bin < N; ++bin )
        for( unsigned j = 2; j <= bin_sizes[bin]; ++j ) factor<-1>(j);
}

#include <iostream>
#include <cmath>
int main(void)
{
    unsigned bin_sizes[] = { 8, 1, 18, 19, 10, 10, 7, 18, 7, 2, 16, 8, 5, 8, 2, 3, 19, 19, 12, 1, 5, 7, 16, 0, 1, 3, 13, 15, 13, 9, 11, 6, 15, 4, 14, 4, 7, 13, 16, 2, 19, 16, 10, 9, 9, 6, 10, 10, 16, 16 };
    calc_combinations(bin_sizes);
    char* sep = "";
    for( unsigned i = 0; i < powers.size(); ++i ) {
        if (powers[i]) {
            std::cout << sep << i;
            sep = " * ";
            if (powers[i] > 1)
                std::cout << "**" << powers[i];
        }
    }
    std::cout << "\n\n";
}

Answer 2

The definition of N choose R is to compute the two products and divide one with the other,

(N * N-1 * N-2 * ... * N-R+1) / (1 * 2 * 3 * ... * R)

However, the multiplications may become too large really quick and overflow existing data type. The implementation trick is to reorder the multiplication and divisions as,

(N)/1 * (N-1)/2 * (N-2)/3 * ... * (N-R+1)/R

It's guaranteed that at each step the results is divisible (for n continuous numbers, one of them must be divisible by n, so is the product of these numbers).

For example, for N choose 3, at least one of the N, N-1, N-2 will be a multiple of 3, and for N choose 4, at least one of N, N-1, N-2, N-3 will be a multiple of 4.

C++ code given below.

int NCR(int n, int r)
{
    if (r == 0) return 1;

    /*
     Extra computation saving for large R,
     using property:
     N choose R = N choose (N-R)
    */
    if (r > n / 2) return NCR(n, n - r); 

    long res = 1; 

    for (int k = 1; k <= r; ++k)
    {
        res *= n - k + 1;
        res /= k;
    }

    return res;
}