bdiv vs redc

[Niels Möller]

Let us review the differences between bdiv and redc, and focus on the
balanced case (since that's what redc uses).

We have an odd modulo M of size n limbs, and a numerator U of size 2n
limbs.

bdiv computes

  Q = U M^-1 (mod B^n)
  R = B^{-n} (U - Q M)

where Q is n limbs, and R is n limbs plus a bit, it's in the range

  -M < R < B^n

The redc functions in gmp used to generate an R of n limbs and a Q of n
limbs and a bit, so we had quite a bit of impedance mismatch. But with
the current code, it seems that the redc functions compute

  Q = - U M^{-1} (mod B^n)
  R = B^{-n} (U + Q M)

where Q is n limbs (but not returned), and R is n limbs plus a bit,

  0 <= R < B^n + M

The redc functions return a carry, while the bdiv functions return a
borrow.

I think the reason for the redc changes was to let the caller decide if
the final conditional subtraction should be done unsing a constant time
method (powm_sec) or faster but with data-dependent timing (regular
powm).

But this means that if we were to implement a bdiv_r consistent with the
other bdiv functions, one could replace

  cy = mpn_redc_* (rp, up, mp, n, ...)
  if (cy != 0)
    mpn_sub_n (rp, rp, mp, n)

by

  cy = mpn_bdiv_r (rp, up, 2*n, mp, n, ...));
  if (cy != 0)
    mpn_add_n (rp, rp, mp, n);

That would be nice cleanup, I think. Now there are a couple of other
differences between the bdiv and redc interfaces:

1. The bdiv interface places the remainder at the high end of the U
   area, while redc arranges to put it at the low end. Can or should we
   change that? I guess it's possible to use a different convention for
   a new bdiv_r than for bdiv_qr. I guess the currrent bdiv_qr
   convention is beneficial for the dc routines, but I haven't looked at
   that code recently.

2. Obviously, the bdiv functions return the quotient. It seems
   mpn_sbpi1_bdiv_qr uses a loop very similar to redc, with addmul_1
   rather than submul_1. And then there's some logic to negate the
   quotient at the end; all that can be omitted for bdiv_r.

3. There's a redc_2 function, but no pi2_bdiv functions.

4. There's a dc_bdiv function, but no redc_dc. After some analysis, I
   think there's a place for redc_dc, even if we assume that the inverse
   computation for redc_n is free.

5. There's a mu_bdiv function, using a smaller inverse than redc_n,
   which is the closest corresponding function on the redc side. I guess
   it makes sense to always use a full inverse for redc, given that it's
   used mainly when the modulo is invariant?

Regards,
/Niels

-- 
Niels Möller. PGP-encrypted email is preferred. Keyid C0B98E26.
Internet email is subject to wholesale government surveillance.