A Computational Theta Problem
Here is a computational problem given to me by Vinayak Vatsal. Given the Riemann theta function
and related Riemann theta function with characteristic $[\alpha, \beta]$
compute $\theta[\sigma,0](0, 2\tau)$ for $\tau = zA, zB$ and various $z \in \mathbb{C}, \text{Im}(z) > 0$ with $\sigma \in {0, 1/2}^3$. The definitions of the matrices $A,B \in \mathbb{C}^{3 \times 3}$ are
Below are some computations done in an iPython notebook using abelfunctions. You can download the notebook itself here, which requires that abelfunctions be installed in order to execute.
The Setup
In [1]:
import numpy
import matplotlib
import matplotlib.pyplot as plt
from mpl_toolkits.axes_grid1 import ImageGrid
import matplotlib.colors
from abelfunctions import RiemannThetaIn [2]:
A = numpy.matrix(
[[2, 1, 1],
[1, 4, 1],
[1, 1, 6]],
dtype = numpy.complex)
B = numpy.matrix(
[[23, -5, -3],
[-5, 11, -1],
[-3, -1, 7]],
dtype = numpy.complex)In [3]:
# construct the funciton we want to compute the values of:
#
# theta[sigma,0](0,zM)
#
# where sigma, z, and M can vary
def vitsal(z, sigma, mat):
# construct the Riemann (period) matrix tau and 2tau
tau = z*mat
tau2 = 2*tau
# compute theta with characteristic (RiemannTheta.characteristic()
# is untested.)
dot = numpy.dot(sigma.T, numpy.dot(tau2, sigma))
scale = numpy.exp(2*numpy.pi*1.0j * dot).item()
arg = numpy.dot(tau2, sigma)
return scale * RiemannTheta(arg, tau2)
# vectorize over z since we want to see the behavior of
# this function for some varying z
vitsal = numpy.vectorize(vitsal, excluded=['sigma', 'mat'])Fixed $z \in \mathbb{C}$, All Characteristics $\sigma \in {0,1/2}$
We (arbitrarily) choose $z = i$ and $z = 2i$. Although, the plots in the following section suggests that the choice of $z$ in the upper half plane doesn’t matter except, possibly, in convergence rate of the implementation of the Riemann theta function to the suggested value of the function. That is, it seems like the truncated sum of the Riemann theta function converges to the actual value as we add more terms faster for $\text{Im}(z)$ large.
In [4]:
def print_values(z, mat):
# constuct list of all half integer characteristics
halfs = [0,0.5]
sigmas = []
for s1 in halfs:
for s2 in halfs:
for s3 in halfs:
sigmas.append(numpy.matrix([s1,s2,s3]).T)
# for each characteristic, print characteristic and
# value of theta function
for sigma in sigmas:
val = vitsal(z, sigma=sigma, mat=mat)
print 'sigma.T: %s\tvalue: %s'%(sigma.T, val)In [5]:
print_values(1.0j, A)Out[5]:
sigma.T: [[0 0 0]] value: (1.00000697473+0j)
sigma.T: [[ 0. 0. 0.5]] value: (6.52457377777e-09+0j)
sigma.T: [[ 0. 0.5 0. ]] value: (6.9877095594e-06+0j)
sigma.T: [[ 0. 0.5 0.5]] value: (2.2795992109e-14+0j)
sigma.T: [[ 0.5 0. 0. ]] value: (0.00373488548776+0j)
sigma.T: [[ 0.5 0. 0.5]] value: (1.21842677994e-11+0j)
sigma.T: [[ 0.5 0.5 0. ]] value: (1.30369858713e-08+0j)
sigma.T: [[ 0.5 0.5 0.5]] value: (1.58402415218e-19+0j)
In [6]:
print_values(2.0j, A)Out[6]:
sigma.T: [[0 0 0]] value: (1.00000000002+0j)
sigma.T: [[ 0. 0. 0.5]] value: (4.24116597336e-17+0j)
sigma.T: [[ 0. 0.5 0. ]] value: (1.21615991209e-11+0j)
sigma.T: [[ 0. 0.5 0.5]] value: (5.15790006267e-28+0j)
sigma.T: [[ 0.5 0. 0. ]] value: (6.97468471242e-06+0j)
sigma.T: [[ 0.5 0. 0.5]] value: (1.47903977386e-22+0j)
sigma.T: [[ 0.5 0.5 0. ]] value: (4.24116597336e-17+0j)
sigma.T: [[ 0.5 0.5 0.5]] value: (6.27280941128e-39+0j)
In [7]:
print_values(1.0j, B)Out[7]:
sigma.T: [[0 0 0]] value: (1+0j)
sigma.T: [[ 0. 0. 0.5]] value: (2.81426845749e-10+0j)
sigma.T: [[ 0. 0.5 0. ]] value: (9.81431759353e-16+0j)
sigma.T: [[ 0. 0.5 0.5]] value: (1.47903461596e-22+0j)
sigma.T: [[ 0.5 0. 0. ]] value: (4.16240046723e-32+0j)
sigma.T: [[ 0.5 0. 0.5]] value: (1.79873633572e-33+0j)
sigma.T: [[ 0.5 0.5 0. ]] value: (1.79873633572e-33+0j)
sigma.T: [[ 0.5 0.5 0.5]] value: (4.16240046723e-32+0j)
In [8]:
print_values(2.0j, B)Out[8]:
sigma.T: [[0 0 0]] value: (1+0j)
sigma.T: [[ 0. 0. 0.5]] value: (7.9201069508e-20+0j)
sigma.T: [[ 0. 0.5 0. ]] value: (9.63208298267e-31+0j)
sigma.T: [[ 0. 0.5 0.5]] value: (2.18754339521e-44+0j)
sigma.T: [[ 0.5 0. 0. ]] value: (1.73255776496e-63+0j)
sigma.T: [[ 0.5 0. 0.5]] value: (3.23545240544e-66+0j)
sigma.T: [[ 0.5 0.5 0. ]] value: (3.23545240544e-66+0j)
sigma.T: [[ 0.5 0.5 0.5]] value: (1.73255776496e-63+0j)
Varying $z \in \mathfrak{h}_1$, $\sigma = [0,0,0]$ and $\sigma = [1/2,0,0]$
We now draw contour plots for the real, imaginary, and absolute value of the theta function for various $z$ in the upper half plane. We chose $\sigma = [0,0,0]$ because it was the only characteristic resulting in a theta value of 1. We chose $\sigma = [1/2,0,0]$ because it seems to be the one most sensitive to changing $z$ and the furthest away from zero for small values of $z$.
In [9]:
# computing the value of the theta function for various z
# set the sigma, z-grid size, and plotting window
sigma = numpy.matrix([0,0,0]).T
N = 64
zre_left, zre_right = -0.2,0.2
zim_left, zim_right = 1,1.2
# compute theta for all of these z
x = numpy.linspace(zre_left, zre_right, N)
y = numpy.linspace(zim_left, zim_right, N)
X,Y = numpy.meshgrid(x,y)
Z = X + 1.0j*Y
V = vitsal(Z, sigma=sigma, mat=A)In [10]:
# a fancy plotting function so we can see the real, imaginary, and
# absolute values with the same colorbar. took more time to write this
# than the rest of the code in this notebook
def fancy_plot(X, Y, V, figwidth=10):
# (... omitted for brevity ...)In [11]:
fancy_plot(X, Y, V)Out[11]:
In [12]:
# same as above, but with sigma = [0.5,0,0]
# set the sigma, z-grid size, and plotting window
sigma = numpy.matrix([0.5,0,0]).T
N = 64
zre_left, zre_right = -1,1
zim_left, zim_right = 1,2
# compute theta for all of these z
x = numpy.linspace(zre_left, zre_right, N)
y = numpy.linspace(zim_left, zim_right, N)
X,Y = numpy.meshgrid(x,y)
Z = X + 1.0j*Y
V = vitsal(Z, sigma=sigma, mat=A)
# plot
fancy_plot(X, Y, V)Out[12]:
