mpmath icon indicating copy to clipboard operation
mpmath copied to clipboard

Full numeric tower

Open fredrik-johansson opened this issue 17 years ago • 9 comments

bc.. 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.

p. Original issue for "#93":https://github.com/fredrik-johansson/mpmath/issues/93: "http://code.google.com/p/mpmath/issues/detail?id=53":http://code.google.com/p/mpmath/issues/detail?id=53

p. Original author: "https://code.google.com/u/111502149103757882156/":https://code.google.com/u/111502149103757882156/

p. Original owner: "https://code.google.com/u/111502149103757882156/":https://code.google.com/u/111502149103757882156/

fredrik-johansson avatar Aug 26 '08 01:08 fredrik-johansson

bc.. 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.

p. Original comment: "http://code.google.com/p/mpmath/issues/detail?id=53#c1":http://code.google.com/p/mpmath/issues/detail?id=53#c1

p. Original author: "https://code.google.com/u/104039945248245758823/":https://code.google.com/u/104039945248245758823/

certik avatar Aug 26 '08 09:08 certik

bc.. http://www.dd.chalmers.se/~frejohl/code/libarith/ Why not add the rationals from there?

Status: Started

p. Original comment: "http://code.google.com/p/mpmath/issues/detail?id=53#c2":http://code.google.com/p/mpmath/issues/detail?id=53#c2

p. Original author: "https://code.google.com/u/[email protected]/":https://code.google.com/u/[email protected]/

vks avatar Nov 05 '08 17:11 vks

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

p. Original comment: "http://code.google.com/p/mpmath/issues/detail?id=53#c3":http://code.google.com/p/mpmath/issues/detail?id=53#c3

p. Original author: "https://code.google.com/u/111502149103757882156/":https://code.google.com/u/111502149103757882156/

fredrik-johansson avatar Nov 06 '08 17:11 fredrik-johansson

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

p. Original comment: "http://code.google.com/p/mpmath/issues/detail?id=53#c4":http://code.google.com/p/mpmath/issues/detail?id=53#c4

p. Original author: "https://code.google.com/u/[email protected]/":https://code.google.com/u/[email protected]/

vks avatar Nov 06 '08 21:11 vks

bc.. 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.

p. Original comment: "http://code.google.com/p/mpmath/issues/detail?id=53#c5":http://code.google.com/p/mpmath/issues/detail?id=53#c5

p. Original author: "https://code.google.com/u/111502149103757882156/":https://code.google.com/u/111502149103757882156/

fredrik-johansson avatar Nov 07 '08 22:11 fredrik-johansson

bc.. 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.

p. Original comment: "http://code.google.com/p/mpmath/issues/detail?id=53#c6":http://code.google.com/p/mpmath/issues/detail?id=53#c6

p. Original author: "https://code.google.com/u/[email protected]/":https://code.google.com/u/[email protected]/

vks avatar Feb 17 '09 12:02 vks

bc.. > 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.

p. Original comment: "http://code.google.com/p/mpmath/issues/detail?id=53#c7":http://code.google.com/p/mpmath/issues/detail?id=53#c7

p. Original author: "https://code.google.com/u/111502149103757882156/":https://code.google.com/u/111502149103757882156/

fredrik-johansson avatar Feb 17 '09 18:02 fredrik-johansson

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

p. Original comment: "http://code.google.com/p/mpmath/issues/detail?id=53#c8":http://code.google.com/p/mpmath/issues/detail?id=53#c8

p. Original author: "https://code.google.com/u/[email protected]/":https://code.google.com/u/[email protected]/

vks avatar Feb 17 '09 21:02 vks

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

p. Original comment: "http://code.google.com/p/mpmath/issues/detail?id=53#c9":http://code.google.com/p/mpmath/issues/detail?id=53#c9

p. Original author: "https://code.google.com/u/111502149103757882156/":https://code.google.com/u/111502149103757882156/

fredrik-johansson avatar Feb 18 '09 18:02 fredrik-johansson