Return to Answer

Your space can be considered as sections of a complex line bundle over the torus. Note that the usual partial derivatives $\partial_1,\partial_2$ do not preserve it, but the operators $$ D_1 = \partial_1 + i\pi x_2,D_2 = \partial_2 - i\pi x_1 $$ do — they define a (unitary) connection on your line bundle, which essentially means that $D_i(fg) = (\partial_i f)g + f(D_i g)$ for $f$ periodic and $g$ in your subspace. Now we have $$ [D_1,D_2] = -2i\pi $$ so these "connection partial derivatives" no longer commute. This calculation immediately shows that there can be no smooth function $g$ of constant norm $1$ in your subspace; if that were the case, we could set $D_ig = A_ig$ for some periodic functions $A_i$, which yields \begin{align*} (2\pi i)g &= D_2D_1 g - D_1D_2 g\\ &= D_2(A_1 g) - D_1(A_2 g)\\ &= (\partial_2 A_1 + A_1A_2 - \partial_1 A_2 -A_1A_2)g\\ \implies 2\pi i &= \partial_2 A_1 - \partial_1 A_2 \end{align*} But this is impossible since the integral of a partial derivative of a periodic function vanishes. (This proof essentially uses that $[D_1,D_2]$, the curvature of the connection, represents a nontrivial de Rham cohomology class on the torus.)

The main advantage of the Fourier basis is that it diagonalizes the partial derivative operators. This doesn't work for the $D_i$ since they satisfy the canonical commutation relations (their commutator is a non-zero multiple of the identity); this is essentially the Heisenberg uncertainty relation. However, the CCR have a unique (irreducible) representation, and from this we obtain the existence of an isometry from your space to $L^2(\mathbb R)$ which carries $D_1$ to multiplication by $2\pi i x$ and $D_2$ to $-\partial_x$.

In fact, one can write this isomorphism down quite explicitly: the operator $-D_1^2 - D_2^2 = (D_1 -iD_2)^\dagger(D_1 -iD_2) +2\pi$ has smallest eigenvalue $2\pi$ with one-dimensional eigenspace spanned by a function $f_0$ such that $D_1 f_0 = iD_2 f_0$. One can then act by the "creation operator" $(D_1 -iD_2)$; setting $f_n = (D_1 -iD_2)^n f_0$, these functions form an orthonormal basis in which the action of the $D_i$ is tridiagonal (the resulting matrix vanishes except for the two diagonals around the main diagonal).

It remains to solve the equation $D_1 f = iD_2 f$ with your boundary conditions. Without boundary conditions, the general solution on $\mathbb R^2\cong \mathbb C$ is of the form $e^{-\pi\lvert z\rvert^2/2}g(z)$ with $g$ holomorphic. The boundary conditions become \begin{align*} g(z+1) & {}= e^{-\pi z +\pi/2}g(z), \\ g(z+i) & {}= e^{i\pi z +\pi/2}g(z). \end{align*} This is the functional equation of the Weierstrass sigma function of the square lattice.

Of course, this basis is probably not that useful for numerics. I would try to discretize the annihilation and creation operators, find the kernel of the former and act on it by the latter; in the resulting basis the operators $D_i$ should have small off-diagonal terms.

Your space can be considered as sections of a complex line bundle over the torus. Note that the usual partial derivatives $\partial_1,\partial_2$ do not preserve it, but the operators $$ D_1 = \partial_1 - i\pi x_2,D_2 = \partial_2 + i\pi x_1 $$ do — they define a (unitary) connection on your line bundle, which essentially means that $D_i(fg) = (\partial_i f)g + f(D_i g)$ for $f$ periodic and $g$ in your subspace. Now we have $$ [D_1,D_2] = 2i\pi $$ so these "connection partial derivatives" no longer commute. This calculation immediately shows that there can be no smooth function $g$ which is nowhere vanishing in your subspace; if that were the case, we could set $D_ig = A_ig$ for some periodic functions $A_i$, which yields \begin{align*} (2\pi i)g &= D_1D_2 g - D_2D_1 g\\ &= D_1(A_2 g) - D_2(A_1 g)\\ &= (\partial_1 A_2 + A_1A_2 - \partial_2 A_1 -A_1A_2)g\\ \implies 2\pi i &= \partial_1 A_2 - \partial_2 A_1 \end{align*} But this is impossible since the integral of a partial derivative of a periodic function vanishes. (This proof essentially uses that $[D_1,D_2]\mathrm dx_1\mathrm dx_2$, the curvature of the connection, represents a nontrivial de Rham cohomology class on the torus, so that the line bundle is nontrivial, while such a function $g$ would determine a trivialization.)

The main advantage of the Fourier basis is that it diagonalizes the partial derivative operators. This doesn't work for the $D_i$ since they satisfy the canonical commutation relations (their commutator is a non-zero multiple of the identity); this is essentially the Heisenberg uncertainty relation. However, the CCR have a unique (irreducible) representation, and from this we obtain the existence of an isometry from your space to $L^2(\mathbb R)$ which carries $D_1$ to multiplication by $2\pi i x$ and $D_2$ to $-\partial_x$.

In fact, one can write this isomorphism down quite explicitly: the operator $-D_1^2 - D_2^2 = (D_1 -iD_2)^\dagger(D_1 -iD_2) +2\pi$ has smallest eigenvalue $2\pi$ with one-dimensional eigenspace spanned by a function $f_0$ such that $D_1 f_0 = iD_2 f_0$. One can then act by the "creation operator" $(D_1 -iD_2)$; setting $f_n = (D_1 -iD_2)^n f_0$, these functions form an orthogonal basis in which the action of the $D_i$ is tridiagonal (the resulting matrix vanishes except for the two diagonals around the main diagonal).

It remains to solve the equation $D_1 f = iD_2 f$ with your boundary conditions. Without boundary conditions, the general solution on $\mathbb R^2\cong \mathbb C$ is of the form $e^{-\pi\lvert z\rvert^2/2}g(z)$ with $g$ holomorphic. The boundary conditions become \begin{align*} g(z+1) & {}= e^{-\pi z +\pi/2}g(z), \\ g(z+i) & {}= e^{i\pi z +\pi/2}g(z). \end{align*} This is the functional equation of the Weierstrass sigma function of the square lattice.

Of course, this basis is probably not that useful for numerics. I would try to discretize the annihilation and creation operators, find the kernel of the former and act on it by the latter; in the resulting basis the operators $D_i$ should have small off-diagonal terms.

Your space can be considered as sections of a complex line bundle over the torus. Note that the usual partial derivatives $\partial_1,\partial_2$ do not preserve it, but the operators $$ D_1 = \partial_1 + i\pi x_2,D_2 = \partial_2 - i\pi x_1 $$ do - they define a (unitary) connection on your line bundle, which essentially means that $D_i(fg) = (\partial_i f)g + f(D_i g)$ for $f$ periodic and $g$ in your subspace. Now we have $$ [D_1,D_2] = -2i\pi $$ so these "connection partial derivatives" no longer commute. This calculation immediately shows that there can be no smooth function $g$ of constant norm $1$ in your subspace; if that was the case, we could set $D_ig = A_ig$ for some periodic functions $A_i$, which yields \begin{align*} (2\pi i)g &= D_2D_1 g - D_1D_2 g\\ &= D_2(A_1 g) - D_1(A_2 g)\\ &= (\partial_2 A_1 + A_1A_2 - \partial_1 A_2 -A_1A_2)g\\ \implies 2\pi i &= \partial_2 A_1 - \partial_1 A_2 \end{align*} But this is impossible since the integral of a partial derivative of a periodic function vanishes. (This proof essentially uses that $[D_1,D_2]$, the curvature of the connetion, represents a nontrivial de Rham cohomology class on the torus.)

The main advantage of the Fourier basis is that it diagonalizes the partial derivative operators. This doesn't work for the $D_i$ since they satisfy the canonical commutation relations (their commutator is proportional to the identity); this is essentially the Heisenberg uncertainty relation. However, the CCR have a unique (irreducible) representation, and from this we obtain the existence of an isometry from your space to $L^2(\mathbb R)$ which carries $D_1$ to multiplication by $2\pi i x$ and $D_2$ to $-\partial_x$.

In fact, one can write this isomorphism down quite explicitly: the operator $-D_1^2 - D_2^2 = (D_1 -iD_2)^\dagger(D_1 -iD_2) +2\pi$ has smallest eigenvalue $2\pi$ with one-dimensional eigenspace spanned by a function $f_0$ such that $D_1 f = iD_2 f$. One can then act by the "creation operator" $(D_1 -iD_2)$; setting $f_n = (D_1 -iD_2)^n f_0$, these functions form an orthonormal basis in which the action of the $D_i$ is tridiagonal (the resulting matrix vanishes except for the two diagonals around the main diagonal).

It remains to solve the equation $D_1 f = iD_2 f$ with your boundary conditions. Without boundary conditions, the general solution on $\mathbb R^2\cong \mathbb C$ is of the form $e^{-\pi|z|^2/2}g(z)$ with $g$ holomorphic. The boundary conditions become \begin{align*} g(z+1) = e^{-\pi z +\pi/2}g(z), g(z+i) = e^{i\pi z +\pi/2}g(z) \end{align*} This is the functional equation of the Weierstrass sigma function of the square lattice, compare https://mathworld.wolfram.com/WeierstrassSigmaFunction.html

Of course, this basis is probably not that useful for numerics. I would try to discretize the annihilation and creation operators, find the kernel of the former and act on it by the latter; in the resulting basis the operators $D_i$ should have small off-diagonal terms.

Your space can be considered as sections of a complex line bundle over the torus. Note that the usual partial derivatives $\partial_1,\partial_2$ do not preserve it, but the operators $$ D_1 = \partial_1 + i\pi x_2,D_2 = \partial_2 - i\pi x_1 $$ do — they define a (unitary) connection on your line bundle, which essentially means that $D_i(fg) = (\partial_i f)g + f(D_i g)$ for $f$ periodic and $g$ in your subspace. Now we have $$ [D_1,D_2] = -2i\pi $$ so these "connection partial derivatives" no longer commute. This calculation immediately shows that there can be no smooth function $g$ of constant norm $1$ in your subspace; if that were the case, we could set $D_ig = A_ig$ for some periodic functions $A_i$, which yields \begin{align*} (2\pi i)g &= D_2D_1 g - D_1D_2 g\\ &= D_2(A_1 g) - D_1(A_2 g)\\ &= (\partial_2 A_1 + A_1A_2 - \partial_1 A_2 -A_1A_2)g\\ \implies 2\pi i &= \partial_2 A_1 - \partial_1 A_2 \end{align*} But this is impossible since the integral of a partial derivative of a periodic function vanishes. (This proof essentially uses that $[D_1,D_2]$, the curvature of the connection, represents a nontrivial de Rham cohomology class on the torus.)

The main advantage of the Fourier basis is that it diagonalizes the partial derivative operators. This doesn't work for the $D_i$ since they satisfy the canonical commutation relations (their commutator is a non-zero multiple of the identity); this is essentially the Heisenberg uncertainty relation. However, the CCR have a unique (irreducible) representation, and from this we obtain the existence of an isometry from your space to $L^2(\mathbb R)$ which carries $D_1$ to multiplication by $2\pi i x$ and $D_2$ to $-\partial_x$.

In fact, one can write this isomorphism down quite explicitly: the operator $-D_1^2 - D_2^2 = (D_1 -iD_2)^\dagger(D_1 -iD_2) +2\pi$ has smallest eigenvalue $2\pi$ with one-dimensional eigenspace spanned by a function $f_0$ such that $D_1 f_0 = iD_2 f_0$. One can then act by the "creation operator" $(D_1 -iD_2)$; setting $f_n = (D_1 -iD_2)^n f_0$, these functions form an orthonormal basis in which the action of the $D_i$ is tridiagonal (the resulting matrix vanishes except for the two diagonals around the main diagonal).

It remains to solve the equation $D_1 f = iD_2 f$ with your boundary conditions. Without boundary conditions, the general solution on $\mathbb R^2\cong \mathbb C$ is of the form $e^{-\pi\lvert z\rvert^2/2}g(z)$ with $g$ holomorphic. The boundary conditions become \begin{align*} g(z+1) & {}= e^{-\pi z +\pi/2}g(z), \\ g(z+i) & {}= e^{i\pi z +\pi/2}g(z). \end{align*} This is the functional equation of the Weierstrass sigma function of the square lattice.

Of course, this basis is probably not that useful for numerics. I would try to discretize the annihilation and creation operators, find the kernel of the former and act on it by the latter; in the resulting basis the operators $D_i$ should have small off-diagonal terms.

Your space can be considered as sections of a complex line bundle over the torus. Note that the usual partial derivatives $\partial_1,\partial_2$ do not preserve it, but the operators $$ D_1 = \partial_1 + i\pi x_2,D_2 = \partial_2 - i\pi x_1 $$ do - they define a (unitary) connection on your line bundle, which essentially means that $D_i(fg) = (\partial_i f)g + f(D_i g)$ for $f$ periodic and $g$ in your subspace. Now we have $$ [D_1,D_2] = -2i\pi $$ so these "connection partial derivatives" no longer commute. This calculation immediately shows that there can be no smooth function $g$ of constant norm $1$ in your subspace; if that was the case, we could set $D_ig = A_ig$ for some periodic functions $A_i$, which yields \begin{align*} (2\pi i)g &= D_2D_1 g - D_1D_2 g\\ &= D_2(A_1 g) - D_1(A_2 g)\\ &= (\partial_2 A_1 + A_1A_2 - \partial_1 A_2 -A_1A_2)g\\ \implies 2\pi i &= \partial_2 A_1 - \partial_1 A_2 \end{align*} But this is impossible since the integral of a partial derivative of a periodic function vanishes. (This proof essentially uses that $[D_1,D_2]$, the curvature of the connetion, represents a nontrivial de Rham cohomology class on the torus.)

The main advantage of the Fourier basis is that it diagonalizes the partial derivative operators. This doesn't work for the $D_i$ since they satisfy the canonical commutation relations (their commutator is proportional to the identity); this is essentially the Heisenberg uncertainty relation. However, the CCR have a unique (irreducible) representation, and from this we obtain the existence of an isometry from your space to $L^2(\mathbb R)$ which carries $D_1$ to multiplication by $2\pi i x$ and $D_2$ to $-\partial_x$.

In fact, one can write this isomorphism down quite explicitly: the operator $-D_1^2 - D_2^2 = (D_1 -iD_2)^\dagger(D_1 -iD_2) +2\pi$ has smallest eigenvalue $2\pi$ with one-dimensional eigenspace spanned by a function $f_0$ such that $D_1 f = iD_2 f$. One can then act by the "creation operator" $(D_1 -iD_2)$; setting $f_n = (D_1 -iD_2)^n f_0$, these functions form an orthonormal basis in which the action of the $D_i$ is tridiagonal (the resulting matrix vanishes except for the two diagonals around the main diagonal).

It remains to solve the equation $D_1 f = iD_2 f$ with your boundary conditions. Without boundary conditions, the general solution on $\mathbb R^2\cong \mathbb C$ is of the form $e^{-\pi|z|^2/2}g(z)$ with $g$ holomorphic. The boundary conditions become \begin{align*} g(z+1) = e^{\pi z +\pi}, g(z+i) = e^{i\pi z +\pi} \end{align*} This is the functional equation of the Weierstrass sigma function of the square lattice, compare https://en.m.wikipedia.org/wiki/Weierstrass_functions (the constants are probably off, I'll try to come back and fix them.)

Of course, this basis is probably not that useful for numerics. I would try to discretize the annihilation and creation operators, find the kernel of the former and act on it by the latter; in the resulting basis the operators $D_i$ should have small off-diagonal terms.

Your space can be considered as sections of a complex line bundle over the torus. Note that the usual partial derivatives $\partial_1,\partial_2$ do not preserve it, but the operators $$ D_1 = \partial_1 + i\pi x_2,D_2 = \partial_2 - i\pi x_1 $$ do - they define a (unitary) connection on your line bundle, which essentially means that $D_i(fg) = (\partial_i f)g + f(D_i g)$ for $f$ periodic and $g$ in your subspace. Now we have $$ [D_1,D_2] = -2i\pi $$ so these "connection partial derivatives" no longer commute. This calculation immediately shows that there can be no smooth function $g$ of constant norm $1$ in your subspace; if that was the case, we could set $D_ig = A_ig$ for some periodic functions $A_i$, which yields \begin{align*} (2\pi i)g &= D_2D_1 g - D_1D_2 g\\ &= D_2(A_1 g) - D_1(A_2 g)\\ &= (\partial_2 A_1 + A_1A_2 - \partial_1 A_2 -A_1A_2)g\\ \implies 2\pi i &= \partial_2 A_1 - \partial_1 A_2 \end{align*} But this is impossible since the integral of a partial derivative of a periodic function vanishes. (This proof essentially uses that $[D_1,D_2]$, the curvature of the connetion, represents a nontrivial de Rham cohomology class on the torus.)

The main advantage of the Fourier basis is that it diagonalizes the partial derivative operators. This doesn't work for the $D_i$ since they satisfy the canonical commutation relations (their commutator is proportional to the identity); this is essentially the Heisenberg uncertainty relation. However, the CCR have a unique (irreducible) representation, and from this we obtain the existence of an isometry from your space to $L^2(\mathbb R)$ which carries $D_1$ to multiplication by $2\pi i x$ and $D_2$ to $-\partial_x$.

In fact, one can write this isomorphism down quite explicitly: the operator $-D_1^2 - D_2^2 = (D_1 -iD_2)^\dagger(D_1 -iD_2) +2\pi$ has smallest eigenvalue $2\pi$ with one-dimensional eigenspace spanned by a function $f_0$ such that $D_1 f = iD_2 f$. One can then act by the "creation operator" $(D_1 -iD_2)$; setting $f_n = (D_1 -iD_2)^n f_0$, these functions form an orthonormal basis in which the action of the $D_i$ is tridiagonal (the resulting matrix vanishes except for the two diagonals around the main diagonal).

It remains to solve the equation $D_1 f = iD_2 f$ with your boundary conditions. Without boundary conditions, the general solution on $\mathbb R^2\cong \mathbb C$ is of the form $e^{-\pi|z|^2/2}g(z)$ with $g$ holomorphic. The boundary conditions become \begin{align*} g(z+1) = e^{-\pi z +\pi/2}g(z), g(z+i) = e^{i\pi z +\pi/2}g(z) \end{align*} This is the functional equation of the Weierstrass sigma function of the square lattice, compare https://mathworld.wolfram.com/WeierstrassSigmaFunction.html

Of course, this basis is probably not that useful for numerics. I would try to discretize the annihilation and creation operators, find the kernel of the former and act on it by the latter; in the resulting basis the operators $D_i$ should have small off-diagonal terms.