Binomial improvements

[Torbjorn Granlund]

nisse at lysator.liu.se (Niels Möller) writes:

  At least a table för k should take little memory. For each k < k', it
  could tabulate not only k! but also a hensel inverse (of the odd part,
  and the shift count) of one or two limbs. Should be useful also for the
  case of large n and small k.
  
I implemented this, see https://gmplib.org/devel/.

The function becomes extremely simple:

void
mpz_bc_bin_uiui (mpz_ptr r, unsigned long int n, unsigned long int k)
{
  mp_limb_t i, j, c;

  c = n - k + 1;
  i = n - k + 2;
  for (j = 1; j < k; j++)
    {
      c = c * i;
      i++;
      c = c * zinv[j] >> ctz[j];
    }

  PTR(r)[0] = c;
  SIZ(r) = 1;
}

There are also tables, currently 32 limbs large, tabulating values also
for even numbers.  That might be a bit silly.  One might as well
tabulate just odd numbers, and ctz[j] to the index, like this
(untested):

  for (j = 1; j < k; j++)
    {
      c = c * i;
      i++;
      c = c * zinv[j >> ctz[j]] >> ctz[j];
    }

Is this the best way to handle really small numbers?  This code cannot
handle k > 32 (with present tables) and n greater then GMP_LIMB_BITS +
2.

One may write a variant that uses mpn_mul_1 instead, to allow any n that
fits in a limb.

Unlike my bdiv_bin_uiui code, this code puts all multiplies and
"divides" on the same dependency chain.  This might hurt performance on
some machines (with good multiply throughput and long latency).

Clearly, pre-generating the actual results for small n,k is faster....

-- 
Torbjörn