text:
# 1 Introduction

## 1 Introduction

Recently, Berenstein, Maldacena and Nastase (BMN)  have discovered a particular limit of the AdS/CFT correspondence – known as the Penrose limit – in which string theory can be solved to all orders in $`\alpha ^{}`$ . The limiting geometry, which is the plane wave in , has very different global properties compared to those of $`AdS_5\times S^5`$. Since the Penrose limit focuses on the geometry in the vicinity of a particular massless particle in $`AdS_5\times S^5`$, which is very massive when measured by an $`AdS_5`$ observer, information about the asymptotic $`AdS_5`$ time-like boundary is lost in the limit, and the conformal boundary of the plane wave degenerates to a null line . This peculiar degeneration makes it challenging to find a precise holographic map between string theory and gauge theory. On the other hand, the limit performed on the geometry can also be realized on the dual gauge theory which leads BMN to conjecture that a certain sector of $`𝒩=4`$ SYM captures the dynamics of string theory in the plane wave background.

A crucial step in making the correspondence precise is to identify the charges of string states in space-time with the charges carried by gauge theory operators. The identification is given by

$$\frac{1}{\mu }H=(\mathrm{\Delta }J),\mu R^2P^+=J,$$
(1.1)
where $`H`$ is the generator of $`x^+`$ translations, $`P^+`$ is the generator of $`x^{}`$ translations, $`\mathrm{\Delta }`$ is the generator of dilatations and $`J`$ is a $`U(1)_R`$ charge. In the limit, the gauge theory truncates to the BMN sector of operators carrying large $`U(1)_R`$ charge $`J\sqrt{N}`$ in the $`N\mathrm{}`$ limit. In this double scaling limit we keep the following parameters fixed

$$\lambda ^{}=\frac{g^2N}{J^2}=\frac{1}{(\mu p^+\alpha ^{})^2},g_2=\frac{J^2}{N}=4\pi g_s(\mu p^+\alpha ^{})^2,$$
(1.2)
which effectively count quantum loops  and non-planar corrections , respectively.

Even in the absence of string interactions, the BMN proposal makes a remarkable prediction about the spectrum of anomalous dimensions of BMN operators in the gauge theory in terms of the free string spectrum on the plane wave. This prediction was analyzed perturbatively in the regime where $`\lambda ^{}1`$ in  and confirmed to all orders using superconformal invariance in .

In the light-cone gauge, in which the plane wave string theory can be readily analyzed, the free Hamiltonian receives string corrections

$$H=H_2+g_2H_3+\mathrm{},$$
(1.3)
which to lowest order result in non-trivial transition amplitudes between one and two-string states. In  a concrete proposal was made which relates these string theory transition amplitudes with three point functions of BMN operators in the gauge theory. Their proposal is therefore the first attempt to construct a holographic map at the interacting level between string theory in a plane wave and gauge theory.

This interesting proposal111Various aspects of this proposal were considered in . was put to an explicit test by Spradlin and Volovich in  where some Hamiltonian matrix elements were computed using the string field theory vertex constructed previously in . Exact agreement with the proposal given in  was reported. Unfortunately both the field theory and string theory computations suffer from errors which when taken into account invalidate the proposal. On the field theory side, operator mixing is more important than initially contemplated in . Three-point functions of single trace operators with order $`\lambda ^{}`$ interactions are not conformally invariant, but the correct form is restored after operator mixing is incorporated, as shown in . On the string field theory side, the prefactor acting on the delta functional overlap of three strings that enters the Hamiltonian $`H_3`$  has a minus sign error, which was first reported by Pankiewicz . Once matrix elements are recomputed with the corrected Hamiltonian $`H_3`$, which we calculate in section $`2`$, agreement with field theory is lost. Therefore, what is the correct holographic map between string theory in the plane wave  and $`𝒩=4`$ SYM at the interacting level?

The most straightforward way to proceed, which was first advocated in a paper by Gross, Mikhailov and Roiban , is to take the identification (1.1) between the string field theory Hamiltonian $`H`$ and the generator of scale transformations $`\mathrm{\Delta }J`$ in $`𝒩=4`$ SYM as the holographic map for all $`g_2`$. This holographic map therefore identifies Hamiltonian matrix elements in string field theory with those of the dilatation operator in gauge theory. In order to test this identification one must find an explicit map between the Hilbert space of string states and the Hilbert space of states in the gauge theory. Moreover, in order for the comparison to be meaningful one must compute the matrix elements of these operators in a basis in which the Hilbert space inner product in gauge theory is the same as that in string field theory. The obvious inner product in the Hilbert space of string theory is the familiar inner product where, for example, the one-string states are orthogonal to two-string states.

In gauge theory, the Hilbert space inner product is induced by the matrix of two-point functions of BMN operators222In this formula the various rows and columns describe single trace, double trace, etc components. A more complete characterization of this matrix is given in section $`3`$.

$$|x|^{2\mathrm{\Delta }_0}O_A\overline{O}_B=G_{AB}+\mathrm{\Gamma }_{AB}\mathrm{ln}(x^2\mathrm{\Lambda }^2)^1,$$
(1.4)
where $`G_{AB}`$ is the Hilbert space inner product and $`\mathrm{\Gamma }_{AB}`$ is the matrix of anomalous dimensions. Unlike with the usual Hilbert space inner product in string field theory, which remains diagonal to all orders in $`g_2`$, perturbative corrections in gauge theory induce operator mixing at each order in $`g_2`$ in perturbation theory and the Hilbert space inner product is no longer diagonal. Direct comparison with string field theory calculations requires correcting for operator mixing systematically, order by order in the $`g_2`$, expansion by making $`G_{AB}`$ orthonormal via a change of basis.

In order to relate string theory to gauge theory calculations via (1.1) we must calculate the matrix elements of the dilatation operator. It is straightforward to show that the matrix elements of the dilatation operator between states created by BMN operators are given by the matrix of anomalous dimensions333Since BMN operators are BPS or nearly BPS the contribution from the action of $`\mathrm{\Delta }`$ coming from the bare dimension of the operator is cancelled by the contribution from the R-charge in (1.5).

$$O_A|(\mathrm{\Delta }J)|O_B=\mathrm{\Gamma }_{AB}.$$
(1.5)
Comparison of matrix elements of $`H`$ and $`\mathrm{\Delta }`$ requires first making the gauge theory inner product orthonormal order by order in perturbation theory. We can accomplish this by finding a new basis of operators $`\stackrel{~}{O}_A=U_{AB}O_B`$ such that they are orthonormal

$$UGU^{}=1.$$
(1.6)
When $`g_2=0`$ the correct identification between string states and gauge theory operators was given by BMN. Namely, an $`n`$-string state is described by an $`n`$-trace operator. Once $`g_2`$ corrections are taken into account this identification has to be modified. Therefore, the precise mapping between string field theory states $`|s_A`$ and gauge theory states $`|\stackrel{~}{O}_A`$ when $`g_20`$ is given by444In light-cone string field theory the canonical normalization of states is the usual delta function normalization $`s_A^{}|s_B^{}=p_A^+\delta (p_A^+p_B^+)=J_A\delta _{J_A,J_B}`$, so that $`|s_A^{}=\sqrt{J_A}|s_A`$. Therefore, when comparing string field theory results with gauge theory results we will have to take into account this normalization factor, since gauge theory states have unit norm .

$$|s_A|\stackrel{~}{O}_A=U_{AB}|O_B,s_A|s_B=\stackrel{~}{O}_A|\stackrel{~}{O}_B=\delta _{AB}.$$
(1.7)
In section $`3`$ we give an expression for the change of basis to order $`g_2^2`$. Once the right basis is found, the holographic map reads

$$\frac{1}{\mu }s_A|H|s_B=\stackrel{~}{O}_A|(\mathrm{\Delta }J)|\stackrel{~}{O}_B=(U\mathrm{\Gamma }U^{})_{AB}.$$
(1.8)
The authors of , where this proposal was first made, concluded that this proposal is equivalent to the one given in  involving three point functions. However, we find that the two holographic maps are not the same and that the one given in (1.8) is the correct one.

An important subtlety in the holographic map (1.8) is that the orthonormalization procedure of the gauge theory inner product is not unique, namely the transformation matrix $`U`$ that makes $`G_{AB}`$ orthonormal is not unique555This ambiguity was first pointed out in .. We uniquely fix the form $`U`$ to order $`g_2`$ by demanding that the dilatation operator matrix elements agree with the matrix elements of the corrected string field theory Hamiltonian that we calculate in section $`2`$. We then evaluate $`U`$ to order $`g_2^2`$ which via (1.1) makes a non-trivial prediction for string field theory Hamiltonian matrix elements to order $`g_2^2`$ which have not yet been evaluated. Remarkably, this purely gauge theory result we present is reproduced by the order $`g_2^2`$ string bit Hamiltonian  matrix elements recently computed in .

The plan of the rest of the paper is as follows. In section $`2`$ we revisit the string field theory Hamiltonian $`H_3`$ and point out that there is an incorrect relative minus sign in 666This sign has been corrected in a recent revision.. We recompute Hamiltonian matrix elements with the corrected sign, which we will use in section $`3`$. In section $`3`$ we fix the form of the change of basis in gauge theory to order $`g_2`$ by demanding that (1.8) holds when the string field theory matrix elements in section $`2`$ are used. We also evaluate the dilatation operator matrix elements to order $`g_2^2`$ by making a particular choice of basis to order $`g_2^2`$ which allows us to make a calculable prediction by using (1.8) about the string field theory Hamiltonian matrix elements to that order. We note that the final answer agrees with the recent result in  found using the string bit formalism . We conclude in section $`4`$.

### Note Added

While this work was being finished, the paper  appeared on the archive. They also reexamine the duality proposal made in  when the corrected string field theory Hamiltonian is used and also conclude that once this is taken into account that the proposal in  no longer holds. Moreover, they also adopt the philosophy in  for comparing gauge theory and string theory and also get agreement.

## 2 SFT computation revisited

Before going into our gauge theory computation, let us perform the correct string field theory calculation that we will compare it to. This also corrects some errors in the previous literature. In string field theory, the Hilbert space is a direct product of $`\mathrm{}`$-string states,

$$=_{\mathrm{}}_{\mathrm{}}.$$
(2.9)
All these states are orthogonal with respect to each other and they have an orthonormal inner product. The full Hamiltonian $`H`$, representing infinitesimal evolution along $`x^+`$, includes the freely propagating part, $`H_2`$, and an interaction part, $`H_3`$ and so on,

$$H=H_2+g_2H_3+\mathrm{}.$$
(2.10)
In the plane wave background the freely propagating part $`H_2`$ is simply the energy of an infinite collection of harmonic oscillators $`\alpha _n^i`$ and the three-string interaction part $`H_3`$ contains the delta functional overlap of three strings with a prefactor acting on it which is necessary to appropriately realize the supersymmetry algebra. The string states, which are dual to a certain class of two impurity BMN operators when $`g_2=0`$, are given by

$`|n`$ $`=\alpha _n^1\alpha _n^2|0,1`$ (2.11)
$`|n,y`$ $`=\alpha _n^1\alpha _n^2|0,y|0,1y`$ (2.12)
$`|y`$ $`=\alpha _0^1|0,y\alpha _0^2|0,1y,`$ (2.13)
where the first state represents a single string and the other two represent two-string states. Here $`|0,y`$ is the one-string vacuum state carrying a fraction $`0<y<1`$ of the total longitudinal momentum $`p^+`$ of the multi-string state. The Hamiltonian matrix elements of these states are given by

$`n|H_3|m,y`$ $`\left(F_{(1)|m|}^+F_{(3)|n|}^++F_{(1)|m|}^{}F_{(3)|n|}^{}\right)\left(\overline{N}_{|m|,|n|}^{(13)}\overline{N}_{|m|,|n|}^{(13)}\right),\text{if }mn>0`$ (2.14)
$`n|H_3|m,y`$ $`\left(F_{(1)|m|}^+F_{(3)|n|}^+F_{(1)|m|}^{}F_{(3)|n|}^{}\right)\left(\overline{N}_{|m|,|n|}^{(13)}+\overline{N}_{|m|,|n|}^{(13)}\right),\text{if }mn<0`$ (2.15)
$`n|H_3|y`$ $`F_{(3)|n|}^+\left(F_{(1)0}^+\overline{N}_{0,|n|}^{(23)}+F_{(2)0}^+\overline{N}_{0,|n|}^{(13)}\right),\text{for }n0`$ (2.16)
where $`F_{m(r)}`$ comes from the prefactor and the Neumann matrices $`\overline{N}_{mn}^{(rs)}`$ come from the delta functional overlap. The negative modes are related to the positive modes by777We would like to thank A. Pankiewicz for informing us of the correct factor of $`i`$ in the second formula. (here $`m,n>0`$)

$`\overline{N}_{m,n}^{(rs)}`$ $`=\left(U_{(r)}\overline{N}^{(rs)}U_{(s)}\right)_{m,n},`$ (2.17)
$`F_{(r)m}^{}`$ $`=i\left(U_{(r)}F_{(r)}^+\right)_m,`$ (2.18)
with $`\left(U_{(r)}\right)_{m,n}=\delta _{m,n}\left(\sqrt{m^2+(\mu p_{(r)}^+\alpha ^{})^2}\mu p_{(r)}^+\alpha ^{}\right)/m`$. Note that in (2.18) we have an extra factor of $`i`$ compared to the original literature 888More explicitly, formulas (3.15) and (3.21) in the original version of  should have an extra factor of $`i`$. This invalidates the evaluation of (4.4) ., which was first pointed out in . Using (2.17) and (2.18), we can show that both (2.14) and (2.15) reduce to the same expression.

In order to compare these answers with perturbative gauge theory we must analyze the $`\mu \mathrm{}`$ limit of (2.14), (2.15) and (2.16). The $`\mu \mathrm{}`$ behavior of both the prefactor and the Neumann matrices was evaluated in . In the $`\mu \mathrm{}`$ limit, (2.14), (2.15) and (2.16) yield999The precise overall numerical factor of the cubic string field theory Hamiltonian is not known. It is fixed by comparing with the gauge theory calculation in the next section.

$`{\displaystyle \frac{1}{\mu }}n|H_3|m,y`$ $`={\displaystyle \frac{\lambda ^{}}{2\pi ^2}}(1y)\mathrm{sin}^2(n\pi y),`$ (2.19)
$`{\displaystyle \frac{1}{\mu }}n|H_3|y`$ $`={\displaystyle \frac{\lambda ^{}}{2\pi ^2}}\sqrt{y(1y)}\mathrm{sin}^2(n\pi y).`$ (2.20)
These corrected results invalidate the agreement previously found in the literature101010Also due to this extra $`i`$, the result proven in  that the prefactor reduces to energy difference between incoming and outgoing states should be modified. Instead it reduces to the energy difference in $`\mathrm{cos}`$ modes minus that in $`\mathrm{sin}`$ modes that appear in the worldsheet Fourier decomposition.. In the next section we will use these results to test the holographic proposal in (1.8) and will show that agreement is found for a particular choice of basis.

## 3 Orthonormalization and comparison with SFT

In this section, we make a change of operator basis such that the gauge theory inner product is orthonormal order by order in $`g_2`$. We then compute the matrix elements of the operator $`\mathrm{\Delta }J`$ in this basis. An important subtlety in this procedure is that the basis change that makes the Hilbert space inner product $`G`$ orthonormal is not unique. The leading term at $`g_2=0`$ is uniquely fixed by the original correspondence explained in  between single/double string states and single/double trace operators. When $`g_20`$ the field theory operators start mixing and the matrix with which we make the inner product orthonormal is not unique, due to the familiar ambiguity when diagonalizing a matrix. We propose to fix this ambiguity in the change of basis by taking seriously the proposal in (1.8) and demanding that exact agreement with the corrected string field theory results computed in the previous section is obtained. Exact agreement is obtained for a unique choice of basis. This particular choice is the unique choice for which the transformation matrix is symmetric and real to order $`g_2`$. We then proceed to calculate to next to leading order, to order $`g_2^2`$. Here there is also an ambiguity in the orthonormalization procedure. If we make the strong assumption that the transformation matrix is still symmetric and real to this order we can uniquely determine the transformation matrix to order111111This choice is compatible with the proposal in . $`g_2^2`$. Given these assumptions we can calculate explicitly the order $`g_2^2`$ matrix element between the orthonormal basis states which reduce when $`g_2=0`$ to single trace operators. This result gives a prediction using the proposal (1.8) for light-cone string field theory matrix elements to order $`g_2^2`$ involving single strings, which has not been constructed. Despite this incompleteness of light-cone string field theory, the matrix elements we predict from the gauge theory computation exactly match the order $`g_2^2`$ matrix elements computed recently in  using the string bit formalism . It would be very desirable to understand more precisely the relation between light-cone string field theory and the string bit formalism.

The particular set of gauge theory operators that we are interested in are the following single trace and double trace operators121212Here we are using the notation in .

$`𝒪^J`$ $`={\displaystyle \frac{1}{\sqrt{JN^J}}}\mathrm{Tr}Z^J,`$ (3.21)
$`𝒪_{(i)}^J`$ $`={\displaystyle \frac{1}{\sqrt{N^{J+1}}}}\mathrm{Tr}\left(\varphi _iZ^J\right),`$ (3.22)
$`𝒪_n^J`$ $`={\displaystyle \frac{1}{\sqrt{JN^{J+2}}}}{\displaystyle \underset{l=0}{\overset{J}{}}}e^{2\pi iln/J}\mathrm{Tr}\left(\varphi _1Z^l\varphi _2Z^{Jl}\right),`$ (3.23)
$`𝒯_p^{J,y}`$ $`=:𝒪_p^{yJ}𝒪^{(1y)J}:,`$ (3.24)
$`𝒯^{J,y}`$ $`=:𝒪_{(1)}^{yJ}𝒪_{(2)}^{(1y)J}:.`$ (3.25)
We need to make the inner product $`G`$ appearing in the matrix of two-point functions in (1.4) orthonormal and eventually compute the matrix elements of $`\mathrm{\Gamma }`$ in the orthonormal basis (1.8). Both matrices $`G_{AB}`$ and $`\mathrm{\Gamma }_{AB}`$ have a systematic expansion in powers of $`g_2`$

$`G`$ $`=\mathrm{𝟏}+g_2G^{(1)}+g_2^2G^{(2)}+𝒪(g_2^3),`$ (3.26)
$`\mathrm{\Gamma }`$ $`=\mathrm{\Gamma }^{(0)}+g_2\mathrm{\Gamma }^{(1)}+g_2^2\mathrm{\Gamma }^{(2)}+𝒪(g_2^3).`$ (3.27)
In the following we split these matrices into $`3\times 3`$ blocks representing matrix elements involving $`𝒪_n^J`$, $`𝒯_p^{J,y}`$ and $`𝒯^{J,y}`$ respectively. The indices $`(n,m,\mathrm{})`$ denote the worldsheet momentum of the single trace BMN operators like, for example, $`𝒪_n^J`$. The double indices $`(py,qz,\mathrm{})`$ represent for example the worldsheet momentum and light-cone momentum fraction of the double trace operators $`𝒯_p^{J,y}`$, while $`(y,z,\mathrm{})`$ represent the fraction of momentum carried by the operator $`𝒯^{J,y}`$. These matrices have been computed in . They are given by131313We would like to thank U. Gursoy for reminding us that the off-diagonal elements of the matrix of two-point functions involving double-trace operators are non-zero at order $`g_2^2`$.:

$$G=\mathrm{𝟏}+g_2\left(\begin{array}{ccc}0& C_{n,qz}& C_{n,z}\\ C_{py,m}& 0& 0\\ C_{y,m}& 0& 0\end{array}\right)+g_2^2\left(\begin{array}{ccc}M_{n,m}^1& 0& 0\\ 0& \mathrm{?}& \mathrm{?}\\ 0& \mathrm{?}& \mathrm{?}\end{array}\right),$$
(3.28)
$$\frac{\mathrm{\Gamma }}{\lambda ^{}}=\left(\begin{array}{ccc}n^2\delta _{n,m}& 0& 0\\ 0& \frac{p^2}{y^2}\delta _{p,q}\delta _{y,z}& 0\\ 0& 0& 0\end{array}\right)+g_2\left(\begin{array}{ccc}0& \mathrm{\Gamma }_{n,qz}^{(1)}& \mathrm{\Gamma }_{n,z}^{(1)}\\ \mathrm{\Gamma }_{py,m}^{(1)}& 0& 0\\ \mathrm{\Gamma }_{y,m}^{(1)}& 0& 0\end{array}\right)+g_2^2\left(\begin{array}{ccc}nmM_{n,m}^1+\frac{1}{8\pi ^2}𝒟_{n,m}^1& 0& 0\\ 0& \mathrm{?}& \mathrm{?}\\ 0& \mathrm{?}& \mathrm{?}\end{array}\right).$$
(3.29)
The explicit form of the matrix elements are summarized in the Appendix and we denote by $`\mathrm{?}`$ matrix elements that have not yet been computed. Luckily, we will not need them for our computations. Finding them is, however, an important enterprise since via the holographic map (1.8) they predict yet unknown matrix elements in string field theory, like the order $`g_2^2`$ Hamiltonian matrix element of a two-string state141414One should also compute, however, the mixing between double and triple trace operators to get this result..

Let us now apply a linear transformation $`U`$ to make the inner product orthonormal. We require that

$$UGU^{}=\mathrm{𝟏}$$
(3.30)
and solve this equation order by order. We can expand $`U`$ in a power series in $`g_2`$ and express it as

$$U=\mathrm{𝟏}+g_2U^{(1)}+g_2^2U^{(2)}+𝒪(g_2^3).$$
(3.31)
As explained in the beginning of this section we restrict our attention to symmetric and real matrices. This is motivated in part by the fact that the symmetric and real choice uniquely leads to exact agreement with string field theory via (1.8) to order $`g_2`$ as we will see below. Clearly, having a better understanding of why this choice works is very desirable. Therefore, by assuming that $`U`$ is a symmetric and real matrix, we need to solve (3.30). Solving this equation order by order we get

$`U^{(1)}={\displaystyle \frac{1}{2}}G^{(1)},`$ (3.32)
$`U^{(2)}={\displaystyle \frac{1}{2}}G^{(2)}+{\displaystyle \frac{3}{8}}\left(G^{(1)}\right)^2.`$ (3.33)
In the new orthonormal basis, $`\mathrm{\Gamma }`$ is transformed to

$$\stackrel{~}{\mathrm{\Gamma }}=U\mathrm{\Gamma }U^{}.$$
(3.34)
We can determine $`\mathrm{\Gamma }`$ order by order in $`g_2`$ by expanding

$$\stackrel{~}{\mathrm{\Gamma }}=\stackrel{~}{\mathrm{\Gamma }}^{(0)}+g_2\stackrel{~}{\mathrm{\Gamma }}^{(1)}+g_2^2\stackrel{~}{\mathrm{\Gamma }}^{(2)}+𝒪(g_2^3).$$
(3.35)
The matrix of anomalous dimensions in the new basis is therefore

$`\stackrel{~}{\mathrm{\Gamma }}^{(0)}=\mathrm{\Gamma }^{(0)},`$ (3.36)
$`\stackrel{~}{\mathrm{\Gamma }}^{(1)}=\mathrm{\Gamma }^{(1)}{\displaystyle \frac{1}{2}}\{G^{(1)},\mathrm{\Gamma }^{(0)}\},`$ (3.37)
$`\stackrel{~}{\mathrm{\Gamma }}^{(2)}=\mathrm{\Gamma }^{(2)}{\displaystyle \frac{1}{2}}\{G^{(2)},\mathrm{\Gamma }^{(0)}\}{\displaystyle \frac{1}{2}}\{G^{(1)},\mathrm{\Gamma }^{(1)}\}+{\displaystyle \frac{3}{8}}\{\left(G^{(1)}\right)^2,\mathrm{\Gamma }^{(0)}\}+{\displaystyle \frac{1}{4}}G^{(1)}\mathrm{\Gamma }^{(0)}G^{(1)}.`$ (3.38)
Using (3.28) and (3.29) we can evaluate (3.37) to be

$$\stackrel{~}{\mathrm{\Gamma }}^{(1)}=\left(\begin{array}{ccc}0& \stackrel{~}{\mathrm{\Gamma }}_{n,qz}^{(1)}& \stackrel{~}{\mathrm{\Gamma }}_{n,z}^{(1)}\\ \stackrel{~}{\mathrm{\Gamma }}_{py,m}^{(1)}& 0& 0\\ \stackrel{~}{\mathrm{\Gamma }}_{y,m}^{(1)}& 0& 0\end{array}\right),$$
(3.39)
where

$`\stackrel{~}{\mathrm{\Gamma }}_{n,py}^{(1)}=\stackrel{~}{\mathrm{\Gamma }}_{py,n}^{(1)}=\lambda ^{}{\displaystyle \frac{\sqrt{1y}}{\sqrt{Jy}}}{\displaystyle \frac{\mathrm{sin}^2(\pi ny)}{2\pi ^2}},`$ (3.40)
$`\stackrel{~}{\mathrm{\Gamma }}_{n,y}^{(1)}=\stackrel{~}{\mathrm{\Gamma }}_{y,n}^{(1)}=\lambda ^{}{\displaystyle \frac{1}{\sqrt{J}}}{\displaystyle \frac{\mathrm{sin}^2(\pi ny)}{2\pi ^2}}.`$ (3.41)
We note that after using the proposed holographic map (1.8) that the gauge theory results (3.40), (3.41) match with the string field theory results (2.19), (2.20)151515As explained in a footnote in page 4, in order compare the string field theory answer with the gauge theory answer, one must divide the string result by $`\sqrt{Jy(1y)}`$ so that both string field theory states and gauge theory states have unit norm..

We can now make a prediction about the order $`g_2^2`$ matrix elements in string field theory by using (1.8). In order to do that we must calculate the matrix of anomalous dimensions in the new basis to order $`g_2^2`$. By using (3.28), (3.29) we can perform the sums in (3.38) to get161616The formulas we need to compute the required sums are summarized in the Appendix.

$`\stackrel{~}{\mathrm{\Gamma }}^{(2)}=\left(\begin{array}{ccc}\stackrel{~}{\mathrm{\Gamma }}_{n,m}^{(2)}& 0& 0\\ 0& \mathrm{?}& \mathrm{?}\\ 0& \mathrm{?}& \mathrm{?}\end{array}\right),`$ (3.42)
where171717Numerically, this expression is the same as the the one in  for the non-nearest neighbor genus 1 single-trace two-point function.

$`\stackrel{~}{\mathrm{\Gamma }}_{n,m}^{(2)}=\{\begin{array}{cc}{\displaystyle \frac{\lambda ^{}}{32\pi ^4}}\left({\displaystyle \frac{3}{nm}}+{\displaystyle \frac{1}{(nm)^2}}\right)\hfill & \text{ if }nm,m\hfill \\ {\displaystyle \frac{\lambda ^{}}{16\pi ^2}}\left({\displaystyle \frac{1}{3}}+{\displaystyle \frac{5}{2\pi ^2n^2}}\right)\hfill & \text{ if }n=m\hfill \\ {\displaystyle \frac{15\lambda ^{}}{128\pi ^4n^2}}\hfill & \text{ if }n=m\hfill \end{array}`$ (3.43)
and $`\mathrm{?}`$ are quantities that we cannot determine since the full matrix of two-point functions has not been computed to order $`g_2^2`$. We can nevertheless make the following prediction

$$\frac{1}{\mu }n|H|m|_{g_2^2}=\stackrel{~}{\mathrm{\Gamma }}_{n,m}^{(2)}.$$
(3.44)
We should note that this quantity has not been yet been computed in light-cone string field theory. However, (3.44) exactly agrees with the recent proposal in  for the Hamiltonian matrix elements of the string bit Hamiltonian.

## 4 Discussion

In this paper we have studied the gauge theory realization of string field theory Hamiltonian matrix elements. The answer, which was already anticipated in , is that these matrix elements correspond to matrix elements of the dilatation operator. Using the corrected string field theory results in section $`2`$ we find a preferred basis of states which yields agreement between gauge theory and string theory calculations. Moreover, we make a prediction using a gauge theory computation for the Hamiltonian matrix elements of single string states to order $`g_2^2`$, which have not yet been computed. We, however, report precise agreement with the recent string bit  Hamiltonian calculation presented in . An outcome of the corrected string field theory calculation in section $`2`$ is that the the proposal of  no longer holds. This paper gives evidence that the correct correspondence is between matrix elements of the string Hamiltonian and the dilatation operator in the gauge theory. This proposal suggests that the only observables that can be holographically computed in the plane wave string theory are gauge theory two-point functions . Nevertheless, operator mixing between multi-trace operators contains information about higher point functions in gauge theory. In string field theory, the Hamiltonian matrix elements compute the matrix of anomalous dimensions in the orthonormal basis via (1.8), which can be read from the gauge theory two-point functions.

The mapping between the string field theory and gauge theory Hilbert spaces is non-unique. By comparing the calculation of the dilatation operator matrix elements with the corrected string field theory Hamiltonian matrix elements in section $`2`$ we fixed the ambiguity, which picks a particular basis of states in the gauge theory to be identified with string states. However, it would be very desirable to have a first principles explanation of why this choice is the correct one. In the recent paper  the authors motivate this choice by proving that the string bit Hamiltonian simplifies in this basis, that is, it truncates at finite order in $`g_2`$. It would be desirable to understand the uniqueness of the basis choice more directly. In this choice of basis the gauge theory computation we present in section $`3`$ exactly agrees with the recent calculation  performed at order $`g_2^2`$ using the string bit  Hamiltonian.

A fascinating open problem is to understand more precisely the relation between light-cone string field theory and the string bit formalism. The advantage of the string bit formalism is, as shown in , that the Hamiltonian truncates at order $`g_2^2`$. This however seems to raise a puzzle. The Hamiltonian matrix elements truncate at order $`g_2^2`$, which via the map (1.8) predicts that matrix elements of the dilatation operator in some basis truncate at order $`g_2^2`$ even though the Hilbert space inner product $`G_{AB}`$ and the matrix of anomalous dimensions $`\mathrm{\Gamma }_{AB}`$ have corrections to all orders in $`g_2`$. It would be very interesting to study this prediction in detail. It still remains very desirable, however, to explicitly construct the light-cone string field theory Hamiltonian at order $`g_2^2`$ and explicitly verify that its matrix elements are those computed in (3.44).

Light cone string field theory and the string bit model have some complementary features. In the string bit formalism it is easier to compute the $`g_2^2`$ corrections, while the same problem is notoriously difficult in light-cone string field theory. On the other hand, in light-cone string field theory we can systematically evaluate $`1/\mu `$ corrections, which are hard to obtain in the string bit model. Computing these corrections is crucial in extending the duality beyond leading order in $`\lambda ^{}`$ in the gauge theory. In particular, the results in (2.14), (2.15) and (2.16) make non-trivial predictions about $`𝒪(\lambda ^{})`$ corrections to the matrix of anomalous dimensions via (1.8). The factorization theorems proven in  will be very useful in computing systematically all these corrections. It would also be very desirable to compute these corrections directly in gauge theory.

### Acknowledgment

We would like to thank Hirosi Ooguri and John Schwarz for useful discussions and comments, and Peter Lee for collaboration in the early stage of this work. We would also like to thank Neil Constable, Umut Gursoy, Daniel Freedman, Matthew Headrick and Radu Roiban for comments on the original version of the paper. J. G. is supported by a Sherman Fairchild Prize Fellowship. This research was supported in part by the DOE grant DE-FG03-92-ER40701. S. M. was supported in part by JSPS Postdoctoral Fellowships for Research Abroad H14-472.

## Appendix

In this appendix, we present the definition of several matrices appearing in the main text and some useful formulae necessary for computing the matrix elements of $`\mathrm{\Gamma }_{AB}`$ in the new basis. They have been obtained from .

Matrix elements ($`|m||n|,m0,n0,p𝐙,0<y<1`$)

$``$ $`C_{n,py}=C_{py,n}={\displaystyle \frac{y^{3/2}\sqrt{1y}}{\sqrt{J}\pi ^2}}{\displaystyle \frac{\mathrm{sin}^2(\pi ny)}{(pny)^2}}`$ (4.46)
$`C_{n,y}=C_{y,n}={\displaystyle \frac{1}{\sqrt{J}\pi ^2}}{\displaystyle \frac{\mathrm{sin}^2(\pi ny)}{n^2}}`$
$``$ $`M_{n,n}^1={\displaystyle \frac{1}{60}}{\displaystyle \frac{1}{24\pi ^2n^2}}+{\displaystyle \frac{7}{16\pi ^4n^4}}`$ (4.49)
$`M_{n,n}^1={\displaystyle \frac{1}{48\pi ^2n^2}}+{\displaystyle \frac{35}{128\pi ^4n^4}}`$
$`M_{n,m}^1={\displaystyle \frac{1}{12\pi ^2(nm)^2}}{\displaystyle \frac{1}{8\pi ^4(nm)^4}}+{\displaystyle \frac{1}{4\pi ^4n^2m^2}}+{\displaystyle \frac{1}{8\pi ^4nm(nm)^2}}`$
$``$ $`\mathrm{\Gamma }_{n,py}^{(1)}=\mathrm{\Gamma }_{py,n}^{(1)}=\left({\displaystyle \frac{p^2}{y^2}}{\displaystyle \frac{pn}{y}}+n^2\right)C_{n,py}`$ (4.51)
$`\mathrm{\Gamma }_{n,y}^{(1)}=\mathrm{\Gamma }_{y,n}^{(1)}=n^2C_{n,y}`$
$``$ $`𝒟_{n,n}^1=𝒟_{n,n}^1={\displaystyle \frac{2}{3}}+{\displaystyle \frac{5}{\pi ^2n^2}}`$ (4.53)
$`𝒟_{n,m}^1={\displaystyle \frac{2}{3}}+{\displaystyle \frac{2}{\pi ^2n^2}}+{\displaystyle \frac{2}{\pi ^2m^2}}`$
Useful summation formulae
When one multiplies two matrices in (3.38), the following formulae are useful181818Similar identities can also be found in the Appendix of .:

$`{\displaystyle \underset{p,y}{}}C_{n,py}C_{py,m}={\displaystyle \frac{1}{J\pi ^4}}J{\displaystyle _0^1}dyy^3(1y)\mathrm{sin}^2(\pi ny)\mathrm{sin}^2(\pi my){\displaystyle \underset{p=\mathrm{}}{\overset{\mathrm{}}{}}}{\displaystyle \frac{1}{(pnr)^2(pmr)^2}}`$
$`=\{\begin{array}{cc}{\displaystyle \frac{1}{6\pi ^2(nm)^2}}+{\displaystyle \frac{1}{4\pi ^4n^2m^2}}+{\displaystyle \frac{1}{\pi ^4nm(nm)^2}}{\displaystyle \frac{1}{4\pi ^4(nm)^4}}\hfill & \text{ if }nm\hfill \\ {\displaystyle \frac{1}{30}}{\displaystyle \frac{1}{12\pi ^2n^2}}+{\displaystyle \frac{1}{2\pi ^4n^4}}\hfill & \text{ if }n=m\hfill \end{array}`$ (4.54)
$`{\displaystyle \underset{y}{}}C_{n,y}C_{y,m}={\displaystyle \frac{1}{J\pi ^4}}J{\displaystyle _0^1}dy{\displaystyle \frac{\mathrm{sin}^2(\pi ny)}{n^2}}{\displaystyle \frac{\mathrm{sin}^2(\pi my)}{m^2}}`$
$`=\{\begin{array}{cc}{\displaystyle \frac{1}{4\pi ^4n^2m^2}}\hfill & \text{ if }nm,m\hfill \\ {\displaystyle \frac{3}{8\pi ^4n^4}}\hfill & \text{ if }n=m,m\hfill \end{array}`$ (4.55)
$`{\displaystyle \underset{p,y}{}}{\displaystyle \frac{p}{y}}C_{n,py}C_{py,m}={\displaystyle \frac{1}{J\pi ^4}}J{\displaystyle _0^1}dyy^2(1y)\mathrm{sin}^2(\pi ny)\mathrm{sin}^2(\pi my){\displaystyle \underset{p=\mathrm{}}{\overset{\mathrm{}}{}}}{\displaystyle \frac{p}{(pnr)^2(pmr)^2}}`$
$`=\{\begin{array}{cc}(n+m)\left\{{\displaystyle \frac{1}{12\pi ^2(nm)^2}}+{\displaystyle \frac{1}{4\pi ^4n^2m^2}}+{\displaystyle \frac{1}{8\pi ^4nm(nm)^2}}{\displaystyle \frac{1}{8\pi ^4(nm)^4}}\right\}\hfill & \text{ if }nm\hfill \\ {\displaystyle \frac{n}{30}}{\displaystyle \frac{1}{12\pi ^2n}}+{\displaystyle \frac{7}{8\pi ^4n^3}}\hfill & \text{ if }n=m\hfill \end{array}`$
(4.56)
$`{\displaystyle \underset{p,y}{}}{\displaystyle \frac{p^2}{y^2}}C_{n,py}C_{py,m}={\displaystyle \frac{1}{J\pi ^4}}J{\displaystyle _0^1}dyy(1y)\mathrm{sin}^2(\pi ny)\mathrm{sin}^2(\pi my){\displaystyle \underset{p=\mathrm{}}{\overset{\mathrm{}}{}}}{\displaystyle \frac{p^2}{(pnr)^2(pmr)^2}}`$
$`=\{\begin{array}{cc}{\displaystyle \frac{n^2+m^2}{12\pi ^2(nm)^2}}+{\displaystyle \frac{n^6+m^62nm(n^4+m^4)+n^3m^3}{4\pi ^4n^2m^2(nm)^4}}\hfill & \text{ if }nm\hfill \\ {\displaystyle \frac{n^2}{30}}+{\displaystyle \frac{3}{2\pi ^4n^2}}\hfill & \text{ if }n=m\hfill \end{array}`$ (4.57)


format: text

text:
# 1 Introduction

## 1 Introduction

Quantum graphs have proved an important model in the semiclassical study of systems whose classical analogues are chaotic . Much work in this field concerns the conjecture of Bohigas, Giannoni and Schmit  which connects statistics of the energy level spectrum to those of random matrices, the ensemble of random matrices depending on the symmetries of the system. For systems with no time-reversal symmetry statistics of the Gaussian unitary ensemble (GUE) are expected whereas those of the Gaussian orthogonal (GOE) or symplectic (GSE) ensemble should apply for systems with time-reversal symmetry, depending on whether the spin is integer or half integer, respectively. Most investigations of this conjecture have concentrated on systems with spin zero and there are few examples of systems whose energy level statistics follow the GSE prediction . In line with this the usual graph quantisation applies the Schrödinger operator to a metric graph, finding energy level statistics of the GOE or GUE .

It is our aim to develop a quantisation of graphs including half integer spin and we expect to observe spectral statistics following the GSE when time-reversal symmetry is present. We choose a Dirac operator that we realise as a self-adjoint operator on an appropriate Hilbert space. The Dirac operator on a graph was considered previously by Bulla and Trenkler  as an alternative model of a simple scattering system, however, without addressing the problem of time-reversal invariance. Instead we take closed graphs and find boundary conditions that ensure a self-adjoint realisation of the Dirac operator such that time-reversal symmetry is preserved. By looking at systems of this type we can compare the spectrum with the results for the Schrödinger operator on graphs and general semiclassical results for systems with spin. One well known advantage of the quantum graph is that the trace formula is exact rather than semiclassical. This provides the opportunity to distinguish features present only in the semiclassical limit from those inherent in systems with spin.

One unusual quality of the Dirac operator in one dimension is the possibility of two component spinors rather than the usual four component spinors required in three dimensions. Physically this appears odd. If the Dirac operator on the graph describes an idealisation of a physical system of wires in the limit that the width of the wire tends to zero we are lead to four component spinors of the type required in three dimensions. If however we set up the mathematical problem in which the graph is a topological entity then it is natural to choose two component spinors. In the context of quantising graphs the apparent contradiction will in fact disappear. The spectrum is independent of the choice of spinors so the physics of the system cannot distinguish the language used to describe it.

Our paper is organised as follows: In section 2 we introduce the necessary terminology of graphs. Section 3 sets out the two approaches to the Dirac operator in one dimension, restricting from the Dirac equation in three dimensions which leads to four component spinors or considering an irreducible representation of the Dirac algebra which requires only spinors with two components. In section 4 we find self-adjoint and time-reversal symmetric realisations of the Dirac operator for two component spinors and in section 5 we do the same in the four component case. Numerical examples confirm that even for simple graphs the energy level statistics follow the GSE prediction. We find that the scattering matrix in both formalisms has the same properties, the spectrum is therefore independent of the approach taken. In section 6 we first derive the exact trace formula for the density of states and then use this to calculate the form factor. Making the diagonal approximation we see that the form factor agrees with the GSE result for low $`\tau `$. Finally we suggest how the recipe for quantising a graph with the Dirac operator is consistent with the model of a graph as a three dimensional system of wires, section 7.

## 2 Graphs

We begin with a few general definitions necessary when considering graphs. A graph consists of $`V`$ vertices connected by $`B`$ bonds. The valency $`v_i`$ of a vertex $`i`$ is the number of bonds meeting at the vertex. The topology of a graph can be described using its connectivity matrix $`C`$. This is a $`V\times V`$ matrix with entries

$$C_{ij}:=\{\begin{array}{cc}1& \text{if vertices }i\text{ and }j\text{ are connected},\hfill \\ 0& \text{otherwise},\hfill \end{array}$$
(2.1)
where we have assumed that the graph has at most one bond connecting any pair of vertices. We also suppose the graphs under consideration to be connected, any vertex can be reached from any other by passing down bonds. A graph is directed if the bonds are assigned directions: A bond running from vertex $`i`$ to vertex $`j`$ will be denoted as $`b=(ij)`$.

On a graph one usually considers probabilistic rather than deterministic classical dynamics. This can be described by a stochastic matrix $`M`$ of transition probabilities between the bonds. The matrix $`M=(M_{bc})`$ defines a Markov chain on the graph propagating vectors $`\rho =(\rho _1,\mathrm{},\rho _{2B})^T`$ of probabilities assigned to the bonds, $`\rho M\rho `$. Note that we allow particles to move in either direction on a bond and so there are $`2B`$ possible states of the particle on the graph. (This convention for counting bonds will prove necessary later as plane-wave solutions for the Dirac operator cannot be made to correspond to the directions on a directed graph as in the Schrödinger case.) The non-negative matrix $`M`$ is called irreducible if for every pair $`(b,c)`$ of states there exists a power $`q`$ such that $`(M^q)_{bc}>0`$. If the power $`q`$ can be chosen independent of $`(b,c)`$, such that all elements of the matrix $`M^q`$ are strictly positive, $`M`$ is called primitive. A Markov chain with an irreducible Markov matrix $`M`$ is ergodic, i.e. any initial distribution $`\rho `$ converges in average to the stable distribution $`\rho _0=\frac{1}{2B}(1,\mathrm{},1)^T`$,

$$\underset{N\mathrm{}}{lim}\frac{1}{N}\underset{n=1}{\overset{N}{}}M^n\rho =\rho _0.$$
(2.2)
A Markov matrix corresponding to an ergodic chain has a non-degenerate eigenvalue one. For primitive stochastic matrices all other eigenvalues are smaller in magnitude such that the associated Markov process is mixing. Then an initial distribution converges to the stable distribution,

$$\underset{n\mathrm{}}{lim}M^n\rho =\rho _0,$$
(2.3)
and, moreover,

$$\underset{n\mathrm{}}{lim}(M^n)_{bc}=\frac{1}{2B}.$$
(2.4)
It is this last condition that will be used in section 6.2 to evaluate the form factor in the diagonal approximation. Details on non-negative matrices and Markov chains can be found in .

If the quantum system on the graph is defined by a matrix $`𝕋=(T_{bc})`$ of complex transition amplitudes between the bonds (we will see that $`𝕋`$ is unitary) then

$$M_{(ij)(ki)}=|T_{(ij)(ki)}|^2.$$
(2.5)
$`M_{(ij)(ki)}`$ is the quantum mechanical probability that a particle traveling from $`k`$ to $`i`$ scatters at $`i`$ to travel towards the vertex $`j`$. The unitarity of $`𝕋`$ implies that $`M`$ is doubly stochastic and therefore defines a Markov chain.

While only topological information is necessary to investigate classical dynamics on a graph quantum mechanics requires us to assign lengths, $`L_{(ij)}`$, to bonds, $`(ij)`$. This defines a metric graph. It is natural to consider a bond with length as directed, one vertex lying at zero on the bond and the other at $`L`$, consequently our metric graphs are also directed graphs. In order to avoid degeneracies in the lengths of periodic orbits on the graph we assume the lengths assigned to the bonds are incommensurate, not related by a rational number.

## 3 The Dirac equation in one dimension

In one spatial dimension the Dirac equation is

$$\text{i}\mathrm{}\frac{}{t}\mathrm{\Psi }(x,t)=\left(\text{i}\mathrm{}c\alpha \frac{}{x}+mc^2\beta \right)\mathrm{\Psi }(x,t),$$
(3.1)
where $`\alpha `$ and $`\beta `$ satisfy the relations $`\alpha ^2=\beta ^2=\text{I}`$ and $`\alpha \beta +\beta \alpha =0`$ that define the Dirac algebra. There are two possible interpretations of this equation: Either one views it as a restriction of the Dirac equation in three dimensions, or one considers it as a problem in one dimension from the outset. In the first case this equation originates from an implementation of the Poincaré space-time symmetries in relativistic quantum mechanics such that the notions of spin and of anti-particles can be carried over. If then, for instance, one restricts to the $`y`$ axis the Dirac matrices $`\alpha =\alpha _y`$ and $`\beta `$ are hermitian $`4\times 4`$ matrices that form a reducible representation of the Dirac algebra in one dimension. The second case is void of the physical interpretations deriving from $`3+1`$-dimensional space-time symmetries. One merely considers the Dirac equation (3.1) with a faithful irreducible representation of the Dirac algebra and thus chooses $`\alpha `$ and $`\beta `$ to be hermitian $`2\times 2`$ matrices. For an extensive discussion of the Dirac equation see .

The two cases are naturally connected as a unitary transformation of $`\alpha _y`$ and $`\beta `$ brings them into block diagonal form preserving the algebraic relations. For example the standard (Dirac) representation of the matrices in three dimensions, restricted to the $`y`$ axis, is

$$\alpha =\left(\begin{array}{cccc}0& 0& 0& \text{i}\\ 0& 0& \text{i}& 0\\ 0& \text{i}& 0& 0\\ \text{i}& 0& 0& 0\end{array}\right),\beta =\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 1& 0& 0\\ 0& 0& 1& 0\\ 0& 0& 0& 1\end{array}\right).$$
(3.2)
Now consider the unitary transformation

$$U=\frac{1}{\sqrt{2}}\left(\begin{array}{cccc}1& 1& 0& 0\\ 0& 0& 1& 1\\ 1& 1& 0& 0\\ 0& 0& 1& 1\end{array}\right).$$
(3.3)
Applying $`U`$ to $`\alpha `$ and $`\beta `$ generates a new representation,

$$U\alpha U^1=\left(\begin{array}{cccc}0& \text{i}& 0& 0\\ \text{i}& 0& 0& 0\\ 0& 0& 0& \text{i}\\ 0& 0& \text{i}& 0\end{array}\right),U\beta U^1=\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 1& 0& 0\\ 0& 0& 1& 0\\ 0& 0& 0& 1\end{array}\right),$$
(3.4)
which is obviously reducible. Either of the $`2\times 2`$ blocks represents the Dirac algebra in one dimension irreducibly. It is, however, important to note that passing from the four dimensional to the two dimensional representation the usual interpretations of spin, particles and anti-particles, and time-reversal are lost.

For both two- and four-dimensional representations of the Dirac algebra the Dirac operator reads

$$𝒟:=\text{i}\mathrm{}c\alpha \frac{\mathrm{d}}{\mathrm{d}x}+mc^2\beta .$$
(3.5)
On the real line $`𝒟`$ is defined as an essentially self-adjoint operator on the Hilbert space $`L^2()^n`$ with domain $`C_0^{\mathrm{}}()^n`$, where $`n=2`$ or $`n=4`$.

In order to obtain self-adjoint realisations of the Dirac operator on a compact interval $`I=[0,L]`$ appropriate boundary conditions have to be specified. These can be classified by extensions of the closed symmetric operator $`𝒟^0`$ given as $`𝒟`$ defined on the domain $`W_{2,1}^0(I)^nL^2(I)^n`$. Here $`W_{2,1}^0(I)`$ denotes the Sobolev space of functions $`\phi `$ on $`I`$ that are, along with their (generalised) first derivatives, in $`L^2(I)`$ and satisfy $`\phi (0)=0=\phi (L)`$. The closed symmetric extensions of $`𝒟^0`$ arise from restrictions of $`𝒟`$ defined on the domain $`W_{2,1}(I)^n`$ of $`n`$-component Sobolev spinors with no specified boundary conditions to subspaces on which the skew-Hermitian quadratic form

$$\mathrm{\Omega }(\varphi ,\psi ):=𝒟\varphi ,\psi \varphi ,𝒟\psi $$
(3.6)
vanishes. The maximal isotropic subspaces of $`W_{2,1}(I)^n`$ with respect to $`\mathrm{\Omega }`$, i.e. the maximal subspaces on which $`\mathrm{\Omega }`$ vanishes, then yield domains on which $`𝒟`$ is self-adjoint. This is the approach adopted by Kostrykin and Schrader for the Schrödinger operator  on graphs. Self-adjoint realisations of the Dirac operator on an interval are also classified using a different technique by Alonso and De Vincenzo .

## 4 The Dirac operator for two component spinors

From the two approaches to the problem we begin with the minimum necessary and consider two component spinors first. According to the above we take

$$\alpha =\left(\begin{array}{cc}0& \text{i}\\ \text{i}& 0\end{array}\right),\beta =\left(\begin{array}{cc}1& 0\\ 0& 1\end{array}\right).$$
(4.1)
(The choice of $`\alpha `$ and $`\beta `$ is designed to simplify the calculations.)

Eigenspinors $`\psi =\left(\genfrac{}{}{0pt}{}{\psi _1}{\psi _2}\right)`$ of $`𝒟`$ are solutions of the time independent Dirac equation

$$𝒟\psi (x)=E\psi (x).$$
(4.2)
For positive energies $`E`$ they are hence plane waves of the form

$$\psi _k(x)=\mu (k)\left(\begin{array}{c}1\\ \text{i}\gamma (k)\end{array}\right)\text{e}^{\text{i}kx}+\widehat{\mu }(k)\left(\begin{array}{c}1\\ \text{i}\gamma (k)\end{array}\right)\text{e}^{\text{i}kx}$$
(4.3)
with $`k>0`$ and

$$\gamma (k):=\frac{Emc^2}{\mathrm{}ck},E=\sqrt{(\mathrm{}ck)^2+m^2c^4}.$$
(4.4)
The choice of the coefficients $`\mu (k)`$ and $`\widehat{\mu }(k)`$ depends on the boundary conditions imposed by the self-adjoint realisation of $`𝒟`$. We see immediately that the plane-wave spinor (4.3) is not invariant changing $`x`$ to $`x`$. As the Dirac operator is first order the direction assigned to bonds on the graph becomes significant.

### 4.1 Self-adjoint realisations on graphs

Self-adjoint realisations of the Dirac operator on graphs were considered by Bulla and Trenkler . However, in order to develop a simple form related to the vertex transition matrix, we continue with the approach of Kostrykin and Schrader  for the Schrödinger operator.

The basic idea is to view a graph as a collection of intervals that are glued together according to the connectivity matrix $`C`$. The Hilbert space on which the Dirac operator acts is therefore the direct sum of the Hilbert spaces for each bond,

$$:=\underset{b=1}{\overset{B}{}}L^2\left([0,L_b]\right)^2,$$
(4.5)
where $`b`$ runs over all bonds. The spinors $`\psi `$ therefore consist of the $`B`$ components $`(\psi ^1,\mathrm{},\psi ^B)`$ attached to the bonds, each of which is a two-spinor. The scalar product of spinors on the graph is then

$$\varphi ,\psi =\underset{b=1}{\overset{B}{}}\varphi ^b,\psi ^b_b,$$
(4.6)
where $`\varphi ^b,\psi ^b_b`$ is the scalar product in $`L^2\left([0,L_b]\right)^2`$. Self-adjoint realisations of $`𝒟`$ are constructed in analogy to the case of a single interval by suitably extending the closed symmetric operator $`𝒟^0`$ defined on the domain

$$\underset{b=1}{\overset{B}{}}W_{2,1}^0\left([0,L_b]\right)^2.$$
(4.7)
The domains of its closed symmetric extensions are the isotropic subspaces of

$$\underset{b=1}{\overset{B}{}}W_{2,1}\left([0,L_b]\right)^2$$
(4.8)
with respect to the skew-Hermitian form $`\mathrm{\Omega }`$ as given in (3.6). Integrating by parts yields

$$\mathrm{\Omega }(\varphi ,\psi )=\mathrm{}c\underset{b=1}{\overset{B}{}}(\psi _1^b(0)\overline{\varphi }_2^b(0)\psi _1^b(L_b)\overline{\varphi }_2^b(L_b)\psi _2^b(0)\overline{\varphi }_1^b(0)+\psi _2^b(L_b)\overline{\varphi }_1^b(L_b)),$$
(4.9)
showing that $`\mathrm{\Omega }`$ only depends on the boundary values of the spinor components. Self-adjoint extensions again correspond to maximally isotropic subspaces and thus can be characterised by boundary conditions.

Following the construction of Kostrykin and Schrader  for Schrödinger operators on graphs, we introduce a map from the infinite-dimensional Hilbert space $``$ to the $`4B`$-dimensional space of boundary values, $`\psi 𝝍=\left(\begin{array}{c}𝝍_1\\ 𝝍_2\end{array}\right)`$, with

$$\begin{array}{cc}\hfill 𝝍_1& :=(\psi _1^1(0),\mathrm{},\psi _1^B(0),\psi _1^1(L_1),\mathrm{},\psi _1^B(L_B))^T,\hfill \\ \hfill 𝝍_2& :=(\psi _2^1(0),\mathrm{},\psi _2^B(0),\psi _2^1(L_1),\mathrm{},\psi _2^B(L_B))^T.\hfill \end{array}$$
(4.10)
By including $`\psi _2(0)`$ with a minus sign we can write $`\mathrm{\Omega }`$ as

$$\mathrm{\Omega }(\varphi ,\psi )=(\begin{array}{cc}\mathit{\varphi }_1^{}& \mathit{\varphi }_2^{}\end{array})\left(\begin{array}{cc}0& \text{I}_{2B}\\ \text{I}_{2B}& 0\end{array}\right)\left(\begin{array}{c}𝝍_1\\ 𝝍_2\end{array}\right).$$
(4.11)
The maximal isotropic subspaces in (4.8) with respect to $`\mathrm{\Omega }`$ are now equivalent to maximal subspaces of vectors $`𝝍^{4B}`$ on which the complex symplectic form (4.11) vanishes. Defining a linear subspace of $`^{4B}`$ as the vectors $`𝝍`$ satisfying

$$𝔸𝝍_1+𝔹𝝍_2=0,$$
(4.12)
with complex $`2B\times 2B`$ matrices $`𝔸`$ and $`𝔹`$, one first observes that in order for this subspace to have the desired dimension $`2B`$ the $`2B\times 4B`$ matrix $`(𝔸,𝔹)`$ must have maximal rank. Kostrykin and Schrader show that such a subspace is maximally isotropic if and only if $`(𝔸,𝔹)`$ has maximal rank and $`𝔸𝔹^{}`$ is hermitian.

Equation (4.12) defines boundary conditions for the Dirac operator on a collection of $`B`$ bonds. These boundary conditions yield a self-adjoint realisation of $`𝒟`$ when

$$\mathrm{rank}(𝔸,𝔹)=2B\text{and}𝔸𝔹^{}=𝔹𝔸^{}.$$
(4.13)
So far the connectivity of the graph has played no role. In order to implement this we consider only boundary conditions, specified through the matrices $`𝔸`$ and $`𝔹`$, that connect bonds according to the connectivity matrix $`C`$. The matrices $`𝔸`$ and $`𝔹`$ then have a block structure, each block defining boundary conditions at a single vertex.

At a given vertex $`i`$ we label the $`v_i`$ bonds such that the first $`w_i`$ bonds are outgoing and the remainder incoming, see figure 1. The vector $`𝝍^{(i)}=\left(\begin{array}{c}𝝍_1^{\left(i\right)}\\ 𝝍_2^{\left(i\right)}\end{array}\right)^{2v_i}`$ of boundary values at $`i`$ then consists of the components

$$\begin{array}{cc}\hfill 𝝍_1^{(i)}& =(\psi _1^1(0),\mathrm{},\psi _1^{w_i}(0),\psi _1^{w_i+1}(L_{w_i+1}),\mathrm{},\psi _1^{v_i}(L_{v_i}))^T,\hfill \\ \hfill 𝝍_2^{(i)}& =(\psi _2^1(0),\mathrm{},\psi _2^{w_i}(0),\psi _2^{w_i+1}(L_{w_i+1}),\mathrm{},\psi _2^{v_i}(L_{v_i}))^T.\hfill \end{array}$$
(4.14)
If we denote the $`v_i\times v_i`$ blocks of the matrices $`𝔸`$ and $`𝔹`$ that define the boundary conditions at $`i`$ by $`𝔸^{(i)}`$ and $`𝔹^{(i)}`$, we have the condition

$$𝔸^{(i)}𝝍_1^{(i)}+𝔹^{(i)}𝝍_2^{(i)}=0\text{ with }𝔸^{(i)}𝔹^{(i)}=𝔹^{(i)}𝔸^{(i)}$$
(4.15)
which guarantees current conservation at the vertex $`i`$.

Eigenspinors $`\psi `$ of the Dirac operator $`𝒟`$ with given boundary conditions are composed of plane waves, equation (4.3), on each of the bonds. Assigning incoming and outgoing plane-wave solutions to each of the bonds we can write down vectors of coefficients $`\stackrel{\mathbf{}}{𝝁}`$, $`\stackrel{\mathbf{}}{𝝁}`$ for the outgoing and incoming waves at the vertex, respectively:

$$\stackrel{\mathbf{}}{𝝁}=\left(\begin{array}{c}\mu ^1\\ \mathrm{}\\ \mu ^{w_i}\\ \widehat{\mu }^{w_i+1}\text{e}^{\text{i}kL_{w_i+1}}\\ \mathrm{}\\ \widehat{\mu }^{v_i}\text{e}^{\text{i}kL_{v_i}}\end{array}\right),\stackrel{\mathbf{}}{𝝁}=\left(\begin{array}{c}\widehat{\mu }^1\\ \mathrm{}\\ \widehat{\mu }^{w_i}\\ \mu ^{w_i+1}\text{e}^{\text{i}kL_{w_i+1}}\\ \mathrm{}\\ \mu ^{v_i}\text{e}^{\text{i}kL_{v_i}}\end{array}\right).$$
(4.16)
Here we omit the index $`i`$ labeling the vertex. Writing $`𝝍_1^{(i)}`$ and $`𝝍_2^{(i)}`$ in terms of $`\stackrel{\mathbf{}}{𝝁}`$ and $`\stackrel{\mathbf{}}{𝝁}`$ and substituting into the definition of the boundary conditions we obtain

$$\stackrel{\mathbf{}}{𝝁}=(𝔸^{(i)}+\text{i}\gamma (k)𝔹^{(i)})^1(𝔸^{(i)}\text{i}\gamma (k)𝔹^{(i)})\stackrel{\mathbf{}}{𝝁}.$$
(4.17)
Notice that under the conditions given the matrices $`(𝔸^{(i)}\pm \text{i}\gamma (k)𝔹^{(i)})`$ are always invertible . This defines the vertex transition matrix $`𝕋^{(i)}`$ which connects the coefficients of incoming and outgoing plane waves at the vertex $`i`$,

$$𝕋^{(i)}:=(𝔸^{(i)}+\text{i}\gamma (k)𝔹^{(i)})^1(𝔸^{(i)}\text{i}\gamma (k)𝔹^{(i)}).$$
(4.18)
The condition that $`𝔸^{(i)}𝔹^{(i)}`$ is hermitian implies that the matrix $`𝕋^{(i)}`$ is unitary,

$$(𝕋^{(i)})^1=(𝕋^{(i)})^{}.$$
(4.19)
Self-adjoint realisations of the Dirac operator on the graph are defined by matrices $`𝔸,𝔹`$ with $`𝔸𝔹^{}`$ hermitian. These prescribe the boundary conditions (4.12) and determine unitary vertex transition matrices (4.18) for plane waves on the graph.

### 4.2 Time-reversal symmetry

Following the conjecture of Bohigas, Giannoni and Schmit, we expect to find energy level statistics corresponding to those of the Gaussian symplectic ensemble (GSE) for quantum systems with half-integer spin $`s`$ and time-reversal symmetry. This is due to the fact that the time-reversal operator is anti-unitary and squares to $`(1)^{2s}`$, imprinting a symplectic symmetry on the Hamiltonian. In the case of a two component Dirac equation in one dimension the usual interpretations of spin and time-reversal inherited from three dimensions are lost. Nevertheless, one can introduce an anti-unitary operator squaring to $`1`$ that is formally identical to the time-reversal operator for two component (Pauli-) spinors in three dimensions,

$$𝒯:=\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right)𝒦,$$
(4.20)
where $`𝒦`$ is complex conjugation. This can be interpreted as a physical time-reversal in the sense that it changes the sign of time and momentum. The matrix part of $`𝒯`$ is only required in order that $`𝒯^2=\text{I}`$ as there is no spin operator for two component spinors. In the subsequent work we will therefore refer to (4.20) as a time-reversal operator. A full discussion of time-reversal symmetry with spin is found in the second chapter of .

We now determine conditions under which the vertex transition matrix is time-reversal symmetric. This will raise questions about when general time-reversal symmetric boundary conditions exist for a Dirac operator on a graph. For the system to be time-reversal symmetric $`𝒯`$ must commute with the Hamiltonian $`𝒟`$. This forces the mass term in $`𝒟`$ to vanish. From now on we consider only the case $`m=0`$ for which $`\gamma (k)=1`$ and time-reversal symmetry is possible with two component spinors.

Applying $`𝒯`$ to the plane-wave solutions on a bond (4.3) we find

$$𝒯\psi _k(x)=\text{i}\overline{\mu }\left(\begin{array}{c}1\\ \text{i}\end{array}\right)\text{e}^{\text{i}kx}+\text{i}\overline{\widehat{\mu }}\left(\begin{array}{c}1\\ \text{i}\end{array}\right)\text{e}^{\text{i}kx}.$$
(4.21)
Time-reversing spinors on all the bonds meeting at the vertex $`i`$ defines new vectors of outgoing and incoming waves $`\stackrel{\mathbf{}}{𝝁}_𝒯,\stackrel{\mathbf{}}{𝝁}_𝒯`$,

$$\stackrel{\mathbf{}}{𝝁}_𝒯=\text{i}\left(\begin{array}{cc}\text{I}_{w_i}& 0\\ 0& \text{I}_{(v_iw_i)}\end{array}\right)\overline{(\stackrel{\mathbf{}}{𝝁})},\stackrel{\mathbf{}}{𝝁}_𝒯=\text{i}\left(\begin{array}{cc}\text{I}_{w_i}& 0\\ 0& \text{I}_{(v_iw_i)}\end{array}\right)\overline{(\stackrel{\mathbf{}}{𝝁})}.$$
(4.22)
For the vertex transition matrix to be time-reversal invariant we require

$$\stackrel{\mathbf{}}{𝝁}_𝒯=𝕋^{(i)}\stackrel{\mathbf{}}{𝝁}_𝒯.$$
(4.23)
Using the definition of the vertex transition matrix, $`\stackrel{\mathbf{}}{𝝁}=𝕋^{(i)}\stackrel{\mathbf{}}{𝝁}`$ for $`𝕋^{(i)}`$ unitary, equation (4.23) implies

$$(𝕋^{(i)})^T=\left(\begin{array}{cc}\text{I}_{w_i}& 0\\ 0& \text{I}_{(v_iw_i)}\end{array}\right)𝕋^{(i)}\left(\begin{array}{cc}\text{I}_{w_i}& 0\\ 0& \text{I}_{(v_iw_i)}\end{array}\right).$$
(4.24)
Equation (4.24) is the condition that all vertex transition matrices must satisfy in order that the system possesses time-reversal symmetry.

Splitting $`𝕋^{(i)}`$ into four blocks equation (4.24) is equivalent to

$$𝕋^{(i)}=\left(\begin{array}{cc}𝕋_1& 𝕋_2\\ 𝕋_3& 𝕋_4\end{array}\right)=\left(\begin{array}{cc}𝕋_1^T& 𝕋_3^T\\ 𝕋_2^T& 𝕋_4^T\end{array}\right).$$
(4.25)
The components $`T_{bc}`$ of the transition matrix $`𝕋^{(i)}`$ are thus either symmetric or antisymmetric according to the alignment of the bonds at the vertex, in the sense that:

| $`T_{bc}=T_{cb}`$ |  |  |  |  |
| --- | --- | --- | --- | --- |
| $`T_{bc}=T_{cb}`$ |  |  | or |  |

This raises a number of problems. Back scattering is not possible in this scheme as $`T_{bb}`$ is identically zero (returning to the same bond always requires antisymmetry). Consequently it is not possible to quantise graphs containing vertices with valency one. It is also clear that the components of $`𝕋^{(i)}`$ cannot be invariant under a permutation of the bonds. But it is this natural physical property imposed by Kottos and Smilansky  that allows them to derive a simple form of the trace formula for the Schrödinger operator on a graph. In fact we may ask if any graphs exist for which the vertex transition matrices are both unitary and satisfy the time reversal condition (4.24). There are no examples of $`𝕋^{(i)}`$ for vertices with valency one and only a trivial solution for one incoming and one outgoing bond, $`v_i=2`$. With $`v_i=3`$ there are again no solutions. The simplest nontrivial example of a transition matrix is for a vertex with valency four. In this case, with two incoming and two outgoing bonds, one possible choice of a time-reversal invariant transition matrix is

$$𝕋^{(i)}=\frac{1}{\sqrt{2}}\left(\begin{array}{cccc}0& 1& 1& 0\\ 1& 0& 0& 1\\ 1& 0& 0& 1\\ 0& 1& 1& 0\end{array}\right).$$
(4.26)
Clearly this is both a unitary matrix and has the symmetry and antisymmetry properties prescribed by (4.25). Physically this would represent strange boundary conditions as a spinor reaching the vertex would be prevented not only from back scattering but also from scattering down another of the bonds. However it is still reasonable to ask how closely the spectrum of a graph quantised with such a vertex transition matrix will follow the GSE prediction. In the following subsection we investigate one example of a binary graph (i.e. where all vertices have two incoming and two outgoing bonds).

### 4.3 Bond scattering matrix

To find the energy levels we must first determine the bond scattering matrix $`𝕊`$ for the whole graph from the transition matrices at the vertices. Let us define $`\mu ^{(ij)}`$ to be the coefficient of the plane wave on the bond $`(ij)`$ traveling in the direction $`ij`$. Then $`\mu ^{(ji)}`$ is $`\widehat{\mu }^{(ij)}`$ in our previous notation. The bond scattering matrix $`𝕊=(S_{(ij)(lm)})`$ is defined by

$$S_{(ij)(lm)}(k)=\delta _{mi}T_{(ij)(lm)}^{(m)}\text{e}^{\text{i}kL_{(lm)}}.$$
(4.27)
To scatter from the spinor with coefficient $`\mu ^{(lm)}`$ to that with $`\mu ^{(ij)}`$ the bonds $`(lm)`$ and $`(ij)`$ must be connected at $`m`$. The transition amplitude is then defined by the transition matrix $`𝕋^{(m)}`$ at $`m`$. Before the transition the spinor collects a phase propagating along $`(lm)`$. Equation (4.27) is equivalent to the description of the scattering matrix as a product $`𝕊(k)=𝕋𝔻(k)`$ where $`𝔻(k)`$ is a diagonal matrix of phases and $`𝕋`$ is a matrix of transition amplitudes for the whole graph. As there are two coefficients for the spinors on each bond $`𝕊`$ is a $`2B\times 2B`$ matrix. The unitarity of the vertex transition matrices and the diagonal matrix of phases $`𝔻(k)`$ ensures that $`𝕊(k)`$ is also unitary. Note that the vertex transition matrices $`𝕋^{(m)}`$ and hence the bond S-matrix still depend on the directions assigned to the bonds.

The bond S-matrix $`𝕊(k)`$ acts on the vector of coefficients $`𝝁`$, which is composed of the coefficients $`\stackrel{\mathbf{}}{𝝁}`$ and $`\stackrel{\mathbf{}}{𝝁}`$ at every vertex and therefore defines the plane wave on the graph. Energy eigenfunctions correspond to vectors $`𝝁`$ where

$$𝕊(k)𝝁=𝝁.$$
(4.28)
The energy eigenvalues then correspond to the values of $`k`$ for which

$$\left|\text{I}_{2B}𝕊(k)\right|=0.$$
(4.29)
This quantisation condition for the graph is the same as that for the Schrödinger operator, see .

Figure 2 shows the nearest neighbour level spacing statistics for a fully connected pentagon calculated using the vertex transition matrices (4.26). Directions were assigned to the bonds to produce two incoming and two outgoing bonds at each vertex. The assignment of $`𝕋^{(i)}`$ to the vertices is still not unique as the elements of $`𝕋^{(i)}`$ also vary between bonds with the same direction. The numerics confirm that, even for this small graph with unusual boundary conditions, the energy level statistics correspond well to those of random matrices from the GSE. We now turn to consider how the quantisation of more general graphs can be realised using two component spinors and whether there exists a trace formula for the Dirac operator.

### 4.4 Paired bonds

We expect it should be possible to put a Dirac operator on any (topological) graph while preserving time-reversal symmetry; and the assignment of elements in the transition matrix should follow a general scheme independent of the particular bonds to which they correspond. We saw that difficulties arise due to the fact that, as opposed to the case of a Schrödinger operator, the Dirac operator is a first order differential operator, thus breaking the invariance under $`xx`$. In this section we show that a quantisation in terms of a Dirac operator can be achieved for two component spinors by replacing each bond of the classical (topological) graph with a pair of metric bonds one running in each direction.

We maintain the realisation of a two component Dirac operator as described above, however with the constraint that each vertex $`i`$ has valency $`2v_i`$ with pairs of directed bonds meeting at $`i`$ as indicated in figure 3. The bonds in a pair replace a single bond of the given topological graph and hence are assigned the same length but run in opposite directions. This restores the physical symmetry of the problem as all the original bonds are treated equivalently. Putting plane-wave solutions, equation (4.3), on the bonds we can write an eigenfunction on a pair $`b`$ of bonds as

$$\begin{array}{cc}\hfill \psi _k^b(x_b)& =\mu _\alpha ^b\left(\begin{array}{c}1\\ \text{i}\end{array}\right)\text{e}^{\text{i}kx_b}+\mu _\beta ^b\left(\begin{array}{c}1\\ \text{i}\end{array}\right)\text{e}^{\text{i}kx_b}\hfill \\ & +\widehat{\mu }_\alpha ^b\left(\begin{array}{c}1\\ \text{i}\end{array}\right)\text{e}^{\text{i}kx_b}+\widehat{\mu }_\beta ^b\left(\begin{array}{c}1\\ \text{i}\end{array}\right)\text{e}^{\text{i}kx_b},\hfill \end{array}$$
(4.30)
where $`x_{v_i+b}`$ has been replaced by $`L_bx_b`$ and the extra phases absorbed into the coefficients. The coefficients have been labeled $`\alpha `$ or $`\beta `$ according to the type of spinor, characterised by a lower component i or $`\text{i}`$, rather than the bond on which the spinor travels. Now let us rearrange the coefficients in the vectors $`\stackrel{\mathbf{}}{𝝁}`$ and $`\stackrel{\mathbf{}}{𝝁}`$ so that $`\alpha `$ coefficients are always listed first,

$$\stackrel{\mathbf{}}{𝝁}=\left(\begin{array}{c}\mu _\alpha ^1\\ \mu _\beta ^1\\ \mathrm{}\end{array}\right),\stackrel{\mathbf{}}{𝝁}=\left(\begin{array}{c}\widehat{\mu }_\alpha ^1\\ \widehat{\mu }_\beta ^1\\ \mathrm{}\end{array}\right).$$
(4.31)
The self-adjoint boundary conditions still define a unitary vertex transition matrix,

$$\stackrel{\mathbf{}}{𝝁}=𝕋^{(i)}\stackrel{\mathbf{}}{𝝁}\text{and}(𝕋^{(i)})^1=(𝕋^{(i)})^{}.$$
(4.32)
Applying the time-reversal operator (4.20) to the wavefunctions (4.30) we find time-reversed vectors $`\stackrel{\mathbf{}}{𝝁}_𝒯,\stackrel{\mathbf{}}{𝝁}_𝒯`$,

$$\stackrel{\mathbf{}}{𝝁}_𝒯=\text{i}\left(\begin{array}{c}\overline{\widehat{\mu }_\beta ^1}\\ \overline{\widehat{\mu }_\alpha ^1}\\ \mathrm{}\end{array}\right),\stackrel{\mathbf{}}{𝝁}_𝒯=\text{i}\left(\begin{array}{c}\overline{\mu _\beta ^1}\\ \overline{\mu _\alpha ^1}\\ \mathrm{}\end{array}\right).$$
(4.33)
Time-reversal invariance again implies

$$\stackrel{\mathbf{}}{𝝁}_𝒯=𝕋^{(i)}\stackrel{\mathbf{}}{𝝁}_𝒯.$$
(4.34)
Using the definition of the transition matrix (4.32) we obtain the following condition for time-reversal invariance,

$$(𝕋^{(i)})^T=\left(\begin{array}{ccc}J^1& & \\ & \mathrm{}& \\ & & J^1\end{array}\right)𝕋^{(i)}\left(\begin{array}{ccc}J& & \\ & \mathrm{}& \\ & & J\end{array}\right),$$
(4.35)
where

$$J:=\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right).$$
(4.36)
The form of equation (4.35) is similar to, but not the same as, the definition of a symplectic matrix.

The vertex transition matrix $`𝕋^{(i)}`$ can be divided into $`2\times 2`$ blocks $`(𝕋^{(i)})^{bc}`$ which relate incoming spinors on the pair of bonds $`b`$ to outgoing spinors on the bond pair $`c`$. Let

$$(𝕋^{(i)})^{bc}=\left(\begin{array}{cc}r& s\\ v& w\end{array}\right),$$
(4.37)
then time-reversal symmetry implies

$$(𝕋^{(i)})^{cb}=\left(\begin{array}{cc}w& s\\ v& r\end{array}\right),$$
(4.38)
or equivalently,

$$(𝕋^{(i)})^{cb}=\left|(𝕋^{(i)})^{bc}\right|\left((𝕋^{(i)})^{bc}\right)^1.$$
(4.39)
This suggests both a method of constructing time-reversal symmetric transition matrices and the form of permutation invariance among the bond pairs to be required. Matrices of this form can be written

$$𝕋^{(i)}=\left(\begin{array}{ccc}u_1& & \\ & \mathrm{}& \\ & & u_{v_i}\end{array}\right)\left\{X\left(\begin{array}{cc}1& 0\\ 0& 1\end{array}\right)\right\}\left(\begin{array}{ccc}u_1^1& & \\ & \mathrm{}& \\ & & u_{v_i}^1\end{array}\right),$$
(4.40)
where $`X=(x_{bc})`$ is a symmetric unitary matrix of dimension $`v_i`$ and $`u_j\mathrm{SU}(2)`$. $`|x_{bc}|^2`$ is the probability that any incoming state on the classical bond $`c`$ scatters to any outgoing state on the classical bond $`b`$. We may now choose to implement a physical permutation symmetry between the bonds of the classical graph by requiring the matrix $`X`$ to be invariant under permutations $`\sigma `$,

$$X=\sigma X\sigma ^1.$$
(4.41)
Together with the condition that $`X`$ is unitary this implies

$$X=\text{e}^{\text{i}\theta }\left(\begin{array}{ccc}q1& & q\\ & \mathrm{}& \\ q& & q1\end{array}\right)\text{with}q:=\frac{1+\text{e}^{\text{i}\omega }}{v_i}.$$
(4.42)
Unsurprisingly this is similar to the vertex transition matrix for a Schrödinger operator  which is also symmetric unitary and permutation invariant. Here the coefficient $`\omega `$, which we can think of parameterising the boundary conditions, is independent of $`k`$. The extra phase $`\text{e}^{\text{i}\theta }`$ is not interesting as it corresponds only to an extra $`k`$ independent phase factor. Adjusting the parameter $`\omega `$ we see that the Dirac equation also allows both Neumann like ($`\omega =0`$) and Dirichlet like ($`\omega =\pi `$) boundary conditions. An element $`u_bu_c^1`$ of $`\mathrm{SU}(2)`$ defines a ‘spinor rotation’ from the pair of linearly independent two component spinors incoming on bond $`c`$ to the outgoing pair on bond $`b`$.

Equations (4.40) and (4.42) define general time-reversal invariant vertex transition matrices for two component spinors. The transition amplitudes have a form of permutation invariance (4.41) and in section 6.1 we see that this allows a simple form for the trace formula to be derived. With the vertex transition matrix the bond scattering matrix for the whole graph can be defined as in equation (4.27) with unpaired bonds. The quantisation condition (4.29) is also unchanged, apart from the fact that due to doubling of the number of bonds $`𝕊(k)`$ now is a $`4B\times 4B`$ matrix.

## 5 The Dirac operator for four component spinors

Before looking at examples for the two component spinors on paired bonds we will compare with the alternative construction of the Dirac operator acting on four component spinors. Therefore, we now consider again directed graphs with $`B`$ unpaired bonds. In this case the time independent Dirac equation reads

$$\text{i}\mathrm{}c\alpha \frac{\mathrm{d}\psi }{\mathrm{d}x}+mc^2\beta \psi =E\psi ,$$
(5.1)
where we choose $`\alpha `$ and $`\beta `$ as in (3.2),

$$\alpha =\left(\begin{array}{cccc}0& 0& 0& \text{i}\\ 0& 0& \text{i}& 0\\ 0& \text{i}& 0& 0\\ \text{i}& 0& 0& 0\end{array}\right),\beta =\left(\begin{array}{cccc}1& 0& 0& 0\\ 0& 1& 0& 0\\ 0& 0& 1& 0\\ 0& 0& 0& 1\end{array}\right).$$
(5.2)
We could alternatively have chosen either of the other $`\alpha `$ matrices of the standard representation for the Dirac operator in three dimensions or used another representation entirely. The choice here simplifies the working.

Eigenspinors of $`𝒟`$ are again plane waves. For positive energy they are of the form

$$\begin{array}{cc}\hfill \psi _k^b(x_b)& =\mu _\alpha ^b\left(\begin{array}{c}1\\ 0\\ 0\\ \text{i}\gamma (k)\end{array}\right)\text{e}^{\text{i}kx_b}+\mu _\beta ^b\left(\begin{array}{c}0\\ 1\\ \text{i}\gamma (k)\\ 0\end{array}\right)\text{e}^{\text{i}kx_b}\hfill \\ & +\widehat{\mu }_\alpha ^b\left(\begin{array}{c}1\\ 0\\ 0\\ \text{i}\gamma (k)\end{array}\right)\text{e}^{\text{i}kx_b}+\widehat{\mu }_\beta ^b\left(\begin{array}{c}0\\ 1\\ \text{i}\gamma (k)\\ 0\end{array}\right)\text{e}^{\text{i}kx_b},\hfill \end{array}$$
(5.3)
with $`k>0`$. The dispersion relation and $`\gamma (k)`$ are as in (4.4). These four plane-wave solutions differ from the two component case although they are labeled with coefficients $`\mu _\alpha ^b`$ etc. in analogy with the two component solutions on pairs of bonds. The wavefunction $`\psi ^b(x_b)`$ is a function on a single directed bond of a metric graph. We notice that, in contrast to the pairs of two component wavefunctions (4.30), $`\psi ^b(x_b)`$ is not invariant under a change in the bond direction $`x_bx_b`$.

### 5.1 Self-adjoint realisations on graphs

The construction of self-adjoint realisations of the Dirac operator on the graph proceeds in complete analogy to the case of two component spinors. We therefore must look for maximal closed subspaces of

$$\underset{b=1}{\overset{B}{}}W_{2,1}\left([0,L_b]\right)^4$$
(5.4)
on which the skew-hermitian form

$$\mathrm{\Omega }(\varphi ,\psi )=𝒟\varphi ,\psi \varphi ,𝒟\psi $$
(5.5)
vanishes. In order to convert the problem to the boundary values of the spinors we integrate by parts and express $`\mathrm{\Omega }`$ as a complex symplectic form,

$$\mathrm{\Omega }(\varphi ,\psi )=(\begin{array}{cc}(\mathit{\varphi }^+)^{}& (\mathit{\varphi }^{})^{}\end{array})\left(\begin{array}{cc}0& \text{I}_{4B}\\ \text{I}_{4B}& 0\end{array}\right)\left(\begin{array}{c}𝝍^+\\ 𝝍^{}\end{array}\right),$$
(5.6)
on $`^{8B}`$. Here we have defined

$$\begin{array}{cc}\hfill 𝝍^+& :=(\psi _1^1(0),\mathrm{},\psi _1^B(0),\psi _2^1(0),\mathrm{},\psi _2^B(0),\hfill \\ & \psi _1^1(L_1),\mathrm{},\psi _1^B(L_B),\psi _2^1(L_1),\mathrm{},\psi _2^B(L_B))^T,\hfill \\ \hfill 𝝍^{}& :=(\psi _4^1(0),\mathrm{},\psi _4^B(0),\psi _3^1(0),\mathrm{},\psi _3^B(0),\hfill \\ & \psi _4^1(L_1),\mathrm{},\psi _4^B(L_B),\psi _3^1(L_1),\mathrm{},\psi _3^B(L_B))^T.\hfill \end{array}$$
(5.7)
The map from the spinors on the graph to the vectors $`𝝍^+`$ and $`𝝍^{}`$ mixes the first and second components of the spinors in $`𝝍^+`$ and the third and fourth ones in $`𝝍^{}`$. Identifying maximally isotropic subspaces of $`^{8B}`$ with respect to $`\mathrm{\Omega }`$, equation (5.6), the argument proceeds as for the two component spinors. Such a linear subspace is defined by

$$𝔸𝝍^++𝔹𝝍^{}=0,$$
(5.8)
with $`4B\times 4B`$ matrices $`𝔸`$ and $`𝔹`$, and is maximal isotropic if $`\mathrm{rank}(𝔸,𝔹)=4B`$ and $`𝔸𝔹^{}=𝔹𝔸^{}`$. Under such conditions the Dirac operator on (5.4) with boundary conditions as prescribed by (5.8) is self-adjoint.

For eigenfunctions of the Dirac operator, i.e. plane-wave solutions (5.3) on each of the bonds, we define vectors of incoming and outgoing coefficients at a single vertex, for example figure 1,

$$\stackrel{\mathbf{}}{𝝁}=\left(\begin{array}{c}\mu _\alpha ^1\\ \mu _\beta ^1\\ \mathrm{}\\ \mu _\alpha ^{w_i}\\ \mu _\beta ^{w_i}\\ \widehat{\mu }_\alpha ^{w_i+1}\text{e}^{\text{i}kL_{w_i+1}}\\ \widehat{\mu }_\beta ^{w_i+1}\text{e}^{\text{i}kL_{w_i+1}}\\ \mathrm{}\\ \widehat{\mu }_\alpha ^{v_i}\text{e}^{\text{i}kL_{v_i}}\\ \widehat{\mu }_\beta ^{v_i}\text{e}^{\text{i}kL_{v_i}}\end{array}\right),\stackrel{\mathbf{}}{𝝁}=\left(\begin{array}{c}\widehat{\mu }_\alpha ^1\\ \widehat{\mu }_\beta ^1\\ \mathrm{}\\ \widehat{\mu }_\alpha ^{w_i}\\ \widehat{\mu }_\beta ^{w_i}\\ \mu _\alpha ^{w_i+1}\text{e}^{\text{i}kL_{w_i+1}}\\ \mu _\beta ^{w_i+1}\text{e}^{\text{i}kL_{w_i+1}}\\ \mathrm{}\\ \mu _\alpha ^{v_i}\text{e}^{\text{i}kL_{v_i}}\\ \mu _\beta ^{v_i}\text{e}^{\text{i}kL_{v_i}}\end{array}\right).$$
(5.9)
Then using boundary conditions defined at the vertex by matrices $`𝔸^{(i)}`$ and $`𝔹^{(i)}`$ we find as in the two component case,

$$\stackrel{\mathbf{}}{𝝁}=(𝔸^{(i)}\text{i}\gamma (k)𝔹^{(i)})^1(𝔸^{(i)}+\text{i}\gamma (k)𝔹^{(i)})\stackrel{\mathbf{}}{𝝁}.$$
(5.10)
The vertex transition matrix

$$𝕋^{(i)}:=(𝔸^{(i)}\text{i}\gamma (k)𝔹^{(i)})^1(𝔸^{(i)}+\text{i}\gamma (k)𝔹^{(i)})$$
(5.11)
is unitary due to the condition $`𝔸^{(i)}𝔹^{(i)}=𝔹^{(i)}𝔸^{(i)}`$ that guarantees current conservation at the vertex.

### 5.2 Time-reversal symmetry

We now turn to see what time-reversal invariance implies for the vertex transition matrix. The time-reversal operator in the standard representation is

$$𝒯=\left(\begin{array}{cccc}0& \text{i}& 0& 0\\ \text{i}& 0& 0& 0\\ 0& 0& 0& \text{i}\\ 0& 0& \text{i}& 0\end{array}\right)𝒦.$$
(5.12)
This commutes with the Hamiltonian on the bonds so there is no requirement that the mass be zero with four component spinors. Applying $`𝒯`$ to the wavefunction at the vertex $`i`$ we define new vectors from the coefficients of outgoing and incoming waves, $`\stackrel{\mathbf{}}{𝝁}_𝒯,\stackrel{\mathbf{}}{𝝁}_𝒯`$,

$$\stackrel{\mathbf{}}{𝝁}_𝒯=\text{i}\left(\begin{array}{ccc}J& & \\ & \mathrm{}& \\ & & J\end{array}\right)\overline{(\stackrel{\mathbf{}}{𝝁})},\stackrel{\mathbf{}}{𝝁}_𝒯=\text{i}\left(\begin{array}{ccc}J& & \\ & \mathrm{}& \\ & & J\end{array}\right)\overline{(\stackrel{\mathbf{}}{𝝁})}.$$
(5.13)
For time-reversal invariance we require

$$\stackrel{\mathbf{}}{𝝁}_𝒯=𝕋^{(i)}\stackrel{\mathbf{}}{𝝁}_𝒯.$$
(5.14)
Using the unitarity of $`𝕋^{(i)}`$ time-reversal symmetry (5.14) is equivalent to demanding the transition matrix satisfy

$$(𝕋^{(i)})^T=\left(\begin{array}{ccc}J^1& & \\ & \mathrm{}& \\ & & J^1\end{array}\right)𝕋^{(i)}\left(\begin{array}{ccc}J& & \\ & \mathrm{}& \\ & & J\end{array}\right).$$
(5.15)
This is the same as condition (4.35) derived for two component spinors on paired bonds. The general form of transition matrices (4.40) proposed for the spinors on paired bonds will hence also hold for four component spinors on a directed graph. Furthermore with four component spinors the mass is no longer required to be zero. Instead any $`𝕋^{(i)}`$ which can be constructed from the boundary conditions according to equation (5.11) and which is also time-reversal invariant is available.

## 6 Energy level statistics

To verify the conjecture of Bohigas, Giannoni and Schmit in the case of the time-reversal invariant Dirac operator on a graph we examine the energy level statistics through both numerical calculations and directly via the trace formula.

From now on we consider only boundary conditions with the form of general vertex transition matrices defined in equations (4.40) and (4.42). This still leaves a range of boundary conditions parameterised by $`\omega `$. Of particular interest are the Neumann like boundary conditions, $`\omega =0`$. These are the boundary conditions most often studied for the Schrödinger operator on the graph. The vertex transition matrix $`𝕋^{(i)}`$ for Neumann boundary conditions is generated by the matrices

$$𝔸^{(i)}=U^{(i)}\left\{\left(\begin{array}{ccccc}1& 1& 0& 0& \mathrm{}\\ 0& 1& 1& 0& \mathrm{}\\ & & \mathrm{}& \mathrm{}& \\ 0& \mathrm{}& 0& 1& 1\\ 0& \mathrm{}& 0& 0& 0\end{array}\right)\left(\begin{array}{cc}1& 0\\ 0& 1\end{array}\right)\right\}(U^{(i)})^1,$$
(6.1)
$$𝔹^{(i)}=U^{(i)}\left\{\left(\begin{array}{cccc}0& 0& \mathrm{}& 0\\ \mathrm{}& \mathrm{}& & \mathrm{}\\ 0& 0& \mathrm{}& 0\\ 1& 1& \mathrm{}& 1\end{array}\right)\left(\begin{array}{cc}1& 0\\ 0& 1\end{array}\right)\right\}(U^{(i)})^1,$$
(6.2)
where $`U^{(i)}`$ is a block diagonal matrix of elements of $`\mathrm{SU}(2)`$,

$$U^{(i)}=\left(\begin{array}{ccc}u_1& & \\ & \mathrm{}& \\ & & u_{v_i}\end{array}\right).$$
(6.3)
Substituting these boundary conditions into equation (5.11) and using some algebra we obtain the vertex transition matrix

$$𝕋^{(i)}=U^{(i)}\left\{\left(\begin{array}{ccc}\frac{2}{v_i}1& & \frac{2}{v_i}\\ & \mathrm{}& \\ \frac{2}{v_i}& & \frac{2}{v_i}1\end{array}\right)\left(\begin{array}{cc}1& 0\\ 0& 1\end{array}\right)\right\}(U^{(i)})^1.$$
(6.4)
$`𝕋^{(i)}`$ is independent of the function $`\gamma (k)`$ determined by the mass. With this form of the Neumann boundary conditions there is no need to specify the mass or whether the system is represented by paired two component spinors or four component spinors on a directed graph.

Figure 4 shows the nearest neighbour spacing statistics for a fully connected square with Neumann boundary conditions, the vertex transition matrix defined in equation (6.4). The three elements $`u_b\mathrm{SU}(2)`$ at a vertex were chosen as

$$u_b=\mathrm{exp}(\text{i}\theta _b\sigma _{n_b}),$$
(6.5)
where $`\sigma _{n_b}`$ is a Pauli matrix. The parameters $`\theta _b`$ were then selected randomly at each vertex. To obtain good agreement with the symplectic random matrix ensemble the parameters were also picked such that the elements $`u_b`$ generated at the vertices were all sufficiently different. Deviations from the GSE behaviour become apparent if the elements of $`\mathrm{SU}(2)`$ at one vertex are similar to those at another, i.e. $`\theta _b^l\theta _c^m`$.

Introducing a magnetic vector potential $`A_b`$ on the bonds breaks time-reversal symmetry. The transition elements for the $`S`$-matrix remain the same but a plane wave propagating down a bond $`(ij)`$ picks up an extra phase $`\mathrm{exp}(\text{i}A_{(ij)}L_{(ij)})`$ where we know that $`L_{(ji)}=L_{(ij)}`$ but $`A_{(ji)}=A_{(ij)}`$. Figure 5 shows that for the Dirac operator on a fully connected square breaking the time-reversal invariance in this way produces energy level statistics like those of random GUE matrices as expected. For the Dirac operator we used Neumann boundary conditions and random elements of $`\mathrm{SU}(2)`$ chosen as previously.

### 6.1 Trace formula

A trace formula for the Laplacian on a graph was first produced by Roth . We derive the trace formula for the density of states of a system with Neumann boundary conditions, following the approach developed by Kottos and Smilansky  for the Schrödinger operator. The vertex transition matrices for such a system were defined in equation (6.4). The main tool in setting up the trace formula is the bond scattering matrix $`𝕊(k)`$ defined in equation (4.27).

Consider the density of states for the wave-number $`k`$. It is an infinite series of delta functions located at $`k=k_n`$ corresponding to eigenvalues $`E_n`$ of the Dirac operator. The eigenvalues are counted according to their multiplicity. Wave-numbers $`k_n`$ are values of $`k`$ for which the function

$$\zeta (k):=|\text{I}_{4B}𝕊(k)|$$
(6.6)
is zero. Diagonalising $`𝕊(k)`$ we can write $`\zeta (k)`$ in terms of the eigenvalues $`\text{e}^{\text{i}\varphi _j(k)}`$ of $`𝕊(k)`$,

$$|\text{I}_{4B}𝕊(k)|=\underset{j=1}{\overset{4B}{}}(1\text{e}^{\text{i}\varphi _j(k)})=2^{4B}|𝕊(k)|^{\frac{1}{2}}\underset{j=1}{\overset{4B}{}}\mathrm{sin}\left(\frac{\varphi _j(k)}{2}\right).$$
(6.7)
$`\zeta (k)`$ is complex with the phase contained in the term $`|𝕊(k)|^{\frac{1}{2}}`$. Multiplying by $`|𝕊(k)|^{\frac{1}{2}}`$ we define a real function whose zeros are also the zeros of $`\zeta (k)`$. Then the density of states can be written

$$d(k)=\underset{n=1}{\overset{\mathrm{}}{}}\delta (kk_n)=\frac{1}{\pi }\underset{ϵ0}{lim}\text{Im}\frac{\mathrm{d}}{\mathrm{d}k}\mathrm{log}(|\text{I}_{4B}𝕊(k+\text{i}ϵ)||𝕊(k+\text{i}ϵ)|^{\frac{1}{2}}).$$
(6.8)
The expression (6.8) factorises into a sum of two terms. The smooth part of the spectral density is

$`d_{\text{smth}}(k)`$ $`:=`$ $`{\displaystyle \frac{1}{\pi }}\underset{ϵ0}{lim}\text{Im}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}k}}\mathrm{log}(|𝕊(k+\text{i}ϵ)|^{\frac{1}{2}})`$ (6.9)
$`=`$ $`{\displaystyle \frac{1}{2\pi }}\text{Im}{\displaystyle \frac{\mathrm{d}}{\mathrm{d}k}}\mathrm{log}(\text{e}^{\text{i}4kL})`$
$`=`$ $`{\displaystyle \frac{2L}{\pi }},`$
where $`L`$ is the total length of the graph,

$$L=\underset{b=1}{\overset{B}{}}L_b.$$

(Note that in our definition of $`L`$ the length of each bond of the classical graph is counted only once although the S-matrix can be considered either as scattering between two component spinors on paired bonds or four component spinors on a directed graph.)

The oscillating part of the spectral density is

$$d_{\text{osc}}(k):=\frac{1}{\pi }\underset{ϵ0}{lim}\text{Im}\frac{\mathrm{d}}{\mathrm{d}k}\mathrm{log}(|\text{I}_{4B}𝕊(k+\text{i}ϵ)|).$$
(6.10)
The logarithm can be expanded as a sum of powers of the trace of $`𝕊(k)`$,

$$\mathrm{log}(|\text{I}_{4B}𝕊(k)|)=\underset{n=1}{\overset{\mathrm{}}{}}\frac{1}{n}\mathrm{tr}𝕊^n(k).$$
(6.11)
As in the case of the Schrödinger operator on a graph we can now write the trace of powers of $`𝕊(k)`$ as a sum over periodic orbits on the graph. However, for the Dirac operator, the transition amplitudes between classical bonds are specified by the $`2\times 2`$ matrices $`𝕋^{bc}`$,

$$\mathrm{tr}𝕊^n(k)=\mathrm{tr}\underset{b_1\mathrm{}b_n}{}\text{e}^{\text{i}kL_{b_1}}𝕋^{b_1b_2}\text{e}^{\text{i}kL_{b_2}}𝕋^{b_2b_3}\mathrm{}\text{e}^{\text{i}kL_{b_n}}𝕋^{b_nb_1}.$$
(6.12)
Thus $`\mathrm{tr}𝕊^n(k)0`$ implies that the sequence $`(b_1b_2\mathrm{}b_n)`$ of bonds labels a periodic orbit on the graph. For a single periodic orbit $`p`$ which passes through $`n`$ vertices,

$$\mathrm{tr}(𝕋^{b_1b_2}\mathrm{}𝕋^{b_nb_1})=A_p\text{tr}(d_p)\text{e}^{\text{i}\pi \mu _p}.$$
(6.13)
In this formula $`A_p`$ is a (positive) stability factor, the product of the absolute values of the transition amplitudes for spinor pairs round the orbit,

$$A_p:=\underset{s=1}{\overset{\mu _p}{}}\left(1\frac{2}{v_s}\right)\underset{t=1}{\overset{\nu _p}{}}\left(\frac{2}{v_t}\right),$$
(6.14)
where $`\mu _p`$ is the number of cases of back scattering from vertices with $`v_s>2`$ and $`n=\mu _p+\nu _p`$. An overall sign is collected in the factor $`\text{e}^{\text{i}\pi \mu _p}`$ and $`d_p`$ is the product of the elements of $`\mathrm{SU}(2)`$ that rotate between the two types of spinors on the periodic orbit,

$$d_p:=u^{b_1b_2}u^{b_2b_3}\mathrm{}u^{b_nb_1}.$$
(6.15)
The spin rotations $`u^{b_ib_{i+1}}:=u_{b_{i+1}}u_{b_i}^1`$ are composed of $`\mathrm{SU}(2)`$-blocks comprising the matrix $`U^{(i)}`$, see equation (6.3), attached to the $`i`$th vertex visited along the periodic orbit.

Combining these results,

$$\mathrm{tr}𝕊^n(k)=\underset{pP_n}{}A_p\text{e}^{\text{i}\pi \mu _p}\mathrm{tr}(d_p)\text{e}^{\text{i}kl_p},$$
(6.16)
where the sum is over the set $`P_n`$ of periodic orbits of $`n`$ bonds and $`l_p`$ is the metric length of the orbit $`p`$. Substituting equations (6.16) and (6.11) into equation (6.10) we find the oscillating part of the density of states. Adding the result for the smooth part (6.9) we obtain

$$d(k)=\frac{2L}{\pi }+\frac{1}{\pi }\underset{p}{}\frac{l_p}{r_p}A_p\text{e}^{\text{i}\pi \mu _p}\mathrm{tr}(d_p)\mathrm{cos}(kl_p).$$
(6.17)
The sum is over all periodic orbits. If an orbit $`p`$ consists of $`r_p`$ repetitions of a shorter orbit then only the length of the primitive orbit $`l_p/r_p`$ is included in the formula. This trace formula provides an exact relation for the density of states of the Dirac operator on a graph with Neumann boundary conditions as a sum over the classical periodic orbits.

Equation (6.17) for the density of states can be compared both to the form derived for the Schrödinger operator on a graph  and to previous semiclassical results for systems with spin . Most terms in our density of states are the same as those derived by Kottos and Smilansky for the Schrödinger operator . The additional term for the Dirac operator is the trace of $`d_p`$. This can be thought of as an additional weight factor which depends on the transformation of the spinor pairs around the orbit.

Comparing our trace formula to the general semiclassical form for the Dirac operator , we see that here the Dirac operator also produces a trace over an element of $`\mathrm{SU}(2)`$ associated with the change in spin around the orbit. This weight factor can also be interpreted in terms of an effective rotation angle of a classical spin vector transported along the periodic orbit according to the equation of Thomas precession. That the results agree is not unexpected. The Dirac operator on a graph however does not require a semiclassical approximation to derive the density of states. Finally we also note the similarity with the trace formula for a different system, the quantum cat map with spin, derived by Keppeler, Marklof and Mezzadri, . In this case the quantum map is an equivalent to the Pauli operator with spin 1/2. In their semiclassical trace formula the presence of spin also appears as a trace of an $`\mathrm{SU}(2)`$ matrix which defines the spin transport around a periodic orbit.

### 6.2 Form factor

From the trace formula we derive the spectral two-point form factor for the Dirac operator on quantum graphs. The form factor is an energy level statistic that is widely studied for classically chaotic quantum systems and in particular for the Schrödinger operator on graphs, see . We see that our form factor agrees with the general semiclassical result for systems with spin , and making the diagonal approximation the form factor for low $`\tau `$ is consistent with the GSE prediction from random matrix theory.

We derive the form factor for the energy level spacings. This is the Fourier transform of the two-point correlation function itself defined as

$$R_2(x):=\left(\frac{\pi }{2L}\right)^2\underset{\mathrm{\Lambda }\mathrm{}}{lim}\frac{1}{\mathrm{\Lambda }}_0^\mathrm{\Lambda }d(k)d\left(k\frac{\pi x}{L}\right)dk1.$$
(6.18)
In defining the two-point correlation function we have rescaled the spectrum dividing by the mean level spacing $`2L/\pi `$, see equation (6.9). For systems with half-integer spin and time-reversal invariance the eigenvalues $`k_n`$ come in pairs (Kramers’ degeneracy). We consider the spectrum $`x_n`$ where this degeneracy has been lifted, each pair of eigenvalues being replaced with a single representative. The mean spacing of the $`x_n`$ distribution is therefore $`L/\pi `$. Substituting the trace formula for the density of states (6.17) into the definition of $`R_2(x)`$ and carrying out the integral we find

$$R_2(x)=\frac{1}{2(2L)^2}\underset{p,q}{}\frac{l_pl_q}{r_pr_q}A_pA_q\text{e}^{\text{i}\pi (\mu _p+\mu _q)}\mathrm{tr}(d_p)\mathrm{tr}(d_q)\mathrm{cos}\left(\frac{\pi xl_q}{L}\right)\delta _{l_p,l_q}.$$
(6.19)
The two-point form factor $`K_2(\tau )`$ is the Fourier transform of the two-point correlation function,

$$K_2(\tau ):=_{\mathrm{}}^{\mathrm{}}R_2(x)\text{e}^{2\pi \text{i}x\tau }dx,$$
(6.20)
which is even in $`\tau `$. Taking the Fourier transform of $`R_2(x)`$ term-wise we obtain

$$K_2(\tau )=\frac{1}{4(2L)^2}\underset{p,q}{}\frac{l_pl_q}{r_pr_q}A_pA_q\text{e}^{\text{i}\pi (\mu _p+\mu _q)}\mathrm{tr}(d_p)\mathrm{tr}(d_q)\delta \left(\tau \frac{l_q}{2L}\right)\delta _{l_p,l_q}$$
(6.21)
for $`\tau `$ positive.

The eigenvalue distribution of random matrices averaged over the GSE has a two-point form factor

$$K_2^{GSE}(\tau )=\{\begin{array}{ccc}\frac{1}{2}|\tau |\frac{1}{4}|\tau |\mathrm{log}|1|\tau ||\hfill & & \text{for }|\tau |2\hfill \\ 1\hfill & & \text{for }|\tau |2\hfill \end{array}.$$
(6.22)
Making similar assumptions to those often used in the Schrödinger case we will see that $`K_2(\tau )|\tau |/2`$ for $`|\tau |0`$. This agrees with the random matrix theory prediction of equation (6.22).

The analysis of the form factor based on the trace formula requires us to take a semiclassical limit. The latter can be characterised by going to some regime with increasing mean spectral density. According to the relation (6.9) this corresponds to increasing the total length $`L`$ of the bonds. We realise this limit in passing to graphs with an increasing number $`B`$ of bonds while keeping the mean bond length $`\overline{L}:=L/B`$ constant. As $`\tau 0`$ the form factor should be dominated by the diagonal form factor introduced by Berry ,

$$K_2^{diag}(\tau ):=\frac{1}{4(2L)^2}\underset{p}{}\left(\frac{l_p}{r_p}\right)^2\left\{A_p^2(\mathrm{tr}(d_p))^2+A_pA_{\overline{p}}\mathrm{tr}(d_p)\mathrm{tr}(d_{\overline{p}})\right\}\delta \left(\tau \frac{l_p}{2L}\right),$$
(6.23)
since in this limit only correlations for an orbit $`p`$ with itself and with its time-reversed partner $`\overline{p}`$ are relevant. Combined with the semiclassical limit we now have to consider $`L\mathrm{}`$, $`\tau 0`$ with $`L\tau \mathrm{}`$.

Comparing $`p`$ and $`\overline{p}`$ we see that $`A_p=A_{\overline{p}}`$ and $`d_{\overline{p}}=d_p^1`$. For elements of $`\mathrm{SU}(2)`$ $`\mathrm{tr}(d^1)=\mathrm{tr}(d)`$,

$$K_2^{diag}(\tau )=\frac{1}{2(2L)^2}\underset{p}{}\left(\frac{l_p}{r_p}\right)^2A_p^2(\mathrm{tr}(d_p))^2\delta \left(\tau \frac{l_p}{2L}\right).$$
(6.24)
In the limit $`L\tau \mathrm{}`$ it is the long orbits which dominate the sum and for these the proportion which are repetitions of shorter orbits tends to zero so we ignore periodic orbits with $`r_p1`$. Breaking the sum on periodic orbits into a sum on sets $`P_n`$ of orbits of $`n`$ bonds and using $`l_pn\overline{L}`$,

$$K_2^{diag}(\tau )\frac{1}{2}\underset{n}{}\left(\frac{n}{2B}\right)^2\delta \left(\tau \frac{n}{2B}\right)\underset{pP_n}{}A_p^2(\mathrm{tr}(d_p))^2.$$
(6.25)
We can regard both $`A_p`$ and $`d_p`$ as functions of $`p`$. The elements of $`\mathrm{SU}(2)`$ in $`d_p`$ were selected randomly, according to Haar measure, so $`A_p`$ and $`d_p`$ are uncorrelated. If the number of periodic orbits of length $`n`$ is $`|P_n|`$ then for large $`n`$

$$\frac{1}{|P_n|}\underset{pP_n}{}A_p^2(\mathrm{tr}(d_p))^2\left(\frac{1}{|P_n|}\underset{pP_n}{}A_p^2\right)\left(\frac{1}{|P_n|}\underset{pP_n}{}(\mathrm{tr}(d_p))^2\right).$$
(6.26)
We now consider the two averages separately. The sum over $`P_n`$ of $`A_p^2`$ can be evaluated in the limit $`n\mathrm{}`$ using the classical ergodicity of the graph. Classical motion on the graph is determined by the doubly stochastic matrix $`M=(M_{bc})`$, with the magnitude squared of the transition amplitudes between the bonds as entries,

$$M_{bc}:=|X_{bc}|^2.$$
(6.27)
$`M`$ defines a Markov chain. For a connected graph a given pair $`(b,c)`$ of bonds can be linked by a path of $`q`$ bonds. Thus $`(M^q)_{bc}>0`$ such that $`M`$ is an irreducible matrix and the Markov chain is ergodic. If, in addition, no eigenvalue other than one lies on the unit circle it is mixing which implies

$$\underset{n\mathrm{}}{lim}(M^n)_{bc}=\frac{1}{2B}.$$
(6.28)
The trace of $`M^n`$ therefore approaches unity, but the trace can also be expressed as a sum over the periodic orbits,

$$\mathrm{tr}M^n=n\underset{pP_n}{}A_p^2.$$
(6.29)
The factor $`n`$ counts the cyclic permutations of the orbit. Thus in the limit of long orbits we may replace $`_{pP_n}A_p^2`$ with $`1/n`$ if the Markov process is mixing. It should be noted that while we have used the mixing property of the classical motion on the graph without spin (defined by the matrices $`X`$) it is equivalent to the same property of a $`4B`$ dimensional matrix $`M`$ defined from $`𝕋`$. In this case $`M_{bc}=|T_{bc}|^2`$ and each spinor is treated separately. We defined the mixing property of the classical graph using $`X`$ to make the connection with the amplitudes $`A_p`$ clear.

In the limit of $`n\mathrm{}`$ we replace the average over $`P_n`$ of $`(\mathrm{tr}(d_p))^2`$ with the integral over $`\mathrm{SU}(2)`$ with Haar measure. Evaluating the integral directly we obtain

$$\underset{n\mathrm{}}{lim}\frac{1}{|P_n|}\underset{pP_n}{}(\mathrm{tr}(d_p))^2=_{\mathrm{SU}(2)}(\mathrm{tr}(u))^2du=1.$$
(6.30)
Substituting these results into the diagonal approximation for the form factor (6.25) we find

$$K_2^{diag}(\tau )\frac{1}{2}\underset{n}{}\frac{1}{2B}\frac{n}{2B}\delta \left(\tau \frac{n}{2B}\right)\frac{1}{2}|\tau |$$
(6.31)
in the combined limit $`L\mathrm{}`$, $`\tau 0`$ with $`L\tau \mathrm{}`$. This is in agreement with the GSE prediction for small $`|\tau |`$. The validity of the assumptions made in the argument is of less consequence than that they are of the same type used in other semiclassical approximations, particularly for the Schrödinger operator on a graph. The key additional assumption made to include the spin evolution was equation (6.26), taking the elements $`d_p`$ of $`\mathrm{SU}(2)`$ to be uncorrelated with the weights $`A_p`$ of the orbit.

Figure 6 compares a numerical calculation of the form factor with the GSE prediction (6.22). As the form factor takes the form of a distribution, see equation (6.21), meaningful results are only obtained when it is averaged with a test function. For simplicity we use a function which is one in the region $`\tau \pm 0.1`$ and zero outside. The theoretical prediction for the GSE form factor can also be averaged in the same way. As the computation of the form factor leads to a delta function at $`\tau =0`$ in the limit of the number of levels going to infinity the first few values of $`\tau `$ are omitted in the numerical calculation.

## 7 Gedanken experiment

It is common when studying graphs, in this context, to describe them as an idealised model of a “typical” quantum system. When compared to quantum maps or systems of constant negative curvature the case is certainly strong. However in introducing the Dirac operator to graphs we may have inadvertently weakened the argument and that is what we seek to correct here.

Let us think of the classical graph that we wish to quantise as a network or web of wires set up in our ideal laboratory, see figure 7. To quantise this system using the Dirac operator we must introduce pairs of spinors propagating in both directions along the wires. At the vertices changes in the composition of the spinor pairs are described by elements $`u_{bc}\mathrm{SU}(2)`$, $`c`$ labels the incoming classical bond and $`b`$ the outgoing bond. These elements $`u_{bc}`$ necessarily depend on the pair of wires under consideration, $`u_{bc}u_{de}`$. This is not a feature of the usual Schrödinger quantisation  where all wires meeting at a vertex are treated equivalently. It is therefore fair to ask whether our model web could really distinguish bonds in such a way.

One simple geometric argument connecting the matrices $`u_{bc}`$ to our experimental network is provided by the double covering map $`f:\mathrm{SU}(2)\mathrm{SO}(3)`$. Let $`R(\theta _{bc})\mathrm{SO}(3)`$ describe the rotation from bond $`c`$ to $`b`$ at a vertex. Setting $`u_{bc}=f^1(R(\theta _{bc}))`$, with some choice of the branch $`f^1(\mathrm{SO}(3))`$, we define a set of elements of $`\mathrm{SU}(2)`$ from the architecture of the network. This construction is time-reversal symmetric as $`u_{cb}=u_{bc}^1`$. (To write this in the form used for the vertex boundary conditions (6.4) choose a reference direction at each vertex and take the $`v`$ elements $`u_c`$ to be $`f^1(R(\theta _c))`$ where $`R(\theta _c)`$ rotates from the reference direction to bond $`c`$.) We have not attempted to provide a physical justification for this transformation of spinor pairs at the vertices. Our formulation rather shows that the extent to which bonds must be distinguished when quantising a graph with the Dirac operator remains consistent with the picture of the graph as a simple ideal quantum system. There is sufficient geometrical information in the structure of a graph embedded in three dimensions to determine rotations of the spinors at the vertices.

### Acknowledgments

We would like to thank Stefan Keppeler for useful discussions and Arnd Bäcker for advice with the numerics. This work has been fully supported by the European Commission under the Research Training Network (Mathematical Aspects of Quantum Chaos) no. HPRN-CT-2000-00103 of the IHP Programme.