stla / EllipticFunctions.jl

Jacobi theta functions and related functions.
https://stla.github.io/EllipticFunctions.jl/dev/
MIT License
16 stars 4 forks source link

kleinj does not work always #24

Open gufrei opened 1 year ago

gufrei commented 1 year ago

I made some test calculations with the kleinj function. It does not work for some numbers near the real axis. e.g.

julia> kleinj(1+0.01im) NaN + NaN*im

while kleinj(0.01im) works (both should be the same value). Another example

julia> kleinj(0.005025125628140614 + 0.001im ) NaN + NaN*im

while kleinj(0.005+ 0.001im ) produces a result. Plotting a color-coded phase plot of the function has visible errors near -1,0, and 1.

stla commented 1 year ago

The Klein-j function is defined on the open upper half-plane, ie the upper half plane without the real axis. It is then "normal" that its behaviour is not nice near the real axis.

Could you show the color-coded phase plot? I did some plots of this function and I didn't notice some problems. But welll, in my plots I keep a minimal distance from the real axis.

gufrei commented 1 year ago

I used the ComplexPortraits package to make the phase plots. I had to redefine the function by 0 (red color) where it is NaN to make it work. The area shown is upper left =-1.1 + 2.2im, lower right = 1.1 + 0.01im.
kleinj

stla commented 1 year ago

Nice pic. Indeed, it's weird that the problems are concentrated about -1 and 1, and not everywhere along the real axis.

stla commented 1 year ago

I have the same implementation in R and I found a way which works better sometimes. Instead of

  x <- lbd * (1 - lbd)
  256 * (1-x)^3 / x^2

you can do

  x <- lbd * (1 - lbd)
  256.0 * (1.0/x - 1)**2 * (1.0 - x)

With the second way, kleinj(1 + 0.01im) = 2.337834e+120-4.23137e+119im, whereas it is Inf - Inf im with the first way.

About the first problem, kleinj(0.005025125628140614 + 0.001im ), you can apply the formula kleinj(z) = kleinj(-1/z) and it works.

It would be interesting to compare with the Python library mpmath.

stla commented 1 year ago

Just checked with mpmath. Results are different. One can also get a different result while we use an equivalent expression of Klein-j. So we cannot trust these results.