How to calculate C(n,k) efficiently
Han de Bruijn
The best algorithm I could find for calculating C(n,k) is based upon
Pascal's triangle:
C(n,k) = C(n-1,k) + C(n-1,k-1) ending with C(n,n) = 1 or C(n,0) = 1 .
Giving the following recursive procedure (in Pascal):
function C(n,r : integer) : integer;
{
Number of Combinations
}
var
u : integer;
begin
u := 1;
if (n > r) and (r > 0) then
u := C(n-1,r) + C(n-1,r-1);
C := u;
end;
But .. is this really efficient? Or is there a better way to do it?
Han de Bruijn
Michael Jørgensen
"Han de Bruijn" <Han.de...@DTO.TUDelft.NL> wrote in message
news:44548$46dfede3$82a1e228$28...@news1.tudelft.nl...
The above is not very efficient. Each call to C() generates two further
calls to C(). Most of the values for C() are computed more than once. This
is highly inefficient.
The best way I know of is a simple loop. Consider C(10,4). The following
manipulations should give you the general idea:
C(10,4) = 10*9*8*7/(1*2*3*4) = 10/1 * 9/2 * 8/3 * 7/4.
The idea is that when you come to the next division, you are guaranteed that
the division is exact. The following example (in C) demonstrates the idea.
int C(int n, int r)
{
int u, i;
u = n;
for (i=2; i<=r; ++i)
{
u *= (n+1-i);
u /= i
}
}
Stephen Montgomery-Smith
Add this line - it might save you some time.
if (2*r>n) r = n-r;
matt271...@yahoo.co.uk
No, I don't think it is efficient. E.g. to calculate even just C(25,
12) generates over ten million calls to C. This number explodes with n
and k. Just doing it the "dumb" way is much quicker.
Han de Bruijn
function cc(n,r : integer) : integer;
{
Number of Combinations
}
var
u,i : integer;
begin
u := n;
for i := 2 to r do
begin
u := u * (n+1-i);
u := u div i;
end;
cc := u;
end;
It works! Thank you!
(Apart from a special case: IMHO your function gives the wrong result
for r = 0, only for r = 0)
Han de Bruijn
matt271...@yahoo.co.uk
> "Han de Bruijn" <Han.deBru...@DTO.TUDelft.NL> wrote in messagenews:44548$46dfede3$82a1e228$28...@news1.tudelft.nl...
If k > n/2 then by this method it's quicker to compute C(n, n - k)
(which is the same as C(n, k)).
Han de Bruijn
> On Sep 6, 1:09 pm, Han de Bruijn <Han.deBru...@DTO.TUDelft.NL> wrote:
>
>>Definition: C(n,k) = n!/((n-k)!k!)
>>
>>The best algorithm I could find for calculating C(n,k) is based upon
>>Pascal's triangle:
>>
>>C(n,k) = C(n-1,k) + C(n-1,k-1) ending with C(n,n) = 1 or C(n,0) = 1 .
>>
>>Giving the following recursive procedure (in Pascal):
>>
>>function C(n,r : integer) : integer;
>>{
>> Number of Combinations}
>>
>>var
>> u : integer;
>>begin
>> u := 1;
>> if (n > r) and (r > 0) then
>> u := C(n-1,r) + C(n-1,r-1);
>> C := u;
>>end;
>>
>>But .. is this really efficient? Or is there a better way to do it?
>
> No, I don't think it is efficient. E.g. to calculate even just C(25,
> 12) generates over ten million calls to C. This number explodes with n
> and k. Just doing it the "dumb" way is much quicker.
Yeah. I observed something like this (by experimenting). That's why ..
Han de Bruijn
Han de Bruijn
Yes. That's the same as Stephen Montgomery-Smith wrote:
> Add this line - it might save you some time.
>
> if (2*r>n) r = n-r;
Final version (?) :
function cc(n,r : integer) : integer;
{
Number of Combinations
}
var
u,i : integer;
begin
u := n;
if (2*r > n) then r := n-r;
for i := 2 to r do
begin
u := u * (n+1-i);
u := u div i;
end;
if r = 0 then u := 1;
cc := u;
end;
Han de Bruijn
matt271...@yahoo.co.uk
Looks good to me. Except, if I were being ultra-picky, for neatness I
would probably write the guts of it like this:
If (2 * r > n) Then r = n - r
If r = 0 Then
cc = 1
Else
cc = n
For i = 2 To r
cc = cc * (n + 1 - i)
cc = cc \ i
Next
End If
... if you'll excuse the change of syntax (I don't like posting code
examples that I haven't actually tested. It's impossibly error-prone!)
Jussi Piitulainen
> Final version (?) :
>
> function cc(n,r : integer) : integer;
> {
> Number of Combinations
> }
> var
> u,i : integer;
> begin
> u := n;
> if (2*r > n) then r := n-r;
> for i := 2 to r do
> begin
> u := u * (n+1-i);
> u := u div i;
> end;
> if r = 0 then u := 1;
> cc := u;
> end;
Starting with u := 1 and stepping from i := 1 to min(r, n - r)
should do it without treating r = 0 specially:
(define (cc n r)
(do ((u 1 (quotient (* u (- n k -1)) k))
(k 1 (+ k 1)))
((> k (min r (- n r))) u)))
(map cc '(4 4 4 4 4) '(0 1 2 3 4))
==> (1 4 6 4 1)
function cc(n,r : integer) : integer;
{
Number of Combinations
}
var
u,i : integer;
begin
if (2*r > n) then r := n-r;
u := 1;
for i := 1 to r do
begin
u := u * (n+1-i);
u := u div i;
end;
cc := u;
end;
Duncan Muirhead
> Definition: C(n,k) = n!/((n-k)!k!)
>
> The best algorithm I could find for calculating C(n,k) is based upon
> Pascal's triangle:
>
> C(n,k) = C(n-1,k) + C(n-1,k-1) ending with C(n,n) = 1 or C(n,0) = 1 .
>
> Giving the following recursive procedure (in Pascal):
>
<snip>
If you are using 64 bit integers, C(n,n/2) for n>67 won't fit.
(For 32 bits its n>34).
So you could have a table of around 2200 (64bit) integers that held
all the values. You could fill this table in an initialisation stage, or
even write a program to generate the table as a source file.
David W. Cantrell
> Han de Bruijn writes:
>
> > Final version (?) :
> >
> > function cc(n,r : integer) : integer;
> > {
> > Number of Combinations
> > }
> > var
> > u,i : integer;
> > begin
> > u := n;
> > if (2*r > n) then r := n-r;
> > for i := 2 to r do
> > begin
> > u := u * (n+1-i);
> > u := u div i;
> > end;
> > if r = 0 then u := 1;
> > cc := u;
> > end;
>
> Starting with u := 1 and stepping from i := 1 to min(r, n - r)
> should do it without treating r = 0 specially:
Thanks! IMO, treating r = 0 specially was quite inelegant.
> (define (cc n r)
> (do ((u 1 (quotient (* u (- n k -1)) k))
> (k 1 (+ k 1)))
> ((> k (min r (- n r))) u)))
>
> (map cc '(4 4 4 4 4) '(0 1 2 3 4))
> ==> (1 4 6 4 1)
>
> function cc(n,r : integer) : integer;
> {
> Number of Combinations
> }
> var
> u,i : integer;
> begin
> if (2*r > n) then r := n-r;
> u := 1;
> for i := 1 to r do
> begin
> u := u * (n+1-i);
> u := u div i;
> end;
> cc := u;
> end;
I don't know C (or whatever language was used above), but I would presume
that "integer" is not restricted to being nonnegative. If that is the case,
then the code is not robust because it fails to take into account the
possibility that n or r is negative. For example, we want cc(3,-1) to be 0.
David W. Cantrell
Han de Bruijn
Of course, one's favorite language essentially doesn't matter here.
Han de Bruijn
Jussi Piitulainen
> I don't know C (or whatever language was used above), but I
> would presume that "integer" is not restricted to being
> nonnegative. If that is the case, then the code is not
> robust because it fails to take into account the
> possibility that n or r is negative. For example, we want
> cc(3,-1) to be 0.
Correct. None of the presented procedures work right if the
arguments are not of the how-many-ways-to-choose-r-from-n
type. Also if r > n.
My Scheme procedure would do the wrong thing if the arguments
weren't integers. I suppose the others have severe overflow
problems; someone already pointed this out.
matt271...@yahoo.co.uk
> Jussi Piitulainen <jpiit...@ling.helsinki.fi> wrote:
> > Han de Bruijn writes:
>
> > > Final version (?) :
>
> > > function cc(n,r : integer) : integer;
> > > {
> > > Number of Combinations
> > > }
> > > var
> > > u,i : integer;
> > > begin
> > > u := n;
> > > if (2*r > n) then r := n-r;
> > > for i := 2 to r do
> > > begin
> > > u := u * (n+1-i);
> > > u := u div i;
> > > end;
> > > if r = 0 then u := 1;
> > > cc := u;
> > > end;
>
> > Starting with u := 1 and stepping from i := 1 to min(r, n - r)
> > should do it without treating r = 0 specially:
>
> Thanks! IMO, treating r = 0 specially was quite inelegant.
Except that, since maximum efficiency is the goal, it's worth noting
that treating r = 0 specially with an if-then-else clause, and
otherwise initialising cc to n, results in slightly quicker execution,
particularly for small r. In my tests the speed-up was about 12% for
C(20, 3).
>
>
>
>
>
> > (define (cc n r)
> > (do ((u 1 (quotient (* u (- n k -1)) k))
> > (k 1 (+ k 1)))
> > ((> k (min r (- n r))) u)))
>
> > (map cc '(4 4 4 4 4) '(0 1 2 3 4))
> > ==> (1 4 6 4 1)
>
> > function cc(n,r : integer) : integer;
> > {
> > Number of Combinations
> > }
> > var
> > u,i : integer;
> > begin
> > if (2*r > n) then r := n-r;
> > u := 1;
> > for i := 1 to r do
> > begin
> > u := u * (n+1-i);
> > u := u div i;
> > end;
> > cc := u;
> > end;
>
> I don't know C (or whatever language was used above), but I would presume
> that "integer" is not restricted to being nonnegative. If that is the case,
> then the code is not robust because it fails to take into account the
> possibility that n or r is negative. For example, we want cc(3,-1) to be 0.
Depends on the application really. Sometimes it makes sense to treat
C(3, -1), or C(2, 5), or whatever, as zero, and other times it means
that something has gone wrong and it's really an error condition.
Virgil
How about this one:
function C(n,r : integer) : integer;
{
Number of Combinations
}
var
u,v : integer;
begin
if n < 2*r
then u := C(n,n-r)
else
begin
u := 1;
for v := 1 to r do u := u*(n-v-1)/v;
end;
c := u;
end;
Jens Kruse Andersen
See e.g. http://www.luschny.de/math/factorial/FastBinomialFunction.html
For each p_i, compute p_i^e_i using repeated squarings and the fastest
algorithm for the operand sizes, FFT multiplication for really large
numbers.
Multiply all the p_i^e_ i into increasingly large numbers with the fastest
multiplication for the operand sizes.
Simple but probably not optimal order: Always replace the two smallest
numbers with their product.
If the final result is in binary then ignore 2^e_1 and add e_1 0s at the
end.
If the final result is in decimal and there is 2^e_1 and 5^e_3 then replace
2^e_1 with 2^(e_1-e_3) and ignore 5^e_3. Add e_3 0s at the end.
--
Jens Kruse Andersen