divappr or approximate fraction (a diversion from sqrt...)

Niels Möller nisse at lysator.liu.se
Sun Oct 18 06:31:26 UTC 2015


tg at gmplib.org (Torbjörn Granlund) writes:

> nisse at lysator.liu.se (Niels Möller) writes:
>
>   When I looked into sqrt some months ago, I found that it would be nice
>   with a divappr function which didn't require the low mostly
>   uninteresting dividend limbs as input, a bit similar to bdiv_q.
>   
>   Maybe one shouldn't think about that as an approximate division, but
>   rather as an approximate fraction. I've implemented one variant (see
>   below), which takes as input two n limb numbers U and D, with D
>   normalized. 
>   
> This is akin to good old mpn_tdiv_qr with its extra-limb argument which
> we never actually allowd to be != 0.  I like the idea to provide the
> fraction thing separately.

And then the fraction function could perhaps require that U < D, so
there's no non-fraction bit.

> Depending on submul_1c is a bit of a practical problem, since it not
> available in general.  We also lag in implementing addmul_1c (and lack
> completely e.g. addmul_2c).

I think it would make sense to make addmul_1c and submul_1c supported
and public (we discussed a while ago how to do the C versions, if we can
get away with only the _1c functions, and let the public _1 calls be via
a macro or function wrapper calling _1c with zero. In effect, _1c would
be the mandatory function, and a cheaper _1 entry point would be
optional). But I understand it's going to be some non-trivial work. I
had a quick look at implementing _1c for coreihwl, but it wasn't obvious
to me which register held the recurrency limb. 

Not sure about _2c functions.

> Long ago I looked into using addmul_1 for plain division (not bdiv);
> that would speed up things and also simplify low-level code.  We might
> want to look into that again.

Say U is n limbs, and U' is the complement, B^n-1 - U. Then R = U -
QD < D is the remainder, we have

  U' + QD = B^n - 1 - R 

Is that the trick you've been considering? Complement U up front, then
add multiples of D in school-book-style to "eliminate" high limbs by making them equal
to B-1. This produces an R', which is the complement of the remainder.

Regards,
/Niels

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


More information about the gmp-devel mailing list