Bug: `jtheta_ab` return incorrect results at certain points

[hammerfunctor]

$z=(-a+1/2)\tau + (-b+1/2)$ should be a zero of $\vartheta_{a,b}(z,\tau)$, but

julia> jtheta_ab(1/3,1/3,pi*( (-1/3+1/2)*(2+im) -1/3+1/2 ), qfromtau(2+im))
-0.5539025646339003 + 0.660115371349597im

[hammerfunctor]

It's getting even more weird for me

julia> veri2(τ,a,b) = EllipticFunctions._jtheta_ab.(promote(a,b, pi*((-a+1/2)*τ-b+1/2), τ)...)
veri2 (generic function with 1 method)

julia> veri3(τ,a,b) = EllipticFunctions.jtheta_ab(a,b, pi*((-a+1/2)*τ-b+1/2), qfromtau(τ))
veri3 (generic function with 1 method)

julia> veri2(2+im,1/3,1/3)
2.554007772509945e-17 + 5.76594535568097e-17im

julia> veri3(2+im,1/3,1/3)
-0.5539025646339003 + 0.660115371349597im

[hammerfunctor]

Alright, I figured it out. For some reason the periodicity in real part of $\tau$ doesn't show up in the result:

julia> taufromq(qfromtau(2+im))
-7.796343665038751e-17 + 1.0im

julia> veri4(τ,a,b) = EllipticFunctions._jtheta_ab.(a,b, pi*((-a+1/2)*τ-b+1/2), taufromq(qfromtau(τ)))
veri4 (generic function with 1 method)

julia> veri4(1.9+im,1/3,1/3)
-0.5065059738916385 + 0.6971456559948964im

There's no such problem in python's mpmath

julia> using PyCall

julia> mp = pyimport("mpmath")
PyObject <module 'mpmath' from '/home/huzf/.julia/conda/3/x86_64/lib/python3.10/site-packages/mpmath/__init__.py'>

julia> thetaab(z,τ,a,b) = exp(im*pi*a^2*τ + 2*pi*im*a*(z+b)) * mp.jtheta(3, pi*(z+a*τ+b),qfromtau(τ))
thetaab (generic function with 1 method)

julia> veri(τ,a,b) = thetaab((-a+1/2)*τ-b+1/2, τ, a, b)
veri (generic function with 1 method)

julia> veri(2+im,1/3,1/3)
PyObject mpc(real='6.2935875993041592e-16', imag='6.9820195867768032e-16')

[stla]

I don't think this relation is correct with our implementation. There should be a factor pi missing somewhere. I'll try tomorrow.

[simonp0420]

Note that τ is not a single-valued function of $q$. Given a value of $q$ as an input to jtheta_ab, there is no way to determine which τ value is desired as the output of taufromq. We use the standard convention, selecting the principal branch of the complex log function to choose a particular value of τ to return. So in @hammerfunctor's original test, if one selects τ using τ = taufromq(qfromtau(2 + im)), it works:

julia> τ = taufromq(qfromtau(2 + im))  # move τ to the principal branch
-7.796343665038751e-17 + 1.0im

julia> jtheta_ab(1/3,1/3,pi*( (-1/3+1/2)*τ -1/3+1/2 ), qfromtau(τ))
6.12184648516288e-17 - 1.5139607520801326e-17im

Since we use $q$ as the fundamental input for all the user-facing functions (and not τ), I believe this situation is fine. The user must select τ using the principal branch of the log function (or by just using taufromq) whenever he/she wants to employ a formula involving τ, such as the zero identity in the first post above.

[hammerfunctor]

I don't think theta function with characteristics is well-defined as a function of $q$. Take the periodicity 1 notation $$\vartheta{a, b} (z, \tau) = e^{i \pi a (a \tau + z + 2 b)} \vartheta (z + a \tau + b, \tau) = \sum{n \in \mathbb{Z}} e^{\pi i (n + a)^2 \tau + 2 \pi i (n + a) (z + b)}$$ It's clearly not invariant if you take $\tau \rightarrow \tau +2$, which leave $q$ unchanged. The problem, if we insist on $q$ input, is that after restricting $\Re \tau$ to be in $(0,2)$, we have no way to know how much to translate $z$ in order to make sure $a \tau + z$ changes integer multiples of $2$.

Shall we change to $\tau$ input for theta function with characteristics?

[simonp0420]

This is @stla's call. In the meantime, I'm really not at all familiar with theta functions with characteristics. Is there a basic introduction available anywhere? Google search didn't turn up much.

[hammerfunctor]

The best introduction to the family of theta functions I've ever seen is still Tata lectures on theta (I-III) by David Mumford. For theta function with characteristics, check out first five sections of Tata I.

[stla]

Here is a good book introducing the theta function with characteristics. bkcorrect.zip

[stla]

I'm ok for changing jtheta_ab to a function of tau.

[simonp0420]

Thanks to you both for providing the references. I'll take a look at these in the next few days. 😄

[simonp0420]

My training in formal mathematics did not get as far as group theory or functional analysis, so I'm afraid both references are impenetrable to me. I defer to both of your judgements as to making jtheta_ab a function of tau rather than q.

[stla]

Hello @hammerfunctor

You wrote (or is it me?):

a, b, z, τ = promote(a, b, z, τ)

And I remember that @simonp0420 told me we should not do that, I don't remember why. For example in CarlsonRF, @simonp0420 did:

xx, yy, zz, _ = promote(x, y, z, 1.0)

@simonp0420 Is there a problem with a, b, z, τ = promote(a, b, z, τ)? Could you remind me why you didn't use x, y, z, _ = promote(x, y, z, 1.0)?

[simonp0420]

Could you remind me why you didn't use x, y, z, _ = promote(x, y, z, 1.0)?

To avoid type instability.

[hammerfunctor]

To impose type stability, we need to promote values with 1.0im in a lot of places due to factor exp(x * im ) everywhere.

[simonp0420]

I actually haven't studied the entire code base. Up to now I've only tried to help out with some of the Jacobi elliptic functions $\theta_n$, where I had (limited) previous experience and where we've verified via testing that they work for extended precision. In those cases I tried to clean up the code to prevent type instability.