Full numeric tower

[GoogleCodeExporter]

The idea of adding rational numbers to mpmath has been brought up before.
Now, I'm a fan of the Unix philosophy "do one thing, and do it well".
Mpmath is intended to do floating-point arithmetic well. But it turns out
that rational numbers help a lot even in a "pure" floating-point system;
just look at the current mess in specfun.py that is intended to handle
rational arguments for hypergeometric series.

It is arguably better to provide core rational and integer operations in
mpmath than only in SymPy. This way SymPy's classes don't get in the way of
the algorithms. And SymPy becomes simpler because it can import the
necessary functions and assume that they just work, and leave it to mpmath
to make sure that the fastest algorithm is used.

Implementing integers and rationals is trivial. The task is to make them
fit in nicely. Consider a full numeric tower:

mpz - mpq - mpf - mpi - mpc (- matrix)

mpi and mpc become generic and can hold any of the preceding types in the
list as parts. However, mpif/mpcf subclasses are created to take the roles
of the existing, optimized pure-float interval and complex types.

Operations propagate up the tower as expected:

Inexact mpz division -> mpq
Inexact mpz or mpq power -> mpf
Complex mpz, mpq or mpf power -> mpc

Beides the existing trappable ComplexResult exception, there will be a
trappable InexactResult for mpz(q) -> mpf promotions.

There will be some differences in behavior of functions. sqrt(4) will
return an mpz instead of an mpf (you'll have to ask for sqrt(4.0)). Ditto
for log(100,10). In principle, exp(0), sin(0), etc could return mpz as
well, although this is perhaps overkill (but arguably gamma(mpz) could
return an mpz).

Inf and nan could theoretically be spun off to an extended number type, but
I'm not sure if this would be useful. They fit fairly well within the
floating-point discipline.

With more types added, the architecture and especially the coercion code
for binary operations needs to be made much more generic (the existing code
for mpf and mpc alone is too complex). I've started looking at how this can
be done while preserving the present speed hacks for e.g. mpf*int, and I
think it can be done fairly smoothly.

Is an mpz class worth the trouble? The mpq class could arguably provide the
same functionality. The slowdown from using mpq to represent integers is
relatively small, at least compared to the slowdown of going from int to
mpz in the first place. The primary benefit of mpz might be for, say,
matrix, which can be optimized for the all-integer case more easily by
checking if the inputs are all int or mpz than if they are mpq with
denominator 1.

Original issue reported on code.google.com by fredrik....@gmail.com on 25 Aug 2008 at 2:52

[GoogleCodeExporter]

Yes, I am +1.

Btw, I think we are trying to solve the same problem as Sage has with it's 
coercion.
I haven't studied how they do it thoroughly yet, but it may be worth learning 
from
their mistakes.

Original comment by ondrej.c...@gmail.com on 26 Aug 2008 at 9:05

[GoogleCodeExporter]

http://www.dd.chalmers.se/~frejohl/code/libarith/

Why not add the rationals from there?

Original comment by Vinzent.Steinberg@gmail.com on 5 Nov 2008 at 5:26

• Changed state: Started

[GoogleCodeExporter]

No reason not to, except for sorting out the interface issues.

Original comment by fredrik....@gmail.com on 6 Nov 2008 at 5:38

[GoogleCodeExporter]

Which interface issues are to be sorted out?
sqrt(mpq) etc.?
sqrt(4)?

Original comment by Vinzent.Steinberg@gmail.com on 6 Nov 2008 at 9:44

[GoogleCodeExporter]

Yes, the meaning of sqrt(4) is one issue.

Whether mpq(3,2)*2 should return an integer (and similar, mostly speed related 
issues).

One alternative would be to keep the fast fraction type somewhere behind the 
scenes
and implement mpq as a wrapper. Then, in gmpy mode, the fast fraction type 
could even
be replaced by gmpy.mpq.

Original comment by fredrik....@gmail.com on 7 Nov 2008 at 10:49

[GoogleCodeExporter]

I have no problems with sqrt(4) being an integer.

I'm not sure about mpq simplifying to integer. Would it affect speed 
dramatically?

Please note that matrix([[1]]) does not simplify to a scalar.

+1 for mpq as a wrapper for the sake of using gmpy.

Original comment by Vinzent.Steinberg@gmail.com on 17 Feb 2009 at 12:41

[GoogleCodeExporter]

> I'm not sure about mpq simplifying to integer. Would it affect speed 
dramatically?

It does make a significant difference when values are likely to be integral. 
Small
Python integers are much faster than gmpy.mpq, even. In particular, there can 
be a
big advantage for a rational complex type where one component is likely to be 0 
or 1.

Original comment by fredrik....@gmail.com on 17 Feb 2009 at 6:44

[GoogleCodeExporter]

Is there actually any disadvantage when simplifying to integer? If not, let's 
just
implement it.

Original comment by Vinzent.Steinberg@gmail.com on 17 Feb 2009 at 9:00

[GoogleCodeExporter]

I might do it in connection with some refactoring of integer functions.

Original comment by fredrik....@gmail.com on 18 Feb 2009 at 6:51