GMP 4.4 remaining tasks

[Torbjorn Granlund]

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