GMP 4.4 remaining tasks

Torbjorn Granlund tg at gmplib.org
Tue Dec 8 08:36:47 CET 2009

```bodrato at mail.dm.unipi.it writes:

We should try the "add a virtual limb" strategy... But if you have the
"exact or one unit too small" result and want to correct it... I suspect
that you need an expensive nxn multiplication...

Perhaps in the worst case, but in almost all cases a O(n) convolution
like product should suffice.  (I'd like to see an example of the latter,
actually.  This relates to the development of digits in a fraction; our
dividend is quite special.)

Now, I am not suggesting that we should necessarily correct the result
at any O(n) price.  If we get the result right almost always, that's
sufficient.  I believe we can get it right almost always by just
including an "extra" D limb during the Newton iterations, or perhaps
including an extra virtual limb at the end of D, or making sure the
penultimate approximation has at least ceil((n+1)/2) limbs.

> Getting mpn_invert to compute a slightly better value should not be
> hard.  For example, including another limb from D (called A in MCA)
> should help.  Also, one may avoid doubling the number of limbs in each

I did some tests, the "try mpn_inverse" test passed without any correction
with a code using "3/4" instead of "1/2" for the last iteration. Problems:
1) this does not mean that the code was correct for _any_ input (some
simulations I did with gp-pari suggest it was not)
2) that code was much _slower_, and this is obvious: the last step
required a (n*3n/4) and a (3n/4*3n/4) product instead of (n*n/2) and
(n/2*n/2)

Not a very good strategy if you ask me...

> iteration.  The basecase will surely be rather large, judging from
> mpn_binvert it will be over 500 limbs, which means an extra limb will
> imply a very little relative cost.

I did not collect timings very carefully, but my first attempt suggest a
smaller threshold, somewhere around 20 limbs, compared to the mpn_invert
currently available in the repository.

That would be *extremely* impressive.  You should swing your wand over
mpn_binvert next.  :-)

--
Torbjörn
```

More information about the gmp-devel mailing list