classification
Title: pow() for complex numbers is rough around the edges
Type: enhancement Stage: test needed
Components: Interpreter Core Versions: Python 3.4
process
Status: open Resolution:
Dependencies: Superseder:
Assigned To: Nosy List: alex, jcea, mark.dickinson, mattip, rhettinger, serhiy.storchaka, terry.reedy
Priority: normal Keywords:

Created on 2012-09-20 22:00 by mattip, last changed 2013-10-13 18:26 by rhettinger.

Files
File name Uploaded Description Edit
rcomplex_testcases2.txt mattip, 2012-09-20 22:00 test cases for complex_power
Messages (9)
msg170856 - (view) Author: mattip (mattip) * Date: 2012-09-20 22:00
complex(1., 0.) ** complex(float('inf'), 0.) raises a ZeroDivisionError. In general, complex_power() needs to handle more corner cases. Barring a clear standard for pow() in C99, the documentation for pow 3 in glibc
http://www.kernel.org/doc/man-pages/online/pages/man3/pow.3.html
seems solid for a start, however it only describes behaviour for float/double values.

Where would be an appropriate place to add tests? I propose adding a test-case file similar to cmath_testcases.txt (attached) and a test runner similar to test_cmath.py
msg170866 - (view) Author: Mark Dickinson (mark.dickinson) * (Python committer) Date: 2012-09-21 07:57
Well, C99 covers pow for *real* numbers just fine;  it's complex numbers where no-one wants to pin down what the behaviour should be.  So I don't think we need the man page reference.

If we're writing tests for complex pow, we might also want to consider adding tests for multiplication and division;  those aren't entirely trivial either for special cases.

I do agree that (for the most part), complex pow applied to arguments with zero imaginary part should behave like regular float pow.  There are some cases where it's clear what the behaviour should be, and others that are murkier.

E.g., for a positive real z and arbitrary complex w, the special cases for z**w should behave in the same way as for exp for z > 0, and with some reflection of that behaviour for 0 < z < 1;  1**w should always be 1.

For nonzero finite values it's straightforward:  we just want to compute the best approximation to exp(w * log(z)), with the branch cut for the log along the negative real axis as usual.

But there are a *lot* of special cases to think about.  Consider that each real or imaginary part of the input is either:

(1) -infinity,
(2) -finite,
(3) -0.0
(4) +0.0
(5) +finite
(6) +infinity
(7) nan

and that we've got 2 complex inputs, or in effect 4 real inputs.  This divides our argument space into 7**4 = 2401 pieces.  With luck we can find rules that cover lots of those pieces at once, but it's still going to be a long job.

It doesn't help that it isn't particularly clear what the underlying mathematical model should be.  For floats, we can think about the two-point compactification of the real line (okay, with a doubled zero, which messes things up a little bit), which is a fairly sane space to work in.
msg170870 - (view) Author: Serhiy Storchaka (serhiy.storchaka) * (Python committer) Date: 2012-09-21 09:07
> Well, C99 covers pow for *real* numbers just fine;  it's complex numbers
> where no-one wants to pin down what the behaviour should be. 

C99 contains cpow. Perhaps we should use conditional compilation?
msg170891 - (view) Author: Mark Dickinson (mark.dickinson) * (Python committer) Date: 2012-09-21 14:08
> C99 contains cpow. Perhaps we should use conditional compilation?

I dread to think what horrors lurk in OS math library implementations of cpow;  I suspect we'd soon find out, if we had used cpow and have any tests at all for special cases.

OS math libraries are bad enough at *float* math, let alone complex;  I'd rather not depend on them unless we have to.  And given that at least on Windows we need our own complex pow implementation anyway, I'd prefer to use the same code on all platforms, so that we have at least some degree of consistency from platform to platform.
msg170951 - (view) Author: Terry J. Reedy (terry.reedy) * (Python committer) Date: 2012-09-22 02:14
Given that
>>> 1.0**float('inf'), 1.0**float('-inf')
(1.0, 1.0)

works,

>>> (1.0+0j)**(float('inf') + 0j)
Traceback ...
ZeroDivisionError: 0.0 to a negative or complex power

(and same for ('-inf') seems like a clear bug in raising an exception, let alone a clearly wrong exception. Clarification of murky cases, if it changes behavior, might be an enhancement.
msg171009 - (view) Author: Mark Dickinson (mark.dickinson) * (Python committer) Date: 2012-09-22 17:09
> (1.0+0j)**(float('inf') + 0j)

Oddly enough, this is nan+nanj on OS X.  I haven't investigated what the difference is due to---probably something to do with the errno results.
msg171048 - (view) Author: Mark Dickinson (mark.dickinson) * (Python committer) Date: 2012-09-23 15:59
Reclassifying this as an enhancement; I don't think it's appropriate to rewrite complex_pow for the bugfix releases.
msg199729 - (view) Author: Mark Dickinson (mark.dickinson) * (Python committer) Date: 2013-10-13 16:33
See also http://stackoverflow.com/q/18243270/270986 , which points out the following inconsistencies:

>>> 1e300 ** 2
OverflowError: (34, 'Result too large')
>>> 1e300j ** 2
OverflowError: complex exponentiation
>>> (1e300 + 1j) ** 2
OverflowError: complex exponentiation
>>> (1e300 + 1e300j) ** 2
(nan+nanj)
msg199752 - (view) Author: Raymond Hettinger (rhettinger) * (Python committer) Date: 2013-10-13 18:26
> OS math libraries are bad enough at *float* math,
> let alone complex;  I'd rather not depend on them unless we have to.

This makes good sense.   We should control how the special cases resolve and not be subject the whims of various C libraries.
History
Date User Action Args
2013-10-13 18:26:20rhettingersetnosy: + rhettinger
messages: + msg199752
2013-10-13 16:33:21mark.dickinsonsetmessages: + msg199729
2012-09-23 20:52:39jceasettype: behavior -> enhancement
stage: test needed
2012-09-23 20:51:50jceasetnosy: + jcea

type: enhancement -> behavior
stage: test needed -> (no value)
2012-09-23 15:59:10mark.dickinsonsettype: behavior -> enhancement
messages: + msg171048
versions: - Python 3.3
2012-09-22 17:09:00mark.dickinsonsetmessages: + msg171009
2012-09-22 02:14:06terry.reedysetversions: + Python 3.3, Python 3.4
nosy: + terry.reedy

messages: + msg170951

stage: test needed
2012-09-21 14:08:59mark.dickinsonsetmessages: + msg170891
2012-09-21 09:07:55serhiy.storchakasetmessages: + msg170870
2012-09-21 07:57:40mark.dickinsonsetmessages: + msg170866
2012-09-21 07:26:39serhiy.storchakasetnosy: + serhiy.storchaka
2012-09-20 22:02:18alexsetnosy: + alex
2012-09-20 22:00:27mattipcreate