Divide and conquer binomial

[bodrato at mail.dm.unipi.it]

Ciao,

Il Gio, 22 Dicembre 2011 6:42 pm, Torbjorn Granlund ha scritto:
> bodrato at mail.dm.unipi.it writes:

>   But the basecase now is quite a stupid one...
>
> What is the branching factor?  I.e., is the basecase invoked a number of
> times linear in the data size, or superlinear (like in Toom multiply)?

The basecase is invoked a power of two times, and the power is linear in k.
For binomial(n,k) the 2^h basecases are: compute (the odd part of) the
product of (k>>h) (maybe + 1) consecutive numbers divided by
factorial(k>>h).
Because of the "maybe +1", it is not always binomial(n-a*k>>h,k>>h), but
similar. In my current code, (the odd part of) factorial(k>>h) is computed
once, and used for all 2^h basecases.

The code is available from my web site at
http://bodrato.it/software/fac_ui.22-dic-2011.c

The file contains all the other combinatorial functions I'm working on,
the function we are speaking about is called BIN_UI; it pre-computes the
number of recursion steps and the divisor for each recursion level (stored
an an array), then calls the recursive mpz_dc_oddbin_ui function.
It can be inserted in the library copying this file in mpz/ and replacing
the file bin_ui.c with the following lines:

#define OPERATION_bin_ui
#include "nameyougavetothefile.ext"

But I prefer to compile it standalone with:
gcc -DMAIN -O2 -Impz -Itune fac_ui.22-dic-2011.c .libs/libgmp.a
tune/.libs/libspeed.a -o bin.test

Then, to test timings for binomial(12345678901234567890123, 654321>>a)
with a={20,19,...,0}, you can invoke: ./bin.test 654321 -U 20

If you want to look at the code and suggest, you are welcome.

Regards,
Marco

-- 
http://bodrato.it/software/combinatorics.html#bin_ui