patch: _mpmath_mpf_mul

[GoogleCodeExporter]

In the attached patch there is _mpmath_mpf_mul;
defining in mpmath rev843 in libmpf.py  
gmpy_mpf_mul = gmpy._mpmath_mpf_mul
mpf_mul is 2x faster at mp.dps = 100

In the patch there is also a slight speed-up for _mpmath_create

Original issue reported on code.google.com by mario.pe...@gmail.com on 19 Jan 2009 at 4:30

Attachments:

• gmpy_diff

[GoogleCodeExporter]

This seems good to me, but I'm uneasy about handling inf and nan in gmpy. I 
consider
the representation and semantics of inf and nan in mpmath subject to change. If
possible, I would suggest handling those via some kind of callback to mpmath.

Original comment by fredrik....@gmail.com on 19 Jan 2009 at 5:41

[GoogleCodeExporter]

I'll have more time to look at this in a few days. If we are going to support
multiply, we should support addition, subtraction, division, etc. Then we could 
have
a correctly rounded replacement for gmp's mpf. Changes like that should be 
reserved
for the next version of gmpy and we would need to agree on an internal 
representation
etc.

Original comment by casevh on 21 Jan 2009 at 6:45

• Added labels: Type-Enhancement
• Removed labels: Type-Defect

[GoogleCodeExporter]

In the attached patch _mpmath_mpf_mul does not depend on the explicit 
representation
of inf and nan, which are taken from libmpf.py adding there a call to
gmpy._mpmath_cache after defining the special numbers.

There is also _mpmath_mpf_add, with which mpf_add is 2x faster at mp.dps = 100

There is an inconsistency in the treatment of man in the two functions; 
in _mpmath_mpf_mul 
sman = (PympzObject *)PyTuple_GET_ITEM(s, 1)
as in _mpmath_normalize
while in _mpmath_mpf_add
sman = anynum2mpz(PyTuple_GET_ITEM(s, 1))
as in _mpmath_create
The reason is that in runtests there are some mpf_tuples with man not an mpz,
which appear as arguments to mpf_add, not to mpf_mul.
If in mpmath it were guaranteed that all mpf tuples have mpz argument, 
it would not be needed to call anynum2mpz in _mpmath_mpf_add .

Original comment by mario.pe...@gmail.com on 22 Jan 2009 at 10:56

Attachments:

• gmpy_diff1

[GoogleCodeExporter]

Let's discuss possibly even more radical approaches.

Note that mpmath doesn't assume that mpf values are tuples (except in a few 
places
that could easily be edited); only that they can be unpacked as such.

If sufficiently many core functions (probably not that many -- arithmetic,
comparison, float/int conversion) are implemented in gmpy, it is likely much 
faster
to implement a custom representation type in gmpy that constructs the component
objects only when required by external code. Such a class could hold and 
manipulate
the sign and bitcount data as machine integers.

In fact, it might not be too much effort to use a mixed representation for the
exponent, where it is represented by a machine integer by default (99% of the 
time)
and becomes an mpz pointer only when the exponent is huge. Since exponents are 
mostly
added and subtracted, the overflow checking should be easy.

Original comment by fredrik....@gmail.com on 26 Jan 2009 at 2:05

[GoogleCodeExporter]

In the attached patch to gmpy I added _mpmath_mpf_div .
The patch to mpmath has the modifications to use the first patch.
In timings_mpmath.py at dps = 100 div is 3x faster.

Original comment by mario.pe...@gmail.com on 2 Feb 2009 at 12:55

Attachments:

• gmpy_diff3
• mpmath_diff3

[GoogleCodeExporter]

In the attached patch to gmpy I added _mpmath_mpf_mul_int.
The patch to mpmath has the modifications to use the first patch.

It seems to me that there is a bug somewhere, so that anynum2mpz gives wrong
results in _mpmath_mpf_mul_int; to avoid this bug I introduced anynum2mpz1,
with which runtests is OK.

In this example
from mpmath.libmpf import *
import gmpy
from gmpy import mpz
s, n, prec, rnd = (0, mpz(118587876497), 0, 37), -2, 212, 'n'
r = mpf_mul_int(s, n, prec, rnd)
print r
r = mpf_mul_int(s, n, prec, rnd)
print r

using anynum2mpz
(1, mpz(118587876497), 1, 37)
(0, mpz(118587876497), 1, 37)

using anynum2mpz1 it is OK
(1, mpz(118587876497), 1, 37)
(1, mpz(118587876497), 1, 37)

I traced the bug to the fact that in two calls to anynum2mpz, first
zconst[0]=-2, then zconst[0]=2
but I didn't go farther. I find this bug on hp Pavilion Q8200 2.33GHz
and on an AMD Athlon XP 2600 .

Original comment by mario.pe...@gmail.com on 4 Feb 2009 at 11:29

Attachments:

• gmpy_diff5
• mpmath_diff5