I spent a couple hours putting this together: https://github.com/omajoshi/complex_math
It does lazy updating, staying in polar form until an addition, then staying in rectangular form until the next multiplication. I implemented Tim Peters' tests (from two emails back): https://github.com/omajoshi/complex_math/blob/main/test.py Here are my results vs his: https://github.com/omajoshi/complex_math/blob/main/test.txt The library may be useless, but as far as I can it is associative (at least on Tim's tests). Om ---- On Wed, 23 Feb 2022 22:36:52 -0600 om <om+pyt...@omajoshi.com> wrote ---- > I stumbled across this: > https://math.stackexchange.com/questions/4172817/complex-floating-point-types, > which discusses implementing complex floating point in polar form (r: > unsigned double, theta: fixed point float), rather than rectangular form (x: > signed double, y: signed double). > > Of course, you lose out on easy addition if you implement C in polar form > (the accepted answer mentions this). But it smooths some of the other issues > mentioned above, as you get complex infinities `(+inf,theta)` and "signed" > zeroes `(epsilon,theta)`. One could even let theta represent the angle > scaled by a factor of 1/pi, to increase the precision of theta, so `i = > (1,0.5)` and `-1 = (1,1)`. > > One could also keep track of a complex number in both rectangular and polar > form (x, y, r, theta). So `i = (0,1,1,0,5)` and `-1 = (-1,0,1,1)`. Then just > use use rectangular coordinates for +- operations and polar coordinates for > */. And then after each +- operation, recalculate the polar coordinates, and > vice versa. > > Though now that I think about it more, that might be far more expensive, and > there might also be some issues with recalculating x,y for a complex > infinity. Thoughts? > > Om > > > > > ---- On Wed, 23 Feb 2022 21:40:58 -0600 Christopher Barker > <python...@gmail.com> wrote ---- > > > Tim, > > > > Does 754 say anything about complex numbers? > > > > It strikes me that what is needed is a defined behavior for complex > numbers, rather than expecting them to magically work for all the edge cases > simply by applying the algebraic rules with two floats. > > > > Your example of a complex infinity that should be one thing is a good > one. Maybe the same for complex zero :-) > > > > Anyway, Python is probably not the place to define a standard like this, > but if there's a good one out there, then using it in Python could be nice. > > > > -CHB > > > > > > > > > > > > On Tue, Feb 22, 2022 at 6:41 PM Tim Peters <tim.pet...@gmail.com> wrote: > > I have to say that signed zeroes and signed infinities don't work well > > in 754 for complex numbers. I know Kahan _wanted_ signed zeroes to > > help sort out branch cuts for complex functions, but he appears to be > > alone in the world in finding them useful for that ;-) > > > > Complex multiplication for values with signed zero components isn't > > even associative. At least not under any implementation I've tried. I > > pointed that out on David Hough's "numeric-interest" mailing list > > while 754 was still being hammered out, but it was "already too late" > > to do anything about it. > > > > Here's a short driver: > > > > from itertools import product > > pz = 0.0 > > cs = [complex(*t) for t in product((pz, -pz), repeat=2)] > > for x, y, z in product(cs, repeat=3): > > t1 = (x*y)*z > > t2 = x*(y*z) > > if repr(t1) != repr(t2): > > print(x, y, z, t1, t2) > > > > On my Windows 3.10.1, associativity fails in 24 (of the 64) cases: > > > > 0j 0j (-0-0j) -0j 0j > > 0j -0j (-0+0j) (-0+0j) 0j > > 0j -0j (-0-0j) -0j (-0+0j) > > 0j (-0+0j) (-0+0j) -0j 0j > > 0j (-0-0j) -0j -0j (-0+0j) > > 0j (-0-0j) (-0-0j) (-0+0j) 0j > > -0j 0j (-0+0j) (-0+0j) 0j > > -0j -0j (-0-0j) (-0+0j) 0j > > -0j (-0+0j) (-0+0j) (-0+0j) -0j > > -0j (-0+0j) (-0-0j) -0j 0j > > -0j (-0-0j) 0j (-0+0j) -0j > > -0j (-0-0j) (-0+0j) -0j 0j > > (-0+0j) 0j -0j 0j (-0+0j) > > (-0+0j) -0j 0j 0j (-0+0j) > > (-0+0j) (-0+0j) 0j 0j -0j > > (-0+0j) (-0+0j) -0j -0j (-0+0j) > > (-0+0j) (-0-0j) -0j 0j -0j > > (-0+0j) (-0-0j) (-0-0j) -0j (-0+0j) > > (-0-0j) 0j 0j 0j -0j > > (-0-0j) -0j 0j (-0+0j) -0j > > (-0-0j) -0j -0j 0j (-0+0j) > > (-0-0j) (-0+0j) -0j 0j -0j > > (-0-0j) (-0-0j) 0j 0j (-0+0j) > > (-0-0j) (-0-0j) (-0+0j) (-0+0j) -0j > > > > That isn't just academic. My employer's FORTRAN compiler failed to > > pass the US govt's validation suite at the time, because of a bizarre > > test failure: the accidental signs of these zeroes can have very > > visible consequences (e.g., as arguments to atan2(), which the > > validation suite called on the real and imag components after raising > > a complex zero to an integer power - and _which_ of the four complex > > zeroes you got depended on the grouping of the multiplies). > > > > Don't try to sell me the idea that any of that is logical ;-) BTW, > > exactly which cases in the above fail can also depend on the rounding > > mode (because x + -x is +0 for any finite x, _except_ under > > to-minus-infinity rounding, where the result changes to -0). > > > > Signed infinities are even sillier here. You want a single projective > > infinity in the complex plane, not a pair of signed infinities on each > > axis. And early drafts of 754 had a control bit to determine which > > flavor of infinity you got. But it's one of the few bits of esoterica > > that got dropped near the end. Partly because signed zeroes seemingly > > demand signed infinities too, like lutefisk demands haggis ;-) > > > > There was much to applaud in 754, but to my eyes signed zeroes, and > > NaN != NaN, caused far more trouble than they helped - and the > > exception model is such a pain it's almost never been implemented as > > intended (Python's `decimal` module being an exception, and Apple's > > long-dead SANE environment). > > > > BTW, complex was added as a full-blown built-in type in 1.4b1, long > > before Python had a major user base (mid-1990s). I think Guido added > > it because he was trained as a mathematician, so felt obligated ;-) > > Seriously, various Python versions of it (and rationals) were long in > > use as test cases for hammering out rules and APIs for coercion (and > > other native language features). > > _______________________________________________ > > Python-ideas mailing list -- python-ideas@python.org > > To unsubscribe send an email to python-ideas-le...@python.org > > https://mail.python.org/mailman3/lists/python-ideas.python.org/ > > Message archived at > https://mail.python.org/archives/list/python-ideas@python.org/message/3IY3NTJRAMSQZYTCVQABTYPM2TGS3UVN/ > > Code of Conduct: http://python.org/psf/codeofconduct/ > > _______________________________________________ > > Python-ideas mailing list -- python-ideas@python.org > > To unsubscribe send an email to python-ideas-le...@python.org > > https://mail.python.org/mailman3/lists/python-ideas.python.org/ > > Message archived at > https://mail.python.org/archives/list/python-ideas@python.org/message/ZXPQ5CUCMKLXBGRM352BHW72QLPIL2DP/ > > > Code of Conduct: http://python.org/psf/codeofconduct/ > > > _______________________________________________ Python-ideas mailing list -- python-ideas@python.org To unsubscribe send an email to python-ideas-le...@python.org https://mail.python.org/mailman3/lists/python-ideas.python.org/ Message archived at https://mail.python.org/archives/list/python-ideas@python.org/message/VVGWSBS5SDOCL3EE4P3ADVPL5QEWWCTW/ Code of Conduct: http://python.org/psf/codeofconduct/
