mpn_divexact_1 comments
Torbjorn Granlund
tg at gmplib.org
Wed Oct 16 15:16:07 CEST 2013
I agree with Niels (don't understand + not appropriate for
low-level...).
We should replace mpn_divexact_1 with code that:
(1) Uses Jebelean's trick with a Euclidean division working left-to-
right and a simultaneous Hensel division working right-to-left.
This is faster in the single-limb divisor case since both division
algorithms (Euclidean and Hensel) are typically limited by multiply
instruction latency.
(2) Uses a 2-limb modular inverse but only uses Hensel norm.
Optimal code without divisor invariance would probably look like this:
if (n < THRESHOLD_1)
binv = 1/d mod B
mpn_bdiv_1_pi1 (...)
else if (n < THRESHOLD_2)
binv = 1/d mod B²
mpn_bdiv_1_pi2 (...)
else
binv = 1/d mod B
inv = (B²-1)/d
mpn_div_1_pi1_bdiv_1_pi1 (...)
(This code assumes the divisor d is odd and the limb base B has no
common factor...)
For CPUs with pipelined multiply insns, I'd expect THRESHOLD_1 to be
small, perhaps 3 to 4. THRESHOLD_2 will be much larger, perhaps around
15. (Without pipelined multiply, mpn_bdiv_1_pi1 will always be faster.)
--
Torbjörn
More information about the gmp-devel
mailing list