[Numpy-discussion] mpfr integration into numpy?

[Charles R Harris]

On Mon, May 6, 2013 at 7:57 AM, Nathaniel Smith <njs at pobox.com> wrote:

> On Mon, May 6, 2013 at 8:52 AM, Funky Dev <funkydevnull at gmail.com> wrote:
> > Hi,
> >
> > I've got a project in which it turns out we need much higher precision
> > than even __float128 and playing around with a few alternatives mpfr
> > seems to be the highest performing possibility.  So I've starting
> > writing a cythonized class mpfr_array which provides array-like
> > functionality but with mpfr_t as a "primitive" type.
> >
> > While this is doable it feels very much like re-inventing the wheel as
> > I'm basically rewriting some part of numpy's functionality.  So I was
> > wondering, as I think this is a very standard use case, is there any
> > interest in adding mpfr_t as a "primitive" dtype to numpy?  I know one
> > can define new dtype's but I don't want it to by a python object since
> > then there will be a large python overhead for all operations (i.e. dot
> > products will happen in python, not C).  mpfr_t is a very natural dtype
> > to add as its fast, C-based and supports general precision.
>
> You actually can define new dtypes implemented in C code -- in fact,
> right now, that's the only way to do it. There are some examples of
> how to do this from third-party code here:
>   https://github.com/numpy/numpy-dtypes
> Since mpfr is GPLed, it couldn't be included in core numpy in any
> case, but I don't see any reason you couldn't implement it this way
> and distribute it as a third-party extension.
>

The rational test class recently merged into numpy is more complete and
fixed
up than the version in numpy-dtypes, although I can't help thinking it
could be streamlined even more.


>
> There would be some subtleties to work out regarding memory management
> -- I guess mpfr_t objects have a fixed width (which is a requirement
> for dtypes), but internally malloc some buffers, so you have to be
> careful with copying them etc.? Should be doable with some care.
>
> > On a related note, if this was done, would it automatically work with
> > functionality such as numpy.linalg.inv(), etc...?  In principle such
> > functions could have been written with macros to be more type-flexible
> > (i.e. an ADD(a,b,c) mapping to e.g. a=b+c for floats or to mpfr_add(a,
> > b, c) for mpfr_t) but I suspect this is not the case.
>
> No, it wouldn't work automatically. But in numpy 1.8 you can register
> new implementations of 'inv' and friends when your library is loaded,
> which np.linalg.inv() will then automatically use if it sees your new
> dtype.
>
>
Chuck
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://mail.python.org/pipermail/numpy-discussion/attachments/20130506/0ef397c3/attachment.html>