Imaginary type and IEC 60559-compatible complex arithmetic

[skirpichev]

"Generally, mixed-mode arithmetic combining real and complex variables should be performed directly, not by first coercing the real to complex, lest the sign of zero be rendered uninformative; the same goes for combinations of pure imaginary quantities with complex variables." (c) Kahan, W: Branch cuts for complex elementary functions.

That's why C standards since C99 introduce imaginary types. This patch proposes similar extension to the Python language:

• Added a new subtype (imaginary) of the complex type. New type has few overloaded methods (conjugate() and __getnewargs__()).
• Complex and imaginary types implement IEC 60559-compatible complex arithmetic (as specified by C11 Annex G).
• Imaginary literals now produce instances of imaginary type.
• cmath.infj/nanj were changed to be of imaginary type.
• Modules ast, code, copy, marshal got support for imaginary type.
• Few tests adapted to use complex, instead of imaginary literals
  • Lib/test/test_fractions.py
  • Lib/test/test_socket.py
  • Lib/test/test_str.py

Lets consider some (actually interrelated) problems, shown for unpatched code, which could be solved on this way.

1. First, new code allows to use complex arithmetic for implementation of mathematical functions without special "corner cases". Take the inverse hyperbolic sine as an example:

>>> z = complex(-0.0, 2)
>>> cmath.asinh(z)
(-1.3169578969248166+1.5707963267948966j)
>>> # naive textbook formula doesn't work:
>>> cmath.log(z + cmath.sqrt(1 + z*z))
(1.3169578969248166+1.5707963267948966j)
>>> # "fixed" version does:
>>> cmath.log(z + cmath.sqrt(complex(1 + (z*z).real, (z*z).imag)))
(-1.3169578969248164+1.5707963267948966j)

2. Previously, we have only unsigned imaginary literals with the following semantics:

±a±bj = complex(±float(a), 0.0) ± complex(0.0, float(b))

While this behaviour was well documented, most users would expect instead here:

±a±bj = complex(±float(a), ±float(b))

i.e. that it follows to the rectangular notation for complex numbers.

Things are worse, because the CPython docs sometimes asserts that the rectangular form is used and that some simple invariants holds. For example, sphinx docs for the complex class says: "complex(real=0, imag=0) ... Return a complex number with the value real + imag*1j ...". But:

>>> complex(0.0, cmath.inf)
infj
>>> 0.0 + cmath.inf*1j
(nan+infj)

3. The eval(repr(x)) == x invariant was broken for the complex type. Below are simple examples with signed zero:

>>> complex(-0.0, 1.0)  # also note funny signed integer zero below
(-0+1j)
>>> -0+1j
1j
>> -(0.0-1j)  # "correct" input for above with Python numeric literals
(-0+1j)
>>> -(0-1j)  # also "correct"; integer 0 matters!
(-0+1j)
>>> complex(1.0, -0.0)
(1-0j)
>>> 1-0j
(1+0j)
>>> -(-1 + 0j)
(1-0j)

Similar was true for complex numbers with other special components:

>>> complex(0.0, -cmath.inf)
-infj
>>> -cmath.infj
(-0-infj)

See also ISO/IEC 9899:2011, Annex G and Rationale for C99 and Augmenting a Programming Language with Complex Arithmetic.

---

📚 Documentation preview 📚: https://cpython-previews--1.org.readthedocs.build/