[patrick.pelissier@gmail.com: Fwd: Optimizing 1/x]

[Paul Zimmermann]

Patrick Pe'lissier forwarded me the following mail (I'm not on gmp-discuss).
Two remarks:

* what GMP needs is a fast inverse at the mpn level, which would compute,
  given a n-limb normalized input d, a n-limb quotient q (plus possibly one
  carry bit), and a n-limb remainder 0 <= r < d, such that

  2^(2*n*mp_bits_per_limb) = q * d + r.

  With such a function, computing 1/x in mpf would become easy.

* Newton's iteration isn't currently implemented in GMP, even for the
  integer division. The main reason is that you need to carefully analyze
  the roundoff errors to get a correct answer (as in the above proposal).
  I've implemented several variants which are quite fast; however the
  fastest known algorithm is currently under a patent from HP, and I failed
  so far to obtain a license for use within GMP...

Paul Zimmermann

- ---------- Forwarded message ----------
From: Ashod Nakashian <saghmos at xter.net>
Date: Sep 11, 2005 11:18 AM
Subject: Optimizing 1/x
To: gmp-discuss at swox.com


Hi,

Since GMP is highly optimized for the most common (and sometimes not so
common) calculations, I assumed that there would be an optimized
function for calculating the multiplicative inverse of a real (1/'mpf_t').

I don't seem to find such a function, so I assumed that mpf_ui_div might
have a specialized path for values of 1 (i.e. inversion), but there
seems to be no such path.

As I understand it, the code currently creates an mpf_t from the 'ui'
value. This value is normalized such that the division is done using the
mpf divisor as is and corrected by the exponent of the divisor.

This seems like a very simplistic approach that is sound (and probably
optimal) for arbitrary numerators. But the dire question is:

Why not use the Newton Iterations to computer 1/x?

And if it's not worth it (i.e. it won't be faster), why not simply
optimize the current algorithm using the fact that the numerator is a
simple 1.

It is hard to imagine that the GMP user base (in deed the developers as
well) didn't come across the need of a faster inversion calculation.

Am I missing something? Is it assumed since Newton's iterative method is
fairly easy to implement, one will do just that? (That still doesn't
sound like th GMP philosophy to me.)

So what do/did you guys do when you need(ed) to find 1/x quite rapidly?

- -Ash
_______________________________________________