GMP 4.4 remaining tasks

[Niels Möller]

bodrato at mail.dm.unipi.it writes:

> I did try to write an invert.c using Newton's iterations... but I was not
> able to get a correct and fast code.

> Comments or suggestions?

I haven't looked closely neither at your code or at left-to-right
inversion in general. But maybe you can get some ideas from Alg 3.5
ApproximateReciprocal (page 121 in version 0.4) in Brent's and
Zimmermann's MCA book?

The algorithm as written is recursive, but it should be fairly
straight-forward to write it iteratively. The iteration is of the form
(for simplicity, I omit the powers of two needed to get the binary point
right).

1.  T = X_k * A_k
2.  E = 1 - T
3.  Small number of O(n) adjustments in the case E < 0
4.  X_{k+1} = X_k + X_k * E

Here A_k is the most significant part of the input, of size appropriate
for the iteration. A is assumed to be "strictly" normalized with high
limb >= B/2, to relax that I guess one would need one more limb in A.

The first multiplication (line 1) is an obvious target for wraparound
logic. The second multiplication (line 4) would be some kind of mulhi,
if I understand this correctly, where some of the lower limbs of the
result are discarded.

This algorithm looks quite clean and simple to me, and for large
inverses I don't think it's a serious drawback to have to do a cheap
adjustment step for each Newton iteration.

Regards,
/Niels

-- 
Niels Möller. PGP-encrypted email is preferred. Keyid C0B98E26.
Internet email is subject to wholesale government surveillance.