PEP239 (Rational Numbers) Reference Implementation and new issues

[Alex Martelli]

Mark Day wrote:

> In article <ank8im$33l$0 at 216.39.172.122>, Bengt Richter <bokr at oz.net>
> wrote:
> 
>> >Use of floating point values in such operations would either produce
>> >errors:
>> >
>> >  1_3 + 2.5 -> TypeError
>> Why? All floating point numbers can be viewed as exact rationals
>> of the form n_<2**e>  with e>=0
> 
> Automatic conversion to rational pretends those float values are exact,
> when they probably weren't meant that way.  I would prefer requiring an
> explicit conversion rather than an implicit one.

In gmpy, at the courteous but adamant insistence of Pearu Peterson,
I introduced a Stern-Brocot algorithm to convert floats to rationals
giving the smallest denominator compatible with the hypothesis that
the float is precise to N bits (N can be any, since gmpy also exposes
via the mpf type the multi-precision floats of the underlying GMP).

The Farey-algorithm, continued-fractions, and other means which have
recently been discussed in python-Dev are presumably more or less
equivalent (I'll of course be quite happy to release the Stern-Brocot
snippet separately from gmpy if anyone needs it with a softer license:
gmpy is LGPL, because I think it may be considered a "derived work"
of GMP, which is itself LGPL).

By default, the number of bits assumed to be accurate for floats
are "all of them" -- usually the DBL_MANT_BITS (or DBL_MANT_DIGS
if FLT_RADIX is 2) from the system's float.h, normally 53 today --
i.e., an imprecision of half an ULP on either side is always
assumed to be possible.

The results seem to be OK for Pearu's symbolic-mathematics purposes,
for casual use, and for every use I've made of the feature myself
(mostly in peculiar combinatorial-arithmetic exercises connected to
card-play probabilities -- but I don't often end up with floats
there except when I'm being lazy and taking shortcuts, admittedly).

[alex at lancelot alex]$ python
Python 2.3a0 (#1, Oct  4 2002, 13:05:05)
[GCC 2.96 20000731 (Mandrake Linux 8.2 2.96-0.76mdk)] on linux2
Type "help", "copyright", "credits" or "license" for more information.
>>> from gmpy import mpq, mpf, f2q
>>> mpq('1/3') + 2.5
mpq(17,6)
>>> mpq(1/3.)
mpq(1,3)
>>> mpq('1/3')
mpq(1,3)
>>> mpq(1,3) + 1.0/3.0
mpq(2,3)
>>>

> Imagine the uproar when repr(R(1.0/3.0)) looks nothing like
> repr(R(1,3)).

That might be bad, I guess, but as you see that's not what
happens with gmpy -- Stern-Brocot make it "recover" the value
of 1.0/3.0 as being best approximated by the rational 1/3.

Of course, if I wrongfully claimed to have more bits of
precision than I actually do...:

>>> f2q(1.0/3.0, -53)
mpq(1,3)
>>> f2q(1.0/3.0, -54)
mpq(1,3)
>>> f2q(1.0/3.0, -55)
mpq(6004799503160661L,18014398509481984L)
>>>

...then at some point the rational "explodes" in terms of
size of numerator and denominator.  f2q accepts the number
of bits as a second negative argument, or the maximum
tolerable error/approximation when the second argument is 
positive, so you don't _have_ to think in binary temrs if 
that's not convenient in your application:

>>> f2q(1.0/3.0, 0.0001)
mpq(1,3)
>>>


Alex