Update: Schönhage 's algorithm for fast integer GCD

[Niels Möller]

Torbjorn Granlund <tege at swox.com> writes:

> You have your own basecase code, I presume?

Yes. There are actually two basecases. The basecase for hgcd is
different from the old code, because it doesn't compute the full gcd
(the basic idea is to compute r = u a + b v where all of r, u and v
are *half* the size of a). And because of the way the rest of the
algorithm works, it seems non-trivial to use a division-free algorithm
for this basecase.

Then there's also a basecase for the outer gcd loop, which is a plain
gcd. It currently uses lehmer/jebelean, but I guess it could be
changed to use a different algorithm.

About the same should apply to gcdext, but I haven't finished that
code yet.

>   Timings, relative to gmp-4.1.2 mpz_gcd, on sparc, are included below.
>   
> I think a venerable sparc v7 is a poor timing platform for GMP.

;-) I'll run the benchmark on my non-ancient x86 machine at work too.
But for the earlier incarnations of the code, the thresholds turned
out to be quite similar for the two machines.

>   Most of the time is spent in computations on two-limb numbers (that's
>   div2 and nhgcd2) and book-keeping of one-limb quotients
> 
> Is that the same div2 as of the existing mpn_gcdext?  If not,
> what does it compute?

div2 is almost the same, I have only hacked it to also return the
remainder. So perhaps it should be called divmod2.

Regards,
/Niels