numpy fails trigonometry with complex numbers, what to do?
I am trying to complete complex numbers in numpypy. Progress is good, I picked up from previous work on the numpypy-complex2 branch. Complex numbers come with extensive tests, it seems all the corner cases are covered. In porting the tests to numpypy, I came across a problem: numpy returns different results than cmath. Some of the differences are due to the fact that numpy does not raise a ValueError for dividing by 0 or other silly input values, but other differences are inexplicable (note the sign of the imaginary part):
numpy.arccos(complex(0.,-0.))
(1.5707963267948966-0j)
cmath.acos(complex(0.,-0.))
(1.5707963267948966+0j)
or this one:
cmath.acos(complex(float('inf'),2.3))
-infj
numpy.arccos(complex(float('inf'),2.3))
(0.78539816339744828-inf*j)
Should I ignore the inconsistencies, or fix the 700 out of 2300 test instance failures? What should pypy's numpypy do - be consistent with numpy or with cmath? cmath is easier and probably faster (no need to mangle results or input args), so I would prefer cmath to trying to understand the logic behind numpy. Matti
2280c2365844
On Wed, Sep 5, 2012 at 3:30 AM, Matti Picus matti.picus@gmail.com wrote:
I am trying to complete complex numbers in numpypy. Progress is good, I picked up from previous work on the numpypy-complex2 branch. Complex numbers come with extensive tests, it seems all the corner cases are covered. In porting the tests to numpypy, I came across a problem: numpy returns different results than cmath. Some of the differences are due to the fact that numpy does not raise a ValueError for dividing by 0 or other silly input values, but other differences are inexplicable (note the sign of the imaginary part):
numpy.arccos(complex(0.,-0.))
(1.5707963267948966-0j)
cmath.acos(complex(0.,-0.))
(1.5707963267948966+0j)
or this one:
cmath.acos(complex(float('inf'),2.3))
-infj
numpy.arccos(complex(float('inf'),2.3))
(0.78539816339744828-inf*j)
Should I ignore the inconsistencies, or fix the 700 out of 2300 test instance failures? What should pypy's numpypy do - be consistent with numpy or with cmath? cmath is easier and probably faster (no need to mangle results or input args), so I would prefer cmath to trying to understand the logic behind numpy. Matti
If you ask me, cmath is correct and numpy just didn't care. Maybe you should submit a bug report to them instead?
Am 05.09.12 03:30, schrieb Matti Picus:
I am trying to complete complex numbers in numpypy. Progress is good, I picked up from previous work on the numpypy-complex2 branch. Complex numbers come with extensive tests, it seems all the corner cases are covered. In porting the tests to numpypy, I came across a problem: numpy returns different results than cmath. Some of the differences are due to the fact that numpy does not raise a ValueError for dividing by 0 or other silly input values, but other differences are inexplicable (note the sign of the imaginary part):
numpy.arccos(complex(0.,-0.))
(1.5707963267948966-0j)
cmath.acos(complex(0.,-0.))
(1.5707963267948966+0j)
or this one:
cmath.acos(complex(float('inf'),2.3))
-infj
numpy.arccos(complex(float('inf'),2.3))
(0.78539816339744828-inf*j)
Should I ignore the inconsistencies, or fix the 700 out of 2300 test instance failures? What should pypy's numpypy do - be consistent with numpy or with cmath? cmath is easier and probably faster (no need to mangle results or input args), so I would prefer cmath to trying to understand the logic behind numpy. Matti
In NumPy you can change how numerical exception are handled: http://docs.scipy.org/doc/numpy/reference/routines.err.html http://docs.scipy.org/doc/numpy/user/misc.html#how-numpy-handles-numerical-e...
import numpy numpy.__version__
'1.6.2'
numpy.arccos(complex(float('inf'),2.3))
-c:1: RuntimeWarning: invalid value encountered in arccos (nan-inf*j) # Warning only once.
numpy.arccos(complex(float('inf'),2.3))
(nan-inf*j)
old_settings = numpy.seterr(all='raise') old_settings
Out[8]: {'divide': 'warn', 'invalid': 'warn', 'over': 'warn', 'under': 'ignore'}
numpy.arccos(complex(float('inf'),2.3))
--------------------------------------------------------------------------- FloatingPointError Traceback (most recent call last) <ipython-input-11-92051afcce38> in <module>() ----> 1 numpy.arccos(complex(float('inf'),2.3))
old_settings = numpy.seterr(all='ignore') numpy.arccos(complex(float('inf'),2.3))
(nan-inf*j)
HTH, Mike
Hi,
On Wed, Sep 5, 2012 at 3:30 AM, Matti Picus matti.picus@gmail.com wrote:
numpy.arccos(complex(0.,-0.))
(1.5707963267948966-0j)
cmath.acos(complex(0.,-0.))
(1.5707963267948966+0j)
cmath.acos(complex(float('inf'),2.3))
-infj
numpy.arccos(complex(float('inf'),2.3))
(0.78539816339744828-inf*j)
According to the C99 standard Annex G (draft, http://www.open-std.org/jtc1/sc22/wg14/www/docs/n1124.pdf), the cmath answer is the correct one in both cases. I don't know if that really means that numpy didn't care about the details. It sounds a bit strange given that it has tests for it; I fear it rather means that numpy implemented a different standard. But maybe that's me being too optimistic/pessimistic (depending on the point of view). I would indeed ask on numpy mailing lists or submit a bug entry and see their reaction.
A bientôt,
Armin.
After more than a year, I am trying to engage with numpy about their non-compliant C99 complex math routines. FWIW, the trigger for this was that we still see failures on our micronumpy -A tests with numpy 1.8. Matti
On 09/06/2012 10:12 AM, Armin Rigo wrote:
Hi,
On Wed, Sep 5, 2012 at 3:30 AM, Matti Picus matti.picus@gmail.com wrote:
numpy.arccos(complex(0.,-0.))
(1.5707963267948966-0j)
cmath.acos(complex(0.,-0.))
(1.5707963267948966+0j)
cmath.acos(complex(float('inf'),2.3))
-infj
numpy.arccos(complex(float('inf'),2.3))
(0.78539816339744828-inf*j)
According to the C99 standard Annex G (draft, http://www.open-std.org/jtc1/sc22/wg14/www/docs/n1124.pdf), the cmath answer is the correct one in both cases. I don't know if that really means that numpy didn't care about the details. It sounds a bit strange given that it has tests for it; I fear it rather means that numpy implemented a different standard. But maybe that's me being too optimistic/pessimistic (depending on the point of view). I would indeed ask on numpy mailing lists or submit a bug entry and see their reaction.
A bientôt,
Armin.
participants (4)
-
Armin Rigo -
Maciej Fijalkowski -
Matti Picus -
Mike Müller