Feature or enhancement
Proposal:
Currently we compute such powers as (1+0j)/z**(-n) (unless absolute value of n is too big to use specialized algorithm for integer exponents). This approach, however, introduce huge accuracy loss due to underflows in division.
For example:
>>> from gmpy2 import *
>>> import math
>>> x, y = map(float.fromhex, ['0x1.47e9c711723f5p+81',
... '0x1.38afd1168e49fp+85'])
>>> z = complex(x, y)
>>> pr = pow(z, -12); pr
0j
>>> pr2 = pow((1/z), 12); pr2
(5.562684646267994e-309+5.56268464626799e-309j)
>>> gr = complex(pow(mpc(z), -12))
>>> abs((pr2-gr).real)/math.ulp(gr.real)
2.0
>>> abs((pr2-gr).imag)/math.ulp(gr.imag)
2.0
Using instead ((1+0j)/z)**(-n) reduced error in this example from ~1e15 ULP to 2ULP. Note that, generic power algorithm is not affected by this issue:
>>> import _testcapi
>>> pr3 = _testcapi._py_c_pow(z, -12)[0]; pr3
(5.56268464626801e-309+5.56268464626799e-309j)
>>> abs((pr3-gr).real)/math.ulp(gr.real)
1.0
>>> abs((pr3-gr).imag)/math.ulp(gr.imag)
2.0
Has this already been discussed elsewhere?
This is a minor feature, which does not need previous discussion elsewhere
Links to previous discussion of this feature:
No response
Linked PRs
Feature or enhancement
Proposal:
Currently we compute such powers as
(1+0j)/z**(-n)(unless absolute value ofnis too big to use specialized algorithm for integer exponents). This approach, however, introduce huge accuracy loss due to underflows in division.For example:
Using instead
((1+0j)/z)**(-n)reduced error in this example from ~1e15 ULP to 2ULP. Note that, generic power algorithm is not affected by this issue:Has this already been discussed elsewhere?
This is a minor feature, which does not need previous discussion elsewhere
Links to previous discussion of this feature:
No response
Linked PRs