\[section\] \[lem\][Theorem]{}

Introduction
============

In the first two papers of this series ([@HSI; @HSII], henceforth referred to as I and II), we initiated a detailed study of spectral properties of operators arising in the context of elliptic Calogero–Moser type systems. The operators at issue are not only the Hamiltonians defining these systems, but also families of Hilbert–Schmidt (HS) operators that seem unrelated at first sight. The relation consists in the kernels of the HS operators being eigenfunctions of differences of the Hamiltonians. As we have first detailed in our lecture notes [@EIS2] and also sketched in I, this relation can be exploited to derive novel information about the Hamiltonians, viewed as Hilbert space operators (as opposed to differential or difference operators).

The present paper is concerned with the Heun Hamiltonian $$\begin{gathered}
\label{Hg}
 H(g;x) \equiv -\frac{d^2}{dx^2}+V(g;x),  \qquad g=(g_0,g_1,g_2,g_3)\in{{\mathbb R}}^4,
\\
\label{Vg}
V(g;x)\equiv  \sum_{t=0}^3g_t(g_t-1)\wp(x+\omega_t;\pi/2r,i\alpha/2),\qquad \ r,\alpha>0,
 \end{gathered}$$ where $$\begin{gathered}
\omega_0\equiv 0,\qquad \omega_1\equiv \pi/2r,\qquad \omega_2\equiv i\alpha/2,\qquad \omega_3\equiv -\omega_1-\omega_2.\end{gathered}$$ (For a lucid account of elliptic and allied functions we recommend [@ww], whose conventions and notation we mostly follow.) It can also be viewed as the $BC_1$ ‘nonrelativistic’ elliptic Calogero–Moser Hamiltonian, the $BC_N$ case being the Inozemtsev system [@inoz]. There exists an extensive lore about the eigenfunctions of $H(g)$ viewed as a differential operator, cf. e.g. the reviews . In recent years, this lore has been significantly enlarged by Takemura, whose work can be traced from his recent paper [@take].

Here, however, we are concerned with spectral properties of $H(g)$, reinterpreted as a self-adjoint operator on the Hilbert space $$\begin{gathered}
\label{cH}
{{\cal H}}\equiv L^2([0,\pi/2r],dx).\end{gathered}$$ A few of the properties we derive (such as the occurrence of solely discrete spectrum for suitably restricted $g$) follow from well-known arguments in the theory of Schrödinger operators [@rs2] (or, alternatively, from Sturm–Liouville theory and the theory of eigenfunction expansions). Our main tool for reobtaining known features and obtaining novel ones is the integral operator whose kernel is given by $$\begin{gathered}
\label{Psi}
 \Psi(g;x,y)\equiv w(g;x)^{1/2}{{\cal S}}(g;x,y)w(g';y)^{1/2},\qquad x,y\in(0,\pi/2r).
 \end{gathered}$$

This kernel and the definition of its ingredients can be found in I, cf. I Subsection 3.2. For ease of reference, however, we repeat the relevant definitions: $$\begin{gathered}
\label{defJN}
g'\equiv J_Ng,\qquad J_N\equiv
 \frac{1}{2} \left( \begin{array}{rrrr}
1  &  1  &  -1  &  1  \\
1  &  1  &  1  &  -1  \\
-1  &  1  &  1  &  1  \\
1  &  -1  &  1  &  1
\end{array} \right) ,
\\
\label{cS}
{{\cal S}}(g;x,y)\equiv \exp(-s_g\ln [R(x+y)R(x-y)]),
\\
\label{sg}
s_g\equiv \frac{1}{2} \sum_{t=0}^3g_t=\frac{1}{2}\sum_{t=0}^3g'_t.
\\
\label{w}
w(g;x)\equiv 1/c(g;x)c(g;-x),
\\
  c(g;x)\equiv R(x+i\alpha/2)^{-g_0}R(x+i\alpha/2-\pi/2r)^{-g_1}R(x)^{-g_2}R(x-\pi/2r)^{-g_3},\nonumber
\\
R(z)\equiv \prod_{l=0}^{\infty}(1-\exp[-(2l+1)r\alpha+2irz])(1-\exp[-(2l+1)r\alpha-2irz]),\qquad z\in{{\mathbb C}}.\nonumber\end{gathered}$$

At this point we would like to mention that a close relative of the kernel ${{\cal S}}(g;x,y)$ has appeared in the literature before, a fact we only recently became aware of. Indeed, it can be tied in with the kernel that can be found on p. 124 of a monograph by Slavyanov and Lay [@slla]. (It seems this kernel was first introduced in a paper by Kazakov and Slavyanov [@kasl].) More precisely, when a $\wp$-function substitution is made in the latter kernel, it can be factorized in terms of the Weierstrass $\sigma$-function; the nontrivial part is then an analytic continuation of ${{\cal S}}(g;x,y)$. Both in [@slla] and in Takemura’s paper [@take2], this kernel is used to relate Heun functions defined on different intervals.

The entire function $R(z)$ plays a crucial role in this paper. It satisfies $$\begin{gathered}
\label{R2}
R(z)=\exp\left( -\sum_{n=1}^{\infty}\frac{\cos(2nrz)}{n\sinh(nr\alpha)}\right),\qquad |\Im z|<\alpha/2,
\\
\label{dupl}
R(2z)=\prod_{{\delta}=+,-}R(z+i{\delta}\alpha/4)R(z+\pi/2r+i{\delta}\alpha/4),
\\
\label{Rs}
R(z+i\alpha/2)=ipe^{-irz}s(z),\qquad p=2r\prod_{k=1}^{\infty}(1-e^{-2kr\alpha})^2.\end{gathered}$$ Here, $s(z)$ is related to the $\sigma$-function by $$\begin{gathered}
\label{ssig}
s(z)=\exp\left(-\eta z^2r/\pi\right)\sigma(z;\pi/2r,i\alpha/2),\end{gathered}$$ so that $$\begin{gathered}
f(z+i\alpha)=-\exp(r\alpha -2irz)f(z),\qquad f=s,R.\end{gathered}$$ (For more details, see [@fo], where the functions $s$ and $R$ were introduced.)

From (\[Rs\]) it follows that the weight function (\[w\]) can also be written $$\begin{gathered}
\label{w2}
  w(g;x)=p^{2g_0+2g_1}s(x)^{2g_0}s(\pi/2r-x)^{2g_1}R(x)^{2g_2}R(\pi/2r-x)^{2g_3}.\end{gathered}$$ It is positive on $(0,\pi/2r)$ for all $g\in{{\mathbb R}}^4$, but integrable on $[0,\pi/2r]$ only for $g_0,g_1>-1/2$.

It now readily follows that for $\Psi(g;x,y)$ to belong to the kernel Hilbert space $$\begin{gathered}
\label{cHK}
 {{\cal H}}_K\equiv L^2\left([0,\pi/2r]^2,dxdy\right),
 \end{gathered}$$ the necessary and sufficient condition is that $g$ belong to the parameter set $$\begin{gathered}
\label{Pi}
 \Pi \equiv \left\{ g\in{{\mathbb R}}^4\mid g_0,g_1,g_0',g_1'>-1/2\right\}.
 \end{gathered}$$ Hence it is clear that the integral operator $$\begin{gathered}
\label{cI}
 ({{\cal I}}(g) f)(x)\equiv {\int_{0}^{\pi/2r}}dy\Psi(g;x,y)f(y),\qquad f\in{{\cal H}},\qquad
g\in\Pi,\end{gathered}$$ is well defined and HS.

It is not obvious, but true that its kernel is a zero-eigenvalue eigenfunction of a difference of Heun Hamiltonians: $$\begin{gathered}
\label{eiB}
(H(g;x)-H(g';y))\Psi(g;x,y)=0,\qquad x,y\in(0,\pi/2r).\end{gathered}$$ (Indeed, this is a corollary of Proposition 3.4 in I, cf. I(3.115).) At face value, this eigenfunction identity seems to have no bearing on the Hilbert space eigenfunctions of the two (generically) distinct Hamiltonians. A principal result of this paper is, however, that for $g$ in the restricted parameter set $$\begin{gathered}
\label{Pir}
\Pi_r\equiv \{ g\in\Pi \mid s_g>0\},\end{gathered}$$ the Hilbert space eigenvector ONBs (orthonormal bases) of $H(g)$ and $H(g')$ are given by functions $e_n(g;x)$ and $e_n(g';x)$, resp., which yield the singular value decomposition of ${{\cal I}}(g)$ (cf. [@rs1]): $$\begin{gathered}
\label{sval}
{{\cal I}}(g)=\sum_{n=0}^{\infty}\nu_n(g) e_n(g)\otimes e_n(g'),\qquad \nu_0\ge \nu_1\ge \cdots > 0,
\qquad \sum_{n=0}^{\infty}\nu_n^2<\infty,
\\
(e_m,e_n)={\delta}_{mn},\qquad m,n\in{{\mathbb N}}.\nonumber\end{gathered}$$ (All singular values are positive, since ${{\cal I}}(g)$ has trivial kernel and dense range for $g\in\Pi_r$, cf. (\[KerI\]). To appreciate the restriction on $s_g$, note that ${{\cal I}}(g)$ has rank 1 for $s_g=0$, cf. (\[cS\]).)

This result entails in particular that $H(g)$ and $H(g')$ have solely discrete spectrum, the two spectra being equal. To be sure, in the case at hand it follows from well-known Schrödinger operator theory that $H(g)$ can be reinterpreted as a self-adjoint operator that has an ONB of eigenvectors, provided $g$ is suitably restricted. The main novelty is that these eigenvectors are just those occurring for ${{\cal I}}(g)$, and that this fixes the self-adjoint extensions of the Schrödinger operator $H(g)$ in case it is not essentially self-adjoint on $C_0^{\infty}((0,\pi/2r))$ (namely, for $g_0$ and/or $g_1$ in $(-1/2,3/2)$). As a bonus, this yields a hidden spectral invariance under permutations of four new parameters, which implies spectral consequences for special cases that are quite startling. Last but not least, the state of affairs for the present case gives rise to a simple paradigm for the ‘relativistic’ (analytic difference operator) Heun case, where the usual self-adjoint extension theory is of no help and no previous spectral results are known, cf. the fourth paper in this series [@HSIV].

We proceed to sketch the organization and results of this paper in more detail. We begin Section \[section2\] by deriving features of functions in the range of ${{\cal I}}(g)$, assuming $g\in\Pi_r$. Then we show that ${{\cal I}}(g),g\in\Pi_r$, has dense range and trivial kernel. These crucial features readily follow from Lemma \[lemma2.1\], which asserts that the integral operator with kernel ${{\cal S}}(g;x,y)$ has kernel $\{ 0\}$ on $L^1([0,\pi/2r])$ whenever $s_g$ is positive. The characteristics of the range of ${{\cal I}}(g)$ are now used as a guide for a choice of domain for $H(g)$, namely the dense subspace $D(g_0,g_1)$ given by (\[D01\]). Lemma \[lemma2.2\] then asserts that $H(g)$ is symmetric on $D(g_0,g_1)$ for all $g$ in the parameter space $$\begin{gathered}
\label{tPi}
\tilde{\Pi}\equiv \left\{ g\in{{\mathbb R}}^4 \mid g_0,g_1>-1/2\right\}.\end{gathered}$$ The proofs of these lemmas involve technicalities that are not enlightening, so we have relegated them to Appendix \[appendix-A\].

Next, using the key identity (\[eiB\]), we show in Lemma \[lemma2.3\] that for all $g\in\tilde{\Pi}$ the subspace ${D(g_0,g_1)}$ is a core (domain of essential self-adjointness [@rs2]). This prepares us for studying the self-adjoint closure (again denoted $H(g)$) by using ${{\cal I}}(g)$. Only a few simple steps are then needed to reach the conclusion that the above ONB $\{e_n(g)\}$ may be chosen so that it is an eigenvector ONB for $H(g)$. Moreover, since the eigenvectors belong to ${D(g_0,g_1)}$, it easily follows from the Frobenius method that the spectrum of $H(g)$ is simple, cf. Theorem \[theorem2.4\].

Whenever the potential $V(g;x)$ is bounded below, it is obvious that $H(g)$ is also bounded below. For $g_0$ and/or $g_1$ in $ (0,1)$, however, $V(g;x)$ is not bounded below, so an extra argument is needed to show that for the pertinent self-adjoint extensions of (\[Hg\]) semi-boundedness is preserved. We prove this by a comparison argument involving the trigonometric $BC_1$ Calogero–Moser Hamiltonian $$\begin{gathered}
H_t(g_0,g_1;x)\equiv -\frac{d^2}{dx^2}+V_t(g_0,g_1;x),\nonumber\\
V_t(g_0,g_1;x)\equiv \frac{r^2g_0(g_0-1)}{\sin^2rx}+\frac{r^2g_1(g_1-1)}{\cos^2rx}.\label{Ht}\end{gathered}$$ The point is that this operator is well defined and essentially self-adjoint on ${D(g_0,g_1)}$, and that it has an ONB of eigenvectors in ${D(g_0,g_1)}$ involving the well-known Jacobi polynomials, cf. Theorem \[theorem2.5\].

We would like to add that the trigonometric Hamiltonian (\[Ht\]) plays a key role in Takemura’s paper [@TakeHeun], where he studies the Heun Hamiltonian $H(g;x)$ (\[Hg\]) with $g_0,g_1\ge 1$ from a Hilbert space point of view. His strategy is based on a perturbation expansion in the parameter $p=\exp(-\alpha r)$. (In this connection it should be pointed out that $H(g;x)$ with $g_2=g_3=0$ reduces to $H_t(g_0,g_1;x)$ for $\alpha\to\infty$.)

The main result of Theorem \[theorem2.6\] is that for $g\in\Pi_r$ the spectra of $H(g)$ and $H(g')$ coincide. Unfortunately, we have not been able to show that the singular values in (\[sval\]) are distinct for all $g\in\Pi_r$, a feature that seems plausible to us. If this is indeed the case, then it would follow that $e_n(g)$ is the $H(g)$-eigenvector with eigenvalue $E_n(g)$ in the obvious eigenvalue list $$\begin{gathered}
E_0(g)<E_1(g)<E_2(g)<\cdots,\qquad g\in\Pi_r.\end{gathered}$$

In Section \[section3\] we fix $g\in\Pi_r$ and address the question for which $\hat{g}\in\tilde{\Pi}$ the spectrum of $H(\hat{g})$ coincides with that of $H(g)$. From Theorem \[theorem2.6\] we know this is true for $\hat{g}=g'$, but there are a few more involutory $g$-transformations for which this is valid. We determine the group $G$ generated by these involutions in Theorem \[theorem3.1\]. It leaves the parameter set $$\begin{gathered}
\label{PiG}
\Pi_G\equiv \{ g\in\Pi_r\mid g_0+g_1+2>g_2+g_3\}\end{gathered}$$ invariant, and is isomorphic to the permutation group $S_4$. (For ease of reference, we collect the four parameter sets (\[PiG\]), (\[Pir\]), (\[Pi\]) and (\[tPi\]) in a chain $$\begin{gathered}
\Pi_G\subset \Pi_r\subset \Pi\subset\tilde{\Pi},\end{gathered}$$ and note that each inclusion is proper.) Indeed, defining new couplings $$\begin{gathered}
c_0\equiv g_0+g_3-1,\qquad  c_1\equiv g_1+g_2-1,\qquad c_2\equiv g_1-g_2,\qquad  c_3\equiv g_0-g_3,\end{gathered}$$ the spectrum of $$\begin{gathered}
\hat{H}(c)\equiv H(g(c))\end{gathered}$$ is invariant under arbitrary permutations of $c_0$, $c_1$, $c_2$, $c_3$.

In Section \[section4\] we consider various special $g$-choices. For the sixteen $g$-values yielding , all of the relevant quantities can be determined in great detail. A striking consequence of our spectral invariance results is that there exist quite a few $g$-choices for which $V(g;x)\ne 0$, yet the spectrum of $H(g)$ equals that for one of the cases where $V(g;x)=0$. We also consider the case $s_g=0$, for which ${{\cal I}}(g)$ has rank 1. In particular, using results from I and [@cadiz], we verify that ${{\cal I}}(g)$ connects the ground state of $H(g')$ with the one of $H(g)$. Finally, we add some observations on the case where the $\pi/r$-periodic potential $V(g;x)$ has no poles for real $x$, a well-known setting that is for instance studied in considerable detail in [@arsc].

The paper is concluded with Appendix \[appendix-A\], in which we present the proofs of Lemmas \[lemma2.1\] and \[lemma2.2\].

The eigenvector ONBs of $\boldsymbol{H(g)}$ and $\boldsymbol{H(g')}$ {#section2}
====================================================================

Unless explicitly mentioned otherwise, we assume from now on that the coupling vector $g$ belongs to the restricted parameter set $\Pi_r$ (\[Pir\]). This entails in particular that the integral operator ${{\cal I}}(g)$ (\[cI\]) is a well-defined HS operator on the Hilbert space ${{\cal H}}$ (\[cH\]). We begin this section by obtaining features of functions in the range of ${{\cal I}}(g)$.

Consider first for $f\in{{\cal H}}$ and $x\in[0,\pi/2r]$ the function $$\begin{gathered}
h(x)     \equiv    \int_0^{\pi/2r}dy{{\cal S}}(g;x,y)w(g';y)^{1/2}f(y)
\nonumber\\
\phantom{h(x)}{}   =       \int_0^{\pi/2r}dy\exp\left( 2s_g\sum_{n=1}^{\infty}\frac{\cos(2nrx)\cos(2nry)}{n\sinh(nr\alpha)}\right) w(g';y)^{1/2}f(y),\label{hdef}
  \end{gathered}$$ where we used (\[R2\]). Since we have $w(g';\cdot)^{1/2}\in{{\cal H}}$, we deduce from the Schwarz inequality $$\begin{gathered}
  w(g';y)^{1/2}f(y)\in L^1([0,\pi/2r],dy).
  \end{gathered}$$ It readily follows that $h(x)$ extends to a function $h(z)$, $\Re z=x$, that is analytic in the strip $|\Im z|<\alpha/2$. Also, $h(z)$ is $\pi/r$-periodic and even. Now $R(z)$ is an entire, $\pi/r$-periodic and even function without zeros for $|\Im z|<\alpha/2$, so the function $$\begin{gathered}
\label{ana}
  s(x)^{-g_0}s(\pi/2r-x)^{-g_1}({{\cal I}}(g) f)(x)=R(x)^{g_2}R(\pi/2r-x)^{g_3}h(x),\qquad x\in (0,\pi/2r),
  \end{gathered}$$ extends to a function that is analytic for $|\Im z|<\alpha/2$, $\pi/r$-periodic and even.

Next, we point out the relations $$\begin{gathered}
  g\in\Pi_r \ \Leftrightarrow \ g'\in\Pi_r, \nonumber
  \\
  \label{cIad}
  {{\cal I}}(g)^{*}={{\cal I}}(g'),
  \end{gathered}$$ which are easily verified. As a consequence, the orthocomplement of the range of ${{\cal I}}(g)$ equals the kernel of ${{\cal I}}(g')$. From (\[hdef\]) it is immediate that ${{\cal I}}(g), g\in\Pi$, has rank 1 when $s_g$ vanishes. By contrast, for $s_g>0$ the integral operator ${{\cal I}}(g)$ has dense range and trivial kernel: $$\begin{gathered}
\label{KerI}
{\rm Ker}( {{\cal I}}(g))=\{ 0\},\qquad {\rm Ran}({{\cal I}}(g))^{-}={{\cal H}},\qquad \forall \, g\in\Pi_r.\end{gathered}$$ This is clear from the above and the following lemma, which is proved in Appendix \[appendix-A\].

\[lemma2.1\] Letting $s_g>0$, suppose $\phi(y) \in L^1([0,\pi/2r],dy)$ satisfies $$\begin{gathered}
 \int_0^{\pi/2r}dy{{\cal S}}(g;x,y)\phi(y)=0.\end{gathered}$$ Then $\phi =0$.

To proceed, we note that the kernel function $\Psi(g;x,y)$ is real-valued. From (\[KerI\]), (\[cIad\]) and the HS property we then deduce that for a fixed $g\in\Pi_r$, there exist two ONBs $\{ e_n(g')\}_{n\in{{\mathbb N}}}$ and $\{ e_n(g)\}_{n\in{{\mathbb N}}}$ of real-valued functions such that $$\begin{gathered}
\label{Psirep}
\Psi(g;x,y)=\sum_{n=0}^{\infty}\nu_n(g)e_n(g;x)e_n(g';y),\qquad \nu_0\ge\nu_1\ge \cdots >0,\end{gathered}$$ with $\{ \nu_n\}\in l^2({{\mathbb N}})$. (The choice of the first ONB uniquely determines the second one; later on, we will return to the specific choice we make.) Moreover, since $e_n(g)$ belongs to the range of ${{\cal I}}(g)$, the above analysis implies that it is of the form $$\begin{gathered}
\label{ed}
e_n(g;x)=w(g;x)^{1/2}d_n(g;x),\end{gathered}$$ where $d_n(g;x)$ extends to a function that is analytic in $|\Im z|<\alpha/2$, even and $\pi/r$-periodic.

We are now prepared to deal with the Heun Hamiltonian $H(g)$. We choose as initial domain for $H(g)$ the space $$\begin{gathered}
\label{D01}
D(g_0,g_1)=s(x)^{g_0}s(\pi/2r-x)^{g_1}{{\cal D}},\qquad g_0,g_1>-1/2,\end{gathered}$$ where $$\begin{gathered}
\label{cD}
{{\cal D}}\equiv \left\{ f\in C^2([0,\pi/2r])\mid f'(0)=f'(\pi/2r)=0\right\}.\end{gathered}$$ This domain choice ensures that the ONB functions $e_n(g;x)$ belong to the domain. Indeed, from the paragraph containing (\[ana\]) we see that the function on the r.h.s. of (\[ana\]) belongs to ${{\cal D}}$, so that we have more generally $$\begin{gathered}
\label{Ran}
{\rm Ran}\, ({{\cal I}}(g))\subset D(g_0,g_1).\end{gathered}$$

From $e_n(g)\in D(g_0,g_1)$ it already follows that $D(g_0,g_1)$ is dense. But this feature can also be derived quite easily without invoking ${{\cal I}}(g)$. Indeed, ${{\cal D}}$ contains the functions $\cos(2nrx)$, $n\in{{\mathbb N}}$, which are well known to yield an orthogonal base of ${{\cal H}}$. Therefore, the assumption $f\perp D(g_0,g_1)$ entails $$\begin{gathered}
\int_0^{\pi/2r}dxs(x)^{g_0}s(\pi/2r-x)^{g_1}\cos(2nrx)f(x)=0,\qquad \forall \, n\in{{\mathbb N}}.\end{gathered}$$ From the orthogonal base property we now deduce $f(x)=0$, so that $D(g_0,g_1)$ is dense.

We are now prepared for our second lemma, whose proof is relegated to Appendix \[appendix-A\].

\[lemma2.2\] Assume $g$ belongs to $\tilde{\Pi}$ . Then $H(g)$ is well defined and symmetric on the domain $D(g_0,g_1)$ given by , .

The proof of this lemma does not involve the HS operator ${{\cal I}}(g)$. In the next lemma we do use ${{\cal I}}(g)$ to obtain a short proof of essential self-adjointness on $D(g_0,g_1)$.

\[lemma2.3\] Assuming $g\in\tilde{\Pi}$, the domain $D(g_0,g_1)$ is a core for $H(g)$.

We begin by noting that the potentials $\wp(x+i\alpha/2)$ and $\wp(x-\pi/2r-i\alpha/2)$ yield bounded self-adjoint operators, so that they do not influence domain issues. Therefore, we need only prove the core property for one value of $g_2$ and $g_3$. We choose $g_2=g_3=1$ (say), the point being that this entails $g\in\Pi_r$. (This is easily checked, cf. (\[defJN\]), (\[Pi\]) and (\[Pir\]).) Hence we can now use ${{\cal I}}(g)$.

Assume $\phi\in{{\cal H}}$ satisfies $$\begin{gathered}
(\phi,(H(g)+i)\psi)=0,\qquad \forall \, \psi\in D(g_0,g_1).\end{gathered}$$ By virtue of (\[Ran\]), this implies in particular $$\begin{gathered}
(\phi,(H(g)+i){{\cal I}}(g) f)=0,\qquad \forall \, f\in{{\cal H}}.\end{gathered}$$ Thanks to (\[eiB\]) and (\[cIad\]), this inner product can be rewritten as $$\begin{gathered}
((H(g')-i){{\cal I}}(g')\phi,f)=0,\qquad \forall \, f\in{{\cal H}}.\end{gathered}$$ This implies that the vector ${{\cal I}}(g')\phi$, which belongs to $D(g_0^{'},g_1^{'})$ by (\[Ran\]), is an eigenvector with eigenvalue $i$ for the symmetric operator $H(g')$. Thus we must have ${{\cal I}}(g')\phi=0$. By (\[KerI\]), this entails $\phi=0$. Therefore ${\rm Ran}\, (H(g)+i)$ is dense in ${{\cal H}}$. Likewise, ${\rm Ran}\, (H(g)-i)$ is dense in ${{\cal H}}$. Hence the core property results.

We would like to point out that the dense subspace $$\begin{gathered}
\label{cC}
{{\cal C}}\equiv C_0^{\infty}((0,\omega_1))\subset D(g_0,g_1),\end{gathered}$$ is already a core for $H(g)$ whenever $g_0,g_1> 3/2$. This is basically known, but for completeness we sketch a proof. Choosing at first $g\in{{\mathbb R}}^4$, consider the differential operator $H(g)$ restricted to ${{\cal C}}$. This is obviously a well-defined symmetric operator on ${{\cal H}}$, and it is easy to verify that the domain of its closure contains the space $$\begin{gathered}
{{\cal C}}'    \equiv     \big\{ \psi\in C^1([0,\omega_1])\mid \psi(0)=\psi(\omega_1)=\psi'(0)=\psi'(\omega_1)=0,
\nonumber \\
\phantom{{{\cal C}}'    \equiv     \big\{}{}  \psi'(x)\in C^1((0,\omega_1)), \psi''(x), x^{-2}\psi(x),(\omega_1-x)^{-2}\psi(x)\in{{\cal H}}\big\}.\end{gathered}$$ Now assume $\phi$ is orthogonal to $ (H(g)+i){{\cal C}}$. Then $\phi$ is a weak solution to the ODE $(H(g)-i)\phi=0$ on $(0,\omega_1)$, and by hypo-ellipticity it is a classical solution.

Choosing next $g_0,g_1> 3/2$, it follows from the Frobenius theory that near 0 we have $$\begin{gathered}
\phi(x)= x^{g_0}h_0(x),\end{gathered}$$ with $h_0(x)$ analytic at $x= 0$. (Indeed, $x^{1-g_0}$ is not in $L^2$ near 0.) Likewise, near $x=\omega_1$, $$\begin{gathered}
 \phi(x)= (\omega_1-x)^{g_1}h_1(x),
 \end{gathered}$$ with $h_1(x)$ analytic at $x=\omega_1$. But then we deduce that $\phi$ has all of the properties defining ${{\cal C}}'$, and since the closure of $H(g)$ is symmetric on ${{\cal C}}'$, it cannot have an eigenvector with nonreal eigenvalue. Thus $\phi=0$ and the core property follows.

After this digression we proceed to obtain our first main result, for which we again need to require $g\in\Pi_r$. We denote the domain of the self-adjoint closure of $H(g)$ by $D(g)$, but keep the notation $H(g)$ for the closure. We now use the key identity (\[eiB\]) as in the proof of Lemma \[lemma2.3\] to deduce $$\begin{gathered}
\label{key2}
(\phi,H(g){{\cal I}}(g)\psi)=(H(g'){{\cal I}}(g')\phi,\psi),\qquad \forall \, \phi,\psi\in{{\cal H}}.\end{gathered}$$ Choosing $\lambda\in{{\mathbb C}}$ with $\Im \lambda >0$, this entails in particular $$\begin{gathered}
\label{ipH}
(\phi,(H(g)+\lambda){{\cal I}}(g)\psi)=((H(g')+\overline{\lambda}){{\cal I}}(g')\phi,\psi),\qquad \forall \, \phi\in D(g),\qquad \forall \, \psi\in D(g').\end{gathered}$$ Now for $\phi\in D(g)$ and $\psi\in D(g')$ there exist $e,f\in{{\cal H}}$ such that $$\begin{gathered}
\phi=(H(g)+\overline{\lambda})^{-1}e,\qquad \psi=(H(g')+\lambda)^{-1}f,\end{gathered}$$ so (\[ipH\]) implies $$\begin{gathered}
(e,{{\cal I}}(g)(H(g')+\lambda)^{-1}f)=({{\cal I}}(g')(H(g)+\overline{\lambda})^{-1}e,f),\qquad \forall \, e,f\in{{\cal H}}.\end{gathered}$$ Hence we obtain the relation $$\begin{gathered}
{{\cal I}}(g)(H(g')+\lambda)^{-1}=(H(g)+\lambda)^{-1}{{\cal I}}(g).\end{gathered}$$ From this we readily deduce $$\begin{gathered}
[(H(g)+\lambda)^{-1},T(g)]=0,\end{gathered}$$ where $$\begin{gathered}
\label{T}
T(g)\equiv {{\cal I}}(g){{\cal I}}(g)^{*}=\sum_{n=0}^{\infty}\nu_n^2(g)e_n(g)\otimes e_n(g).\end{gathered}$$

The upshot is that the self-adjoint trace class operator $T(g)$ commutes with the $H(g)$-resolvent. Since $T(g)$ has trivial kernel by (\[KerI\]), all of its eigenspaces are finite-dimensional, and the $H(g)$-resolvent leaves these finite-dimensional spaces invariant. Thus we can choose the eigenvectors of $T(g)$ to be eigenvectors of the $H(g)$-resolvent, hence of $H(g)$ as well. We now summarize and extend these findings.

\[theorem2.4\] Let $g\in\Pi_r$. Then the positive trace class operator $T(g)$  has an ONB $\{ e_n(g)\}_{n\in{{\mathbb N}}}$, of eigenvectors that are also eigenvectors of $H(g)$. The functions $e_n(g;x)$ are real-valued and belong to ${D(g_0,g_1)}$. They are of the form , where $d_n(g;x)$ extends to a function that is analytic in $|\Im z|<\alpha/2$, even and $\pi/r$-periodic. Moreover, we have $$\begin{gathered}
\label{dnnz}
d_n(g;0)\ne 0,\qquad d_n(g;\pi/2r)\ne 0,\end{gathered}$$ and the spectrum of $H(g)$ is nondegenerate.

We have just proved the first assertion. The second and third one have been shown earlier too, cf. the paragraphs containing (\[ed\]) and (\[cD\]). To prove the last one, we first note that it follows from $e_n(g)\in{D(g_0,g_1)}$ that $e_n(g;x)$ is a classical solution to the ODE $H(g)e_n=\lambda_ne_n$ for a certain $\lambda_n\in{{\mathbb R}}$. From the Frobenius method it then follows that $d_n(g;x)$ satisfies (\[dnnz\]). Moreover, any solution to the same ODE that is linearly independent of $e_n(g;x)$ does not have the features encoded in (\[ed\]), so the eigenvalue $\lambda_n$ is simple.

Quite likely, for $g\in\tilde{\Pi}\setminus \Pi_r$ the $H(g)$-eigenvectors still belong to ${D(g_0,g_1)}$ and have the same features as for $g\in\Pi_r$. (If so, it follows as before that $\sigma(H(g))$ is simple.) However, we cannot prove this via the above arguments, since they hinge on exploiting the HS operator ${{\cal I}}(g)$ with $g$ restricted to $\Pi_r$.

It is a corollary of Theorem \[theorem2.4\] that there exist distinct real numbers $\lambda_n(g)$ such that $$\begin{gathered}
\label{Hlam}
H(g)e_n(g)=\lambda_n(g)e_n(g), \qquad n\in{{\mathbb N}},\qquad g\in\Pi_r.\end{gathered}$$ At this stage, however, we have no other information on the $H(g)$-eigenvalues. In the following theorem, we obtain further results on the $H(g)$-spectrum that do not involve ${{\cal I}}(g)$, and which are valid for all $g\in\tilde{\Pi}$. As a preparation, note that $H(g)$ is obviously bounded below for $g_0,g_1\ge 1$. (Indeed, this entails that $V(g;x)$ is bounded below, cf. (\[Vg\]).) For $g_0$ and/or $g_1$ in $(0,1)$, however, $V(g;x)$ diverges to $-\infty$ as $x\downarrow 0$ and/or $x\uparrow \pi/2r$. Even so, our next theorem implies that $H(g)$ is still bounded below. Its proof is based on a comparison argument involving the trigonometric $BC_1$ Calogero–Moser Hamiltonian (\[Ht\]).

\[theorem2.5\] For all $g\in\tilde{\Pi}$, the spectrum of the operator $H(g)$ is a discrete set that is bounded below. Moreover, each eigenvalue has finite multiplicity.

The trigonometric and elliptic differential operators (\[Ht\]) and (\[Hg\]) differ by a real-valued potential $$\begin{gathered}
V_d(g;x)\equiv V_t(g_0,g_1;x)-V(g;x),\end{gathered}$$ that is bounded on $[0,\pi/2r]$. Hence we can associate to (\[Ht\]) the operator on ${D(g_0,g_1)}$ given by $$\begin{gathered}
\label{Vd}
H_t(g_0,g_1)\equiv H(g)+V_d(g;\cdot).\end{gathered}$$ This operator is essentially self-adjoint on ${D(g_0,g_1)}$, since $H(g)$ is e.s.a. on ${D(g_0,g_1)}$ by virtue of Lemma \[lemma2.3\] and the multiplication operator $V_d(g;\cdot)$ is bounded and self-adjoint. We now introduce the trigonometric weight function $$\begin{gathered}
w_t(g_0,g_1;x)\equiv (\sin rx)^{2g_0}(\cos rx)^{2g_1},\end{gathered}$$ and observe that we may rewrite ${D(g_0,g_1)}$ as $$\begin{gathered}
{D(g_0,g_1)}=w_t(g_0,g_1;x)^{1/2}{{\cal D}}.\end{gathered}$$ (Indeed, ${{\cal D}}$ (\[cD\]) is closed under products and we have $$\begin{gathered}
\left(\frac{s(x)}{\sin rx}\right)^{\lambda},\ \ \left(\frac{s(\pi/2r-x)}{\cos rx}\right)^{\mu}\in{{\cal D}},\qquad \lambda,\mu\in{{\mathbb R}},\end{gathered}$$ as is easily checked.)

Next, we calculate the similarity-transformed differential operator $$\begin{gathered}
A_t(g_0,g_1;x)     \equiv     w_t(g_0,g_1;x)^{-1/2}H_t(g_0,g_1;x)w_t(g_0,g_1;x)^{1/2}
\nonumber \\
 \phantom{A_t(g_0,g_1;x)}{}  =     -\frac{d^2}{dx^2}-2r(g_0\cot rx -g_1\tan rx)\frac{d}{dx}+r^2(g_0+g_1)^2.\end{gathered}$$ This yields an essentially self-adjoint operator $A_t(g_0,g_1)$ on the dense subspace ${{\cal D}}$ (given by (\[cD\])) of the Hilbert space $$\begin{gathered}
{{\cal H}}_t\equiv L^2([0,\pi/2r],w_t(g_0,g_1;x)dx).\end{gathered}$$ The crux is now that $A_t(g_0,g_1)$ has an orthogonal base of eigenvectors belonging to ${{\cal D}}$ with nonnegative eigenvalues, explicitly given by $$\begin{gathered}
\label{dmt}
d_{m,t}(g_0,g_1;x)= {}_2F_1(-m,g_0+g_1+m;g_0+1/2;\sin^2rx),\qquad m\in{{\mathbb N}},
\\
\label{teig}
A_t(g_0,g_1)d_{m,t}(g_0,g_1)=r^2(g_0+g_1+2m)^2d_{m,t}(g_0,g_1),\qquad m\in{{\mathbb N}}.\end{gathered}$$ Indeed, the r.h.s. of (\[dmt\]) yields polynomials of degree $m$ in $\sin^2rx$, and these polynomials amount to the Jacobi polynomials $$\begin{gathered}
P_m^{(\alpha,\beta)}(y),\qquad \alpha\equiv g_0-1/2,\qquad \beta\equiv g_1-1/2,\qquad y\equiv \cos 2rx,\qquad  m\in{{\mathbb N}},\end{gathered}$$ which form an orthogonal base in $$\begin{gathered}
{{\cal H}}_P\equiv L^2\big([-1,1],(1-y)^{\alpha}(1+y)^{\beta}dy\big),\end{gathered}$$ see for instance [@kosw].

From the above it follows that the spectrum of $H_t(g_0,g_1)$, $g_0,g_1>-1/2$, is bounded below, so by (\[Vd\]) the same is true for $H(g)$, $g\in\tilde{\Pi}$. Furthermore, it is clear from the eigenvalues (\[teig\]) that the resolvent of $H_t(g_0,g_1)$ is compact. By the second resolvent identity this also holds true for the resolvent of $H(g)$, so the theorem follows.

This theorem also reveals the special character of the (at most four) self-adjoint extensions of $H(g)$ restricted to ${{\cal C}}$ (\[cC\]) that arise via ${{\cal I}}(g)$: They are the extensions whose trigonometric limits yield the (at most four) self-adjoint extensions of $H_t(g_0,g_1)$ restricted to ${{\cal C}}$ that are associated with the Jacobi polynomials. Of course, this is particularly clear when $V(g;x)$ and $V_t(g_0,g_1;x)$ vanish identically, cf. Section \[section4\].

For the remainder of this section we assume $g\in\Pi_r$. It follows from Theorems \[theorem2.4\] and \[theorem2.5\] that $H(g)$ has nondegenerate eigenvalues $$\begin{gathered}
E_0(g)<E_1(g)<E_2(g)<\cdots.\end{gathered}$$ Also, since $H(g)$ and $H_t(g_0,g_1)$ differ by a bounded self-adjoint operator (cf. (\[Vd\])), it follows from (\[teig\]) that the eigenvalues have asymptotics $$\begin{gathered}
E_m(g)\sim 4r^2m^2,\qquad m\to\infty.\end{gathered}$$ By Theorem \[theorem2.4\] the associated eigenvectors $u_m(g)$ are also eigenvectors of $T(g)$, and they are just the scalar multiples of $e_n(g)$ for some $n\in{{\mathbb N}}$. Hence we may and will fix the ordering choice of the $T(g)$-eigenvectors for a given degenerate eigenvalue (if any occur), attained for example for $\nu_j=\nu_{j+1}=\cdots =\nu_{j+k}$, by requiring that the associated $H(g)$-eigenvalues satisfy $\lambda_j<\lambda_{j+1}<\cdots <\lambda_{j+k}$ . With this ordering of the $e_n(g)$ understood from now on, there is a uniquely determined permutation $\tau_g$ of the nonnegative integers such that $$\begin{gathered}
\label{perm}
\lambda_{\tau_g(m)}(g)=E_m(g),\qquad m\in{{\mathbb N}},\end{gathered}$$ cf. (\[Hlam\]).

Next, we observe that we still have not completely fixed the choice of the basis vectors $e_n(g)$. Indeed, we have required real-valuedness of $e_n(g;x)$, but this leaves a sign undetermined. We may and will fix this ambiguity by requiring $$\begin{gathered}
d_n(g;0)>0,\qquad \forall \, n\in{{\mathbb N}},\end{gathered}$$ cf. (\[dnnz\]).

At this point the alert reader might (or rather should) object that the sign choice just detailed may be in conflict with the decomposition (\[Psirep\]). Indeed, once the ONB-functions $e_n(g;x)$ have been fixed, the ONB-functions $e_n(g';x)$ are uniquely determined, and the signs of the latter might differ from the above sign convention.

It follows from the last theorem of this section that this contingency does not arise. Its main result is, however, that $H(g)$ and $H(g')$ have the same spectrum.

\[theorem2.6\] The integral operator ${{\cal I}}(g)$ is real-analytic in $g$ on $\Pi_r$ in the Hilbert–Schmidt norm. With the above sign convention for $e_n(g;x)$ understood, it has singular value decomposition . Furthermore, we have an equality of spectra $$\begin{gathered}
\sigma(H(g))=\sigma(H(g')),\qquad g\in\Pi_r.\end{gathered}$$

Analyticity of ${{\cal I}}(g)$ in $g$ with respect to the HS norm amounts to analyticity of its kernel $\Psi(g;x,y)$ in $g$ in the strong topology of the Hilbert space ${{\cal H}}_K$ (\[cHK\]). From the explicit formula (\[Psi\]) it is readily verified that the kernel is strongly analytic in $g$ in a complex neighborhood of $\Pi_r$, so the first assertion follows.

Next, we fix $g\in\Pi_r$ and note that the vectors $$\begin{gathered}
b_m(g;x)\equiv e_{\tau_g(m)}(g;x),\qquad m\in{{\mathbb N}},\end{gathered}$$ yield an ONB with $H(g)$-eigenvalues $E_m(g)$, cf. (\[perm\]). Recalling (\[key2\]) and (\[cIad\]), we now deduce $$\begin{gathered}
H(g'){{\cal I}}(g')b_m(g)   =  {{\cal I}}(g')H(g)b_m(g)
 =  E_m(g){{\cal I}}(g')b_m(g).$$ In view of (\[KerI\]), this entails that the vectors ${{\cal I}}(g')b_m(g)$ are eigenvectors of $H(g')$ with eigenvalues $E_m(g)$, and that they are total in ${{\cal H}}$. Hence, $$\begin{gathered}
E_m(g)=E_m(g'),\qquad \forall \, m\in{{\mathbb N}},\end{gathered}$$ and so the last assertion follows. Moreover, it follows that we have $$\begin{gathered}
\tau_g=\tau_{g'},\end{gathered}$$ and $$\begin{gathered}
{{\cal I}}(g')b_m(g)=\mu_m(g')b_m(g'),\qquad \mu_m(g')\in{{\mathbb R}}^{*}.\end{gathered}$$

To prove the theorem, it remains to show that all of the numbers $\mu_m(g)$ are positive. Now they do not vanish and are continuous on $\Pi_r$, since their squares yield eigenvalues of the self-adjoint trace class operator $T(g)$, which (by the first assertion) is continuous in $g$ on $\Pi_r$ in the trace norm topology. As $\Pi_r$ is a connected set, it suffices to check positivity for one $g$ in $\Pi_r$. For $g=(0,0,1,1)$ we have $g'=g$ and ${{\cal I}}(g)$ is self-adjoint. Thus the numbers $\mu_m(g)$ are eigenvalues of ${{\cal I}}(g)$. These are easily determined explicitly, and they are all positive, cf. (\[nun\]).

For all of the special $g\in\Pi_r$ where we have calculated the numbers $\mu_m(g)$ explicitly (cf. Section \[section4\]), they are distinct and decreasing in $m$: $$\begin{gathered}
\mu_0(g)>\mu_1(g)>\mu_2(g)>\cdots >0.\end{gathered}$$ This entails in particular that the permutation $\tau_g$ is the identity permutation. It is an open problem whether these features hold true on all of $\Pi_r$.

A hidden $\boldsymbol{S_4}$ symmetry {#section3}
====================================

In this section our purpose is to study, for a given coupling vector $g\in \Pi_r$, the couplings $\hat{g}$ for which the spectra of $H(g)$ and $H(\hat{g})$ are equal. To this end we note that we can establish equality of spectra by exploiting the core property of ${D(g_0,g_1)}$, cf. Lemma \[lemma2.3\]. First, since the action of the differential operator $H(g)$ on the domain ${D(g_0,g_1)}$ is invariant under taking $g_t\to 1-g_t$, and the domain does not depend on $g_2$ and $g_3$, we get equality of spectra (and even equality of operators) under the transformations $$\begin{gathered}
\label{triv}
g_2\to 1-g_2,\qquad g_3\to 1-g_3.\end{gathered}$$ (A caveat is in order: these maps can take $g$ out of $\Pi_r$.) Less trivially, consider the permutation $$\begin{gathered}
\label{refl}
(g_0,g_1,g_2,g_3)\to (g_1,g_0,g_3,g_2).\end{gathered}$$ Defining the ‘mirror’ unitary $$\begin{gathered}
\label{mir}
(Mf)(x)\equiv f(\pi/2r-x),\qquad x\in[0,\pi/2r],\qquad f\in{{\cal H}},\end{gathered}$$ we clearly get $$\begin{gathered}
M{D(g_0,g_1)}=D(g_1,g_0).\end{gathered}$$ This implies $$\begin{gathered}
MH(g_0,g_1,g_2,g_3)M=H(g_1,g_0,g_3,g_2),\end{gathered}$$ so that (\[refl\]) yields operators with equal spectra; in this case the two operators are different for generic $g$.

Next, we note that the matrix $J_N$ (\[defJN\]) can be viewed as the reflection with respect to the unit vector $(1,-1,1,-1)/2$: $$\begin{gathered}
J_N={\bf 1}_4-2u\otimes u,\qquad u\equiv (1,-1,1,-1)/2.\end{gathered}$$ Thus the self-dual couplings are given by $$\begin{gathered}
g'=g \ \Leftrightarrow \ g_0+g_2=g_1+g_3,\end{gathered}$$ and on this 3-dimensional submanifold of $\Pi_r$ the HS operator ${{\cal I}}(g)$ is self-adjoint. The maps (\[triv\]) and (\[refl\]) just considered do not involve $J_N$. But it is clear that even when we start from a self-dual $g$, the maps (\[triv\]) give rise to new couplings $\hat{g}$ that are generically not self-dual. When also $\hat{g}\in\Pi_r$, we can act with $J_N$ to obtain $\hat{g}'\in\Pi_r$ for which $H(\hat{g}')$ has the same spectrum as $H(g)$ (by virtue of Theorem \[theorem2.6\]).

All of the maps studied thus far are involutions. The question now arises which group is generated by these involutions. We proceed to answer this question.

For this purpose it is expedient to switch to different couplings. We do this in two steps. First, we set $$\begin{gathered}
\label{lg}
\lambda =g-\zeta/2,\qquad \zeta\equiv (1,1,1,1),\end{gathered}$$ so that $$\begin{gathered}
\label{Hl}
H(g)\to -\frac{d^2}{dx^2}+\sum_t\left(\lambda_t^2-\frac{1}{4}\right)\wp(x+\omega_t).\end{gathered}$$ Also, (\[triv\]) and (\[refl\]) become $$\begin{gathered}
\label{ltr}
\lambda_2\to -\lambda_2,\qquad \lambda_3\to -\lambda_3,
\\
\label{lre}
(\lambda_0,\lambda_1,\lambda_2,\lambda_3)\to (\lambda_1,\lambda_0,\lambda_3,\lambda_2),\end{gathered}$$ and we have $$\begin{gathered}
\lambda = \lambda' \ \Leftrightarrow  \ \lambda_0 +\lambda_2 =\lambda_1 + \lambda_3,\end{gathered}$$ since $\zeta$ is self-dual.

The second step may seem unmotivated at first sight. It consists in defining couplings $$\begin{gathered}
\label{cl}
c_0=\lambda_0 +\lambda_3,\qquad c_3=\lambda_0 -\lambda_3, \qquad c_1=\lambda_1 +\lambda_2,\qquad c_2=\lambda_1 -\lambda_2.\end{gathered}$$ The point of this definition is that the maps (\[ltr\]) become permutations $$\begin{gathered}
\label{ctr}
c_0\leftrightarrow c_3,\qquad c_1\leftrightarrow c_2,\end{gathered}$$ while (\[lre\]) turns into the permutation $$\begin{gathered}
\label{cre}
(c_0,c_1,c_2,c_3)\to (c_1,c_0,c_3,c_2).\end{gathered}$$ Moreover, the $J_N$-action transforms into the permutation $$\begin{gathered}
\label{Jc}
c_2\leftrightarrow c_3,\end{gathered}$$ as is easily checked.

Since the three transpositions (\[ctr\]) and (\[Jc\]) already generate $S_4$, it follows that (\[cre\]) may be viewed as a suitable product of these maps. The symmetries (\[ctr\]), (\[cre\]) generate an 8-element $S_4$-subgroup isomorphic to the dihedral group $I_2(4)$. This subgroup gives rise to a decomposition of $S_4$ into three cosets. Specifically, as coset-representants we may choose the three permutations $$\begin{gathered}
\label{3per}
(c_0,c_1,c_2,c_3)\ \to \ (c_0,c_1,c_2,c_3),\   (c_0,c_1,c_3,c_2),\  (c_0,c_2,c_3,c_1).\end{gathered}$$

Now that we have clarified the character of the symmetry group by exploiting the new couplings $c$, we may and will transform back to $g$. In particular, for a given $g$ the three maps (\[3per\]) are $$\begin{gathered}
g\to g,g',\tilde{g},\end{gathered}$$ where $$\begin{gathered}
  \tilde{g}_0\equiv \frac{1}{2}(g_0+g_1+g_2+g_3-1),\qquad \tilde{g}_1\equiv \frac{1}{2}(g_0+g_1-g_2-g_3+1),
  \nonumber\\
 \tilde{g}_2\equiv \frac{1}{2}(-g_0+g_1-g_2+g_3+1),\qquad \tilde{g}_3\equiv \frac{1}{2}(g_0-g_1-g_2+g_3+1).$$ We now summarize and extend the above analysis in the following theorem.

\[theorem3.1\] The group $G$ generated by the maps $$\begin{gathered}
\label{g1}
(g_0,g_1,g_2,g_3)   \to (g_0,g_1,1-g_2,g_3),
\\
\label{g2}
(g_0,g_1,g_2,g_3)   \to (g_0,g_1,g_2,1-g_3),\end{gathered}$$ and $$\begin{gathered}
\label{g3}
g\to J_Ng,\end{gathered}$$ with $J_N$ given by , is isomorphic to $S_4$ and contains the mirror map $$\begin{gathered}
(g_0,g_1,g_2,g_3)   \to (g_1,g_0,g_3,g_2).\end{gathered}$$ Now assume $g$ belongs to the set $\Pi_G$, cf. . Then the $G$-orbit of $g$ belongs to $\Pi_G$ and one has $$\begin{gathered}
\label{sprel}
\sigma (H(g))= \sigma(H(w(g))),\qquad \forall \, w\in G.\end{gathered}$$

We have already shown the first assertion. Next, we recall that $\Pi_r$ is defined by the five linear constraints $$\begin{gathered}
2g_0>-1,\qquad 2g_1>-1,\qquad g_0+g_1+g_2+g_3>0,\qquad g_0+g_1-g_2+g_3>-1,\\
 g_0+g_1+g_2-g_3>-1.\end{gathered}$$ Together with the extra constraint defining $\Pi_G$, i.e., $$\begin{gathered}
g_0+g_1-g_2-g_3>-2,\end{gathered}$$ this yields six constraints that transform into $$\begin{gathered}
\label{six}
c_{\mu}+c_{\nu}>-2,\qquad \mu,\nu\in \{ 0,1,2,3\},\qquad \mu\ne\nu,\end{gathered}$$ as is readily verified from (\[cl\]) and (\[lg\]). Since the latter are manifestly permutation invariant, the $G$-orbit assertion follows. Finally, since (\[sprel\]) holds for the generators (\[g1\])–(\[g3\]) and $\Pi_G$ is invariant under $G$, it holds for all $w\in G$.

Rewriting (\[Hl\]) in terms of $c$, we obtain $$\begin{gathered}
\hat{H}(c)       =       -\frac{d^2}{dx^2}-\frac{1}{4}\sum_{t=0}^3\wp(x+\omega_t)
 +\frac{1}{4}\left(c_0^2+c_3^2\right)\sum_{t=0,3}\wp(x+\omega_t)+\frac{1}{2}c_0c_3(\wp(x)-\wp(x+\omega_3))
\nonumber\\
\phantom{\hat{H}(c)  = }{}  +\frac{1}{4}\left(c_1^2+c_2^2\right)\sum_{t=1,2}\wp(x+\omega_t)+\frac{1}{2}c_1c_2(\wp(x+\omega_1)-\wp(x+\omega_2)).   \end{gathered}$$ Requiring (\[six\]), the spectral invariance under the $I_2(4)$ subgroup is readily understood in this guise as well. But the invariance under any other $c$-permutation looks bizarre.

Special cases {#section4}
=============

It is plain from (\[Hg\]), (\[Vg\]) that the Heun differential operator reduces to $-d^2/dx^2$ for the sixteen $g$’s obtained by choosing $g_t\in \{ 0,1\}$, $t=0,1,2,3$. Viewing $-d^2/dx^2$ as an operator on ${{\cal H}}$ with domain ${D(g_0,g_1)}$, it is also clear which ONB $\{ b_n\}_{n\in{{\mathbb N}}}$ of eigenvectors in ${D(g_0,g_1)}$ arises. Specifically, with the normalization constant $$\begin{gathered}
{{\cal N}}\equiv (4r/\pi)^{1/2},\end{gathered}$$ the following ONBs are involved: $$\begin{gathered}
{\rm (i)}\ \ g_0=g_1=0\ \Rightarrow\ b_0={{\cal N}}/\sqrt{2},\ \ \ b_n={{\cal N}}\cos 2nrx,\ \ \ n>0,
\\
{\rm (ii)}\ \ g_0=g_1=1\ \Rightarrow\ b_n={{\cal N}}\sin (2n+2)rx,
\\
{\rm (iii)}\ \ g_0=1,\ g_1=0\ \Rightarrow\ b_n={{\cal N}}\sin (2n+1)rx,
\\
{\rm (iv)}\ \ g_0=0,\ g_1=1\ \Rightarrow\ b_n={{\cal N}}\cos (2n+1)rx.\end{gathered}$$

We proceed to study the above $g$-choices in relation to ${{\cal I}}(g)$. For case (i) we must choose $g_2=g_3=1$ to obtain $g\in\Pi_r$. (For $g_2=g_3=0$ we get $s_g=0$, and for $g_2=0$, $g_3=1$ and $g_2=1$, $g_3=0$ we get $g_1'=-1/2$ and $g_0'=-1/2$, resp.) Then we obtain from (\[Psi\])–(\[cS\]) and (\[w2\]) $$\begin{gathered}
\Psi((0,0,1,1);x,y)=R(x)R(\omega_1-x)R(y)R(\omega_1-y)/R(x+y)R(x-y).\end{gathered}$$ It now follows from Theorem \[theorem2.4\] that the ONB (i) consists of eigenvectors for ${{\cal I}}((0,0,1,1))$. In fact, this can be readily verified by a contour integration, which also yields the eigenvalues explicitly. (To simplify the residues, one needs the duplication formula for the $R$-function, cf. (\[dupl\]).) Specifically, we obtain $$\begin{gathered}
\label{nun}
\nu_n=\pi/p\cosh (nr\alpha),\qquad g=(0,0,1,1),\qquad c=(0,0,-1,-1),\end{gathered}$$ with $p$ the infinite product in (\[Rs\]).

Turning to case (ii), we get two self-dual subcases. For the first one, $$\begin{gathered}
\Psi((1,1,0,0);x,y)=p^4 s(x)s(\omega_1-x)s(y)s(\omega_1-y)/R(x+y)R(x-y),\end{gathered}$$ the eigenvalues of the corresponding integral operator ${{\cal I}}(g)$ on the ONB (ii) can be calculated in the same way as before. The result is $$\begin{gathered}
\nu_n=\pi e^{r\alpha}/p\cosh ((n+1) r\alpha),\qquad g=(1,1,0,0),\qquad c=(0,0,1,1).\end{gathered}$$ Using (\[dupl\]) and (\[Rs\]), the second subcase can be rewritten as $$\begin{gathered}
\Psi((1,1,1,1);x,y)=p^2s(2x)s(2y)/R^2(x+y)R^2(x-y).\end{gathered}$$ With due effort (involving in particular a suitable use of (\[Rs\])), the eigenvalues for this $g$-choice can again be determined explicitly, yielding $$\begin{gathered}
\nu_n=\frac{2\pi(n+1)re^{r\alpha}}{p^2\sinh ((n+1)r\alpha)},\qquad g=(1,1,1,1),\qquad c=(1,1,0,0).\end{gathered}$$ It should be noted that the ${{\cal I}}(g)$-spectra of these two subcases do not coincide, even though they are related by a $c$-permutation.

The remaining two subcases are not self-dual, but they are related by mirror symmetry, cf. (\[mir\]). Hence we only study $$\begin{gathered}
\label{remc}
g=(1,1,0,1),\ c=(1,0,1,0) \ \Leftrightarrow \ g'=\frac{1}{2}(3,1,1,1),\ c'=(1,0,0,1).\end{gathered}$$ For this choice it is clear that $H(g')$ has a nontrivial potential. Its eigenvector ONB is related to the ONB (ii) via ${{\cal I}}(g')$. Specifically, we must have $$\begin{gathered}
b_n(g';x)= {{\cal N}}_n {\int_{0}^{\pi/2r}}dy\frac{[w(g';x)w(g;y)]^{1/2}}{[R(x+y)R(x-y)]^{3/2}}\sin (2n+2)ry,\qquad g=(1,1,0,1),\end{gathered}$$ where the normalization constant ${{\cal N}}_n$ ensures that $b_n(g';x)$ has norm 1. Also, the $H(g')$-eigenvalues are given by $$\begin{gathered}
\label{En1}
E_n((3,1,1,1)/2)=(2n+2)^2r^2.\end{gathered}$$

The cases (iii) and (iv) are related by mirror symmetry, so we only consider case (iii). First we choose $$\begin{gathered}
g=(1,0,1,1),\ c=(1,0,-1,0)\ \Leftrightarrow \ g'=\frac{1}{2}(1,1,1,3),\ c'=(1,0,0,-1).\end{gathered}$$ Thus $H(g')$ has a nontrivial potential. Just as for (\[remc\]), its eigenvector ONB is given by $$\begin{gathered}
b_n(g';x)= {{\cal N}}_n {\int_{0}^{\pi/2r}}dy\frac{[w(g';x)w(g;y)]^{1/2}}{[R(x+y)R(x-y)]^{3/2}}\sin (2n+1)ry,\qquad g=(1,0,1,1).\end{gathered}$$ The analog of (\[En1\]) is $$\begin{gathered}
E_n((1,1,1,3)/2)=(2n+1)^2r^2.\end{gathered}$$

We continue with $$\begin{gathered}
g=(1,0,0,0),\ c=(0,-1,0,1) \ \Leftrightarrow \ g'=\frac{1}{2}(1,1,-1,1),\ c'=(0,-1,1,0).\end{gathered}$$ Here we obtain the $H(g')$-ONB $$\begin{gathered}
b_n(g';x)={{\cal N}}_n  {\int_{0}^{\pi/2r}}dy\frac{[w(g';x)w(g;y)]^{1/2}}{[R(x+y)R(x-y)]^{1/2}}\sin (2n+1)ry,\qquad g=(1,0,0,0),\end{gathered}$$ with eigenvalues $$\begin{gathered}
E_n((1,1,-1,1)/2)=(2n+1)^2r^2.\end{gathered}$$ Next we consider the self-dual choice $g_2=0$, $g_3=1$. Then we get $$\begin{gathered}
\Psi((1,0,0,1);x,y)=p^2s(x)R(\omega_1-x)s(y)R(\omega_1-y)/R(x+y)R(x-y).\end{gathered}$$ As before, the eigenvalues of ${{\cal I}}(g)$ on the ONB (iii) can be explicitly determined, yielding $$\begin{gathered}
\nu_n=\pi e^{r\alpha/2}/p\cosh ((n+1/2)r\alpha),\qquad g=(1,0,0,1),\qquad c=(1,-1,0,0).\end{gathered}$$

Finally, we study the choice $$\begin{gathered}
g=(1,0,1,0),\ c=(0,0,-1,1) \ \Leftrightarrow \ g'=(0,1,0,1),\ c'=(0,0,1,-1),\end{gathered}$$ for which $$\begin{gathered}
\Psi((1,0,1,0);x,y)=p^2s(x)R(x)s(\omega_1-y)R(\omega_1-y)/R(x+y)R(x-y).\end{gathered}$$ Once more, a contour integration yields the expected result $$\begin{gathered}
\sin (2n+1)rx =
\nu_n  {\int_{0}^{\pi/2r}}dy \Psi((1,0,1,0);x,y)\cos (2n+1)ry,\end{gathered}$$ with the singular values $\nu_n$ given by $$\begin{gathered}
\nu_n = \pi e^{r\alpha/2}/p\sinh((n+1/2)r\alpha),\qquad g=(1,0,1,0),\qquad c=(0,0,-1,1).\end{gathered}$$

We conclude this section with some remarks on special cases of a different character. We recall first that ${{\cal I}}(g)$ reduces to a rank-one operator for $g\in\Pi$ with $s_g=0$, cf. (\[Psi\]) and (\[cS\]). By continuity in $g$, it therefore follows that the functions $w(g;x)^{1/2}$ and $w(g';x)^{1/2}$ are eigenfunctions of $H(g)$ and $H(g')$, resp., with the same eigenvalue. We point out that this is not at all obvious. With due effort, however, this can be directly verified.

Specifically, it follows from I (3.79), (3.113) and (3.114) that the eigenfunction property is equivalent to constancy of the function $$\begin{gathered}
\sum_{t=0}^3g_t^2\wp(x+\omega_t)-\left(\sum_{t=0}^1\left( g_t\frac{s'(x+\omega_t)}{s(x+\omega_t)}+
g_{t+2}\frac{R'(x+\omega_t)}{R(x+\omega_t)}\right)\right)^2,\qquad \sum_{t=0}^3g_t=0.\end{gathered}$$ We encountered this functional identity before on p. 253 of [@cadiz], where we supplied a (short) proof. Using next I (3.105)–(3.106), the constant (which amounts to the ground state eigenvalue) can be calculated explicitly: $$\begin{gathered}
E_0(g)=\sum_{j=1}^3g_je_j(g_j+2g_0),\qquad e_j\equiv \wp(\omega_j),\qquad j=1,2,3, \qquad \sum_{t=0}^3g_t=0.\end{gathered}$$ From this formula one readily verifies equality to $E_0(g')$.

Finally, we comment on the case where $V(g;x)$ is real-analytic on ${{\mathbb R}}$. For this we need to choose $g_0,g_1\in\{ 0,1\} $, cf. (\[Vg\]). Then one need only combine the ${{\cal H}}$-ONB’s $b_n((0,0,g_2,g_3);x)$ and $b_n((1,1,g_2,g_3);x)$ (consisting of real-analytic $\pi/r$-periodic functions that are even and odd, resp.) to obtain an ONB of $\pi/r$-periodic functions for $$\begin{gathered}
{{\cal H}}_b\equiv L^2([-\pi/2r,\pi/2r],dx).\end{gathered}$$ Likewise, an ONB of $\pi/r$-antiperiodic functions for ${{\cal H}}_b$ results from combining the ${{\cal H}}$-ONB $b_n((0,1,g_2,g_3);x)$ of even functions and the ${{\cal H}}$-ONB of odd functions $b_n((1,0,g_2,g_3);x)$. In terms of Floquet theory, therefore, one can only arrive at the ${{\cal H}}_b$-ONB’s associated with the multipliers 1 and $-1$. Whenever the four $g$-values are all in $\Pi_r$, one can invoke the spectral invariance under taking $g\to g'$. Note that $V(g';x)$ may be singular at $x=0$ and/or $x=\pi/2r$, as we have already seen in the above special cases where $V(g;x)$ vanishes identically.

Proofs of Lemmas \[lemma2.1\] and \[lemma2.2\] {#appendix-A}
==============================================

We recall that the function $$\begin{gathered}
h(x)\equiv {\int_{0}^{\pi/2r}}dy{{\cal S}}(g;x,y)\phi(y),\qquad \phi \in L^1([0,\pi/2r]),\qquad x\in [0,\pi/2r],\end{gathered}$$ extends to an even analytic function in $|{{\rm Im}\, }z|<\alpha/2$, cf. the paragraph containing (\[hdef\]). Thus its vanishing is equivalent to vanishing of the coefficients in its Taylor expansion $$\begin{gathered}
\label{hcn}
h(x)=\sum_{k=0}^{\infty}c_kx^{2k},\end{gathered}$$ at $x=0$. We proceed to study these coefficients.

Obviously, we need information on the coefficients in the expansion of ${{\cal S}}(g;x,y)$ at $x=0$. To this end we first consider the expansion $$\begin{gathered}
F(x)    \equiv   -\frac{1}{2}( \ln R(x+y)+\ln R(x-y))
    =     \sum_{j=0}^{\infty}\frac{F^{(2j)}(0)}{(2j)!}x^{2j}.
  \end{gathered}$$ Using $$\begin{gathered}
  -\frac{d^2}{dz^2}\ln R(z)={\tilde{\wp}}(z)+2\eta r/\pi,\qquad {\tilde{\wp}}(z)\equiv \wp(z+i\alpha/2),
  \end{gathered}$$ (cf. (\[Rs\]), (\[ssig\])), we obtain $$\begin{gathered}
  F(0)=-\ln R(y),\qquad F^{(2)}(0)={\tilde{\wp}}(y)+2\eta r/\pi,
\\
  F^{(2k)}(0)={\tilde{\wp}}^{(2k-2)}(y),\qquad k>1.
  \end{gathered}$$ Now from standard elliptic function lore [@ww] we deduce $$\begin{gathered}
  \wp^{(2k-2)}(z)=(2k-1)!\wp(z)^k +  {\rm l.d.},\qquad k>1,
  \end{gathered}$$ where l.d. stands for a polynomial in $\wp(z)$ of degree $<k$. Thus we have $$\begin{gathered}
  F(x)=-\ln R(y)+\sum_{j=1}^{\infty}a_j x^{2j},
  \end{gathered}$$ with $$\begin{gathered}
\label{ak}
  a_k=\frac{1}{2k}{\tilde{\wp}}(y)^k +  {\rm l.d.}
  \end{gathered}$$

We are now prepared to consider the expansion of ${{\cal S}}$. It is given by $$\begin{gathered}
  {{\cal S}}(g;x,y)     =     R(y)^{-2s_g}\exp\left( 2s_g\sum_{j=1}^{\infty}a_jx^{2j}\right)
  =    R(y)^{-2s_g}\sum_{k=0}^{\infty}b_kx^{2k},\label{Sex}\end{gathered}$$ so that the first few $b_k$ read $$\begin{gathered}
   b_0=1,\qquad b_1=2s_ga_1,\qquad b_2=2s_ga_2+2s_g^2a_1^2,\qquad b_3=2s_ga_3+4s_g^2 a_2a_1 + 4s_g^3 a_1^3/3.
   \end{gathered}$$ More generally, it is easily seen that $b_k$ is the solution to the recurrence $$\begin{gathered}
   kb_k=2s_g\sum_{j=1}^k ja_j b_{k-j},\qquad  k>0,\qquad b_0=1.
   \end{gathered}$$ In view of (\[ak\]), this implies that $b_k$ is given by $$\begin{gathered}
\label{bQ}
 b_k=Q^{(k)}(s_g){\tilde{\wp}}(y)^k+  {\rm l.d.},
 \end{gathered}$$ where $Q^{(k)}(u)$ is a polynomial of the form $$\begin{gathered}
 Q^{(k)}(u)=\sum_{l=1}^kp_l^{(k)}u^l,\qquad p^{(k)}_l\in(0,\infty).
 \end{gathered}$$

We are now in the position to use our assumption $s_g>0$. It entails that the coefficient $Q^{(k)}(s_g)$ in (\[bQ\]) is positive: $$\begin{gathered}
\label{bkp}
 b_k=q_k{\tilde{\wp}}(y)^k+ {\rm l.d.},\qquad q_k\in(0,\infty),\qquad \forall \, k\in{{\mathbb N}}.
 \end{gathered}$$ Returning to (\[hcn\]), we first note that (\[Sex\]) entails $$\begin{gathered}
 c_k={\int_{0}^{\pi/2r}}dyb_k(y)\psi(y),
 \end{gathered}$$ where $$\begin{gathered}
\label{psiR}
 \psi(y)\equiv R(y)^{-2s_g}\phi(y).
 \end{gathered}$$ Since the coefficients $c_k$ vanish by assumption, it follows recursively from (\[bkp\]) that we have $$\begin{gathered}
 {\int_{0}^{\pi/2r}}dy {\tilde{\wp}}(y)^k\psi(y)=0,\qquad \forall \, k\in{{\mathbb N}}.
 \end{gathered}$$

Now ${\tilde{\wp}}(y)$ is monotone increasing on $[0,\pi/2r]$. Thus we can invoke the Stone-Weierstrass theorem [@rs1] to conclude that the span of the functions $1,{\tilde{\wp}}(y),{\tilde{\wp}}(y)^2,{\tilde{\wp}}(y)^3,\ldots$ is dense (in the supremum norm) in the space of real-valued continuous functions on $[0,\pi/2r]$. Since $\psi(y)$ is an $L^1$-function, it follows by dominated convergence that we have $$\begin{gathered}
 {\int_{0}^{\pi/2r}}dy C(y)\psi(y)=0,\qquad \forall \, C\in C_{{{\mathbb R}}}([0,\pi/2r]).
 \end{gathered}$$ Hence $\psi=0$, so by (\[psiR\]) $\phi=0$.

Letting $e(x)\in {D(g_0,g_1)}$, we first show that the action of the differential operator $H(g)$ (\[Hg\]) on $e(x)$ yields a function in ${{\cal H}}$. To this end we use the relations $$\begin{gathered}
 e(x)=s(x)^{g_0}s(x-\omega_1)^{g_1}d(x),\qquad d\in{{\cal D}},
 \end{gathered}$$ and $$\begin{gathered}
 \frac{s'(z)^2}{s(z)^2}=\frac{s''(z)}{s(z)}+\wp(z)+\eta \omega_1,
 \end{gathered}$$ to calculate $$\begin{gathered}
e''(x)     =    e(x)\sum_{t=0,1}\left(
[g_t^2-g_t][\wp(x+\omega_t)+\eta \omega_1]
+g_t^2\frac{s''(x+\omega_t)}{s(x+\omega_t)}
\right)
\nonumber\\
\phantom{e''(x)     =}{}  +e(x)\left(
  Q(x)+2\frac{d'(x)}{d(x)}\sum_{t=0,1}g_t\frac{s'(x+\omega_t)}{s(x+\omega_t)}
   +\frac{d''(x)}{d(x)}\right),  \end{gathered}$$ where $$\begin{gathered}
\label{defQ}
Q(x)\equiv 2g_0g_1 \frac{s'(x)}{s(x)} \frac{s'(x+\omega_1)}{s(x+\omega_1)}.\end{gathered}$$ Subtracting the potential term $$\begin{gathered}
e(x)\sum_{t=0}^3 g_t(g_t-1)\wp(x+\omega_t),\end{gathered}$$ we obtain a sum of terms that are manifestly in ${{\cal H}}$, save for the terms $$\begin{gathered}
\label{rem}
e(x)\left( Q(x)+\sum_{t=0,1} g_t^2\frac{s''(x+\omega_t)}{s(x+\omega_t)}\right)
+2s(x)^{g_0}s(x-\omega_1)^{g_1}
\sum_{t=0,1}g_td'(x)\frac{s'(x+\omega_t)}{s(x+\omega_t)}.\end{gathered}$$ Now $s(z)$ is odd and has a first order zero for $z=0$, so $s''(z)$ vanishes at $z=0$ at least to first order. Thus the two terms in the first sum are smooth for $x\in{{\mathbb R}}$. Next, $s'(z)$ vanishes at $z=\omega_1$, since $s(z)$ is odd and $2\omega_1$-antiperiodic. Therefore $Q(x)$ (\[defQ\]) is also smooth. Finally, since $d'(0)=d'(\omega_1)=0$, the two terms in the second sum are continuous on $[0,\omega_1]$. Hence the function given by (\[rem\]) is in ${{\cal H}}$, so that $H(g)$ is well defined on ${D(g_0,g_1)}$.

It remains to show that $H(g)$ is symmetric on ${D(g_0,g_1)}$. To this end, let $d_1,d_2\in{{\cal D}}$ and consider $$\begin{gathered}
(s(\cdot)^{g_0}s(\cdot +\omega_1)^{g_1}d_1,H(g)s(\cdot)^{g_0}s(\cdot +\omega_1)^{g_1}d_2)\nonumber\\
\qquad{}-
(H(g)s(\cdot)^{g_0}s(\cdot +\omega_1)^{g_1}d_1,s(\cdot)^{g_0}s(\cdot +\omega_1)^{g_1}d_2).\label{ipd}\end{gathered}$$ From the above calculation we see that this equals $$\begin{gathered}
  \int_0^{\omega_1} dx s(x)^{2g_0}s(x+\omega_1)^{2g_1}\biggl( \overline{d_1^{''}(x)}d_2(x)-\overline{d_1(x)}
d_2^{''}(x)
   \nonumber\\
\qquad{} +2[\overline{d_1^{'}(x)}d_2(x)-
\overline{d_1(x)}d_2^{'}(x)]\sum_{t=0,1}g_t\frac{s'(x+\omega_t)}{s(x+\omega_t)}\biggr).\end{gathered}$$ Clearly, this can be rewritten as $$\begin{gathered}
\label{dd}
  \int_0^{\omega_1}\! dx \left(\!-\overline{d_1(x)}\frac{d}{dx}[s(x)^{2g_0}s(x+\omega_1)^{2g_1}d_2'(x)]
  +d_2(x)\frac{d}{dx}[s(x)^{2g_0}s(x+\omega_1)^{2g_1}\overline{d_1'(x)}]\!\right)\!.\!\!\!
  \end{gathered}$$ Now for all $d\in{{\cal D}}$ we have $$\begin{gathered}
d'(x)=\int_0^x d''(t)dt =O(x),\quad x\downarrow 0,\qquad d'(x)=-\int_x^{\omega_1}d''(t)dt=O(x-\omega_1),\quad x\uparrow \omega_1.\end{gathered}$$ Integrating by parts in (\[dd\]), it follows that the boundary terms vanish, whereas the new integrand vanishes identically. As a result, the difference (\[ipd\]) vanishes, so that $H(g)$ is indeed symmetric on ${D(g_0,g_1)}$.

Acknowledgments {#acknowledgments .unnumbered}
---------------

We would like to thank F. Nijhoff, B. Sleeman and K. Takemura for illuminating discussions and for supplying information about related literature.

[99]{}

Ruijsenaars S.N.M., Hilbert–Schmidt operators vs. integrable systems of elliptic Calogero–Moser type I. The eigenfunction identities, [*Comm. Math. Phys.*]{} [**286**]{} (2009), 629–657.

Ruijsenaars S.N.M., Hilbert–Schmidt operators vs. integrable systems of elliptic Calogero–Moser type II. The $A_{N-1}$ case: First steps, [*Comm. Math. Phys.*]{} [**286**]{} (2009), 659–680.

Ruijsenaars S.N.M., Elliptic integrable systems of Calogero–Moser type: Some new results on joint eigenfunctions, in Proceedings of the 2004 Kyoto Workshop on Elliptic integrable systems, Editors M. Noumi and K. Takasaki, [*Rokko Lectures in Math.*]{}, Dept. of Math., Kobe Univ., 2004, 223–240.

Whittaker E.T., Watson G.N., A course of modern analysis, Cambridge University Press, Cambridge, 1973.

Inozemtsev V.I., Lax representation with spectral parameter on a torus for integrable particle systems, [*Lett. Math. Phys.*]{} [**17**]{} (1989), 11–17.

Ronveaux A. (Editor), Heun’s differential equations, Oxford University Press, Oxford, 1995.

Sleeman B.D., Kuznetsov V.B., Heun functions, Chapter HE in Handbook of mathematical functions, Editor F. Olver, to appear.

Maier R., The 192 solutions of the Heun equation, [*Math. Comp.*]{} [**76**]{} (2007), 811–843, [math.CA/0408317](http://arxiv.org/abs/math.CA/0408317).

Takemura K., Middle convolution and Heun’s equation, [*SIGMA*]{} [**5**]{} (2009), 040, 22 pages, [arXiv:0810.3112](http://arxiv.org/abs/0810.3112).

Reed M., Simon B., Methods of modern mathematical physics. II. Fourier analysis, self-adjointness, Academic Press, New York, 1975.

Slavyanov S.Y., Lay W., Special functions. A unified theory based on singularities, Oxford University Press, Oxford, 2000.

Kazakov A.Y., Slavyanov S.Y., Integral relations for Heun-class special functions, [*Theoret. and Math. Phys.*]{} [**107**]{} (1996), 733–739.

Takemura K., Integral representation of solutions to Fuchsian system and Heun’s equation, [*J. Math. Anal. Appl.*]{} [**342**]{} (2008), 52–69, [arXiv:0705.3358](http://arxiv.org/abs/0705.3358).

Ruijsenaars S.N.M., First order analytic difference equations and integrable quantum systems, [*J. Math. Phys.*]{} [**38**]{} (1997), 1069–1146.

Reed M., Simon B., Methods of modern mathematical physics. I. Functional analysis, Academic Press, New York, 1972.

Ruijsenaars S.N.M., Hilbert–Schmidt operators vs. integrable systems of elliptic Calogero–Moser type IV. The relativistic Heun (van Diejen) case, in preparation.

Takemura K., The Heun equation and the Calogero–Moser–Sutherland system. II. Perturbation and algebraic solution, [*Electron. J. Differential Equations*]{} [**2004**]{} (2004), no. 15, 30 pages, [math.CA/0112179](http://arxiv.org/abs/math.CA/0112179).

Ruijsenaars S.N.M., Integrable $BC_N$ analytic difference operators: Hidden parameter symmetries and eigenfunctions, in Proceedings Cadiz 2002 NATO Advanced Research Workshop “New trends in integrability and partial solvability”, [*NATO Science Series*]{}, Vol. 132, Editors A.B. Shabat, A. González-López, M. Mañas, L. Martínez Alonso and M.A. Rodríguez, Kluwer, Dordrecht, 2004, 217–261.

Arscott F.M., Periodic differential equations, Pergamon Press, London, 1964.

Koekoek R., Swarttouw R.W., The Askey-scheme of hypergeometric orthogonal polynomials and its $q$-analogue, Delft University of Technology, Faculty of Information Technology and Systems, Department of Technical Mathematics and Informatics, Report no. 98-17, 1998, [math.CA/9602214](http://arxiv.org/abs/math.CA/9602214).
