two questions on GMP & C float

[Sylvain Pion]

Marc Glisse a écrit :
>> In the case of a float, there is a big chance that lots of zeros will end
>> up in the denominator or the numerator, because one of them will have
>> a large power of 2 factor.  Just consider 3.1415 as mpq.
> 
> That is quite a different issue though. These 0s are useful.

I think it's the only issue here.  Maybe I was just reading between the
lines of the OP's mail, because I know the issue.  But, just to make it clear,
I was talking of 3.1415, as a floating-point rounded value of the 3.1415
exact decimal, of course.  Then, once that is converted as an mpq, the denominator
is a power of 2.  Of course these zeroes are useful, but there exists a more
compact representation of a power of 2 than its binary string...

>> I have the same application as Vojtěch, geometric predicates, which are
>> usually signs of polynomial expressions with machine floats as input.
>> The ideal way to compute these exactly is IMO to use MPFR instead of 
>> mpq_t.
>>
>> Of course, once you use a division, you will need a rational type, so 
>> MPFR
>> is not enough, but maintaining a quotient of 2 "exact" MPFR will at least
>> compress the zeroes, which an mpq does not.
> 
> mpfr/mpz, so you have only one exponent? Or does using mpfr/mpfr reduce 
> the length of the code?

mpfr/mpz might be enough, although it might require the addition of
more interoperable operations, I don't know what's more practical.
Anyway, wasting one redundant exponent is probably not a severe issue,
so mpfr/mpfr might be enough, simpler, and more interoperable.

> I was just thinking:
> 
> assume the parameters are all integers and your predicate is a 
> polynomial with integer coefficients. Since you used interval 
> arithmetics first, in the few cases where you are doing the exact 
> computation, you know the result is going to be close to 0 (and you have 
> a bound on that). You could thus do the computations modulo some 
> slightly larger number (a power of 2^64 for instance, an unsigned long 
> may even be sufficient in some cases). I don't know how hard that would 
> be in your case, as the parameters are not integers, and there are often 
> divisions since noone is very fond of the homogeneous kernel :-(

Modular arithmetic here is cool, but it tends to be trickier to put
in place due to its constraints.


>> This is why I recently raised on the MPFR list the idea to add functions
>> dedicated to exact "polynomial" operations to MPFR (not that those are
>> hard to implement outside MPFR, of course).
> 
> What you are asking for is an exact ring and an exact field based on mpz 
> and mpq multiplied by a power of 2 (and stored in canonical form with 
> only odd integers, except possibly for 0).

An mpz_t multiplied by a power of 2 is nothing else than an mpfr_t
(apart from the inf/nan refinements) in terms of representable values.
This is my point.

> I don't think this fits with the goals of mpfr (where is the 
> *R*ounding?).

I would not characterize these operations as not fitting.
They are the exact math operations that you can do with the
values represented by mpfr_t, and its a mathematically
coherent set (from an algebraic structure point of view).
What else fits best...
To me, this is the set of most intrinsic operations that
you can do with this data type.  How could they not fit?
The thing which differs is that the current mpfr operations
are targeted with a fixed precision, hence the rounding.
I don't see the problem with adding some simple operations
which you can do exactly.

> If anything, it would fit better within gmp (than it would 
> fit with mpfr or even than mpf fits with gmp)...

The boundary between gmp and mpfr is very fuzzy anyway.  As far
as data types are concerned, mpfr_t can represent scaled integers,
just as mpf_t (and just as an std::pair<mpz_t, int> or whatever).

<historical context>
At some point in the past, there was a goal to have mpfr_t replace mpf_t.
This was given up, even though some things have changed which might
make it worth reconsidering, like the integration of MPFR as part of
the GNU project.  (BTW, what's the goal of mpf_t besides staying here
for backward compatibility?  Does it address areas which mpfr_t does
not, or not nicely at all?).
Unfortunately the outcome of this failed merge is that there are now
2 multiple precision floating-point types in the GNU project, with
the possible confusion and interoperability issues that can arise.
</historical context>

Do you suggest to add a third one in yet another separate library,
to deal with the operations in question (after all, inf/nan is not
needed there, and my operations are way different from mpf_t's since
they are exact, so "it does not fit"...)?  Or do you suggest to add
operations on mpfr_t in GMP?  Or do you suggest a "toy library"
re-using mpfr_t but not integrated in MPFR, so that it's going
to cost much more in terms of packaging&such as it is to maintain
it as part of MPFR?

In my view, those operations are either gathered in MPFR, or just
duplicated in whatever project where they are needed.  I think it
would strengthen MPFR to provide them, but it's just my user opinion
of course.

-- 
Sylvain Pion
INRIA Sophia-Antipolis
Geometrica Project-Team
CGAL, http://cgal.org/