gcdext_1 questions
Torbjorn Granlund
tg at gmplib.org
Sat Dec 5 23:06:51 CET 2009
Since you're talking gcd for small operands, let me contribute a
variation of Euclid:
mp_limb_t
gcd (mp_limb_t a, mp_limb_t b)
{
mp_limb_t q, d, r;
ASSERT (a + b/2 >= a);
do
{
q = (a + b/2) / b;
d = q * b;
if (d > a)
r = d - a;
else
r = a - d;
a = b;
b = r;
}
while (r != 0);
return a;
}
The quotient is rounded to nearest (approximately) by adding +b/2. That
means the remainder a - qb would be negative half the time; the
remainder is infact centered around 0. Since we don't want to deal with
negative numbers, we compute immediately the absolute value of the
remainder.
This is a simple trick of making Euclid's algorithm converge a good bit
faster.
(The if statement should of course be organised such as there will be no
branch, since such a branch's outcome pattern will be completely
unpredictable.)
--
Torbjörn
More information about the gmp-devel
mailing list