Nonlinear sigma models with compact hyperbolic target spaces

We explore the phase structure of nonlinear sigma models with target spaces corresponding to compact quotients of hyperbolic space, focusing on the case of a hyperbolic genus-2 Riemann surface. The continuum theory of these models can be approximated by a lattice spin system which we simulate using Monte Carlo methods. The target space possesses interesting geometric and topological properties which are reflected in novel features of the sigma model. In particular, we observe a topological phase transition at a critical temperature, above which vortices proliferate, reminiscent of the Kosterlitz-Thouless phase transition in the O(2) model [1, 2]. Unlike in the O(2) case, there are many different types of vortices, suggesting a possible analogy to the Hagedorn treatment of statistical mechanics of a proliferating number of hadron species. Below the critical temperature the spins cluster around six special points in the target space known as Weierstrass points. The diversity of compact hyperbolic manifolds suggests that our model is only the simplest example of a broad class of statistical mechanical models whose main features can be understood essentially in geometric terms.


Introduction
The breaking of continuous symmetries is accompanied by the appearance of massless Goldstone modes. Fluctuations of these modes destroy long-range order at any finite temperature in dimensions d ≤ 2; this is the statement of the Mermin-Wagner theorem. This does not exclude the possibility of "quasi long-range order", however, in which correlators exhibit power-law, rather exponential decay. This behavior is associated with a continuous phase transition driven by the proliferation of topological defects [3].
We can expect topological transitions to occur in a spin system whenever the fundamental group of the spin space is non-trivial. The XY model is the simplest example in which the spin space of S 1 has fundamental group π 1 (S 1 ) = Z. This is to be compared with the O(3)-invariant Heisenberg model which has π 1 (S 2 ) = 0 and thus does not exhibit quasi-long-range order in d ≤ 2. The first non-trivial example with two-dimensional spin space is the torus T 2 = S 1 × S 1 which is essentially two decoupled copies of the XY model. A more interesting possibility is to consider the two-handled double torus which we study in this paper.
The novelty of the double torus is that, as a consequence of Gauss-Bonnet theorem, it admits a metric with constant negative curvature, JHEP06(2016)145 these models and their g = 0 and g = 1 counterparts is the existence of preferred points in the manifold, an example of which are the Weierstrass points. To define a Weierstrass point, first consider the set of all meromorphic functions which are holomorphic away from a specified point P . P is a Weierstrass point if this set contains more meromorphic functions with poles at P of some specific order than are guaranteed by the Riemann-Roch theorem [4]. We will work with a particularly simple double torus corresponding to a tiling of the hyperbolic plane by regular octagons in which opposite sides are identified. The Weierstrass points corresponding to this particular identification are known [4] to be the center of the octagon together with the midpoints of the sides, as well as the point determined by the eight identified vertices. We will refer to this manifold as the regular double torus.
In this paper we will perform Monte-Carlo simulations of the double torus model on a two dimensional square lattice with periodic boundary conditions. We will study the impact of the non-standard topological and geometric properties of this model on the phase transition. Our results can be summarized as follows: • There is a phase transition at a finite temperature T c ≈ 0.58. Numerical results are consistent with a second order phase transition, but do not exclude the possibility that the phase transition is of infinite order as in the XY model. Our numerical results also do not exclude the possibility of additional phase transitions.
• For temperatures slightly below T c , the spins cluster around one of six special points on the regular double torus, which can be defined as fixed points of the discrete automorphism group. This is quite unlike the XY model, in which there are no preferred points in the target space.
• For temperatures somewhat above T c , vortices of many different topological types appear, their numbers following thermal distributions as one would predict from treating them as free, independent excitations.
The organization of the rest of this paper is as follows. In section 2, we explain the general framework of lattice models with hyperbolic quotients as target spaces and specify the precise model we are interested in. In section 3 we briefly describe our numerical methods and then explain our results, focusing on the points just summarized. We end with a discussion including future directions in section 4.

Hyperbolic tilings and quotient spaces
The n-dimensional hyperbolic space H n is a maximally symmetric Euclidean manifold of constant negative curvature. It can be embedded in R 1,n in a manifestly SO(1, n) symmetric manner, More precisely, H n is the upper sheet of the two-sheeted hyperboloid described by eq. (2.1). We consider non-linear sigma models with target space given by the quotient space H n /Γ

JHEP06(2016)145
where Γ is a discrete subgroup of the orientation-preserving isometries SO + (1, n). In the 'upstairs' picture we can think of the covering space H n as being tiled or tessellated by cells, each of which is related to the fundamental cell by the action of a particular element of the group Γ. We only consider orientation-preserving isometries so that the resulting topological space is orientable. Moreover, we assume that Γ is freely acting so that the resulting quotient space does not have fixed points. The Poincaré ball model maps the hyperboloid H n to the unit ball where the induced metric is given by We will mostly focus on the n = 2 case where the isometries SO + (1, 2) are realized as fractional linear transformations acting on the unit disc in C, with a, b ∈ C such that |a| 2 − |c| 2 = 1. The geodesics in the Poincaré coordinates are either diameters or arcs of circles intersecting orthogonally with the boundary of the disc. These geodesics form the edges of the regular hyperbolic polygons which tile H 2 . A tiling by regular p-gons with q polygons meeting at each vertex exists provided that 1/p+1/q < 1/2. We only consider tilings with p even so that the sides of each polygon can be paired. The subset γ ⊂ Γ of group elements which pair the sides of the polygon are the generators of Γ. Note that γ is not a group and the group Γ is obtained by multiplying the elements of γ in all possible ways, Γ = {g 1 · · · g n | g i ∈ γ} . (2.4) In this way we construct a compact orientable surface with constant negative curvature. Note that this procedure will not always lead to a smooth surface. For example in the case of the {8, 8} tessellation we obtain the genus-2 hyperbolic Riemann surface, i.e. the regular double torus, but the surface obtained from the {4, 5} tessellation must necessarily have cusps in order to be consistent with the Gauss-Bonnet theorem.

Lattice Hamiltonian
The general structure of the Hamiltonians we consider is Here x and y label sites on the square lattice, and the sum is over nearest neighbors. s x and s y are the 'spins,' in other words s x ∈ M for each point x on the lattice. For the XY model M = S 1 and the standard approach is to replace it with a real variable θ x , with θ x and θ x + 2π identified. The function h(s x , s y ) is a map M × M → R which is bounded below, usually with its lower bound attained precisely when s x = s y . For the XY model, the standard choice is h(θ x , θ y ) = 1 − cos(θ x − θ y ) . A convenient choice for some purposes is the so-called 'Villain approximation' [5], where the sum over n enforces 2π periodicity of θ x and θ y . The natural generalization of eq. (2.7) to the quotient with s x and s y points on H 2 ⊂ R 1,2 . The right-hand side is non-negative because any two points on the hyperboloid H 2 are spacelike separated. An important feature of the Villain energy function (2.7) is that it is continuous but only piecewise smooth: there is a discontinuity in its first derivative along the locus where θ x − θ y ≡ π mod 2π. Likewise, the generalization (2.8) is continuous but only piecewise smooth: for example, if s y is at the origin of the fundamental octagon, then h(s x , s y ) has discontinuities in one of its first derivatives at the boundaries of all images of that octagon. An intrinsic coordinate system with periodically defined coordinates is not known, so it is non-trivial to give an explicit, smooth map analogous to (2.6). We will therefore work strictly with the 'Villain' form (2.8), which we may equivalently define as where the dot product is in the standard mostly plus flat metric on R 1,2 . The equivalent forms (2.8) and (2.9) are not suited to computation unless we can efficiently restrict the minimization to a small subset of the elements in Γ. For the XY model, this is easy to do: one requires θ x , θ y ∈ (−π, π), and then the only images one needs to consider are θ x and θ x ± 2π. We must ask: if O ⊂ H 2 is the fundamental octagon, and JHEP06(2016)145 we require s x , s y ∈ O, then what is the analogous subset of images γ(s x ) that we must minimize over to be sure of finding the minimum value of −γ(s x ) · s y ?
A sufficiently large subset of Γ for the regular double torus is the identity element together with elements γ such that γ(O) touches the fundamental octagon either along a side or at a corner. This subset, call it Γ 49 , has 49 such elements, which can be constructed as follows. A standard basis for Γ satisfies the identity 7 i=0 γ i = 1. The 49 group elements of interest are 1 together with where j < k. Counting the distinct γ jk is straightforward if we consider how they move us along the dual graph to the octagonal tiling of H 2 , a subgraph of which is shown in figure 2. A simple topological way to define Γ 49 is that it is the minimal set of generators such that an open set S ⊂ H 2 can be found satisfying A simulation of the torus model amounts to designing a Markov chain process which results in random sampling of configurations from a probability distribution proportional to e −H/T , where H is given by (2.5) and h is given equivalently by (2.8) or (2.9). To build such a process, one must be able to choose spins s x in the fundamental octagon with uniform probability with respect to the natural measure inherited from H 2 and one must be able to evaluate all instances of h(s x , s y ), which in practice is done by restricting the minimization in (2.8) or (2.9) to γ ∈ Γ 49 .

Simulation results
For temperatures T slightly smaller than the critical temperature where the specific heat is maximized (figure 3), we notice a surprising clustering of spins around one of the six Weierstrass points (figure 4). 1 Infrared fluctuations must eventually cause the system to explore all possible regions of phase space, so for sufficiently large lattices we would expect to see a distribution of spins which is democratic among the Weierstrass points at any fixed temperature. Thus the interesting point is that a spatial correlation length is large enough near the critical temperature so that essentially our whole lattice clusters in the vicinity of one Weierstrass point. Above the critical temperature a topological phase transition occurs as a result of of vortex proliferation. Unlike the XY-model where the winding of a loop of spins is characterized by an integer, the winding of spins in the double torus model is characterized by an operator, which can be determined in the following way: for any consecutive pair of spins s 1 and s 2 in the loop we find the image of s 2 that is closest to s 1 . This amounts to finding the operator O s,t ∈ Γ 49 that minimizes the distance s · (Ot). If the loop consists of spins s 1 , . . . , s n , the operator corresponding to the loop is We map this operator to a number by taking the matrix trace. The fact that the winding is characterized by a non-abelian matrix rather than a number suggests that we should consider loops which are larger than the elementary 1 × 1 plaquette. If we define vortices with respect to a single plaquette, then only finitely many 1 We observe only one pronounced local maximum for the specific heat. Thus we only have evidence for one critical temperature where a phase transition may occur. It is interesting, however, that there are many types of vortices, and it is not impossible that there could be many phase transitions as a result. possible group elements can appear in O loop but we may see a bigger set of group elements if we consider larger loops. Indeed, this intuition is corroborated by simulations (see figure 5) in which we see that the number of vortices of fixed trace increases as we increase the loop size from 1 × 1 to 3 × 3. The logarithmic scale on these plots indicates that vortices are proliferating approximately according to Boltzmann statistics; that is,

JHEP06(2016)145
for some N 0 and E. Moreover, the vortices appear for all temperatures, not just for T > T c . A closer examination of the occurrence of vortices is worthwhile. We refer to figure 5 and table 1. The main features to note are: • The exponential law (3.2) persists up to temperatures comparable to T c . At larger temperatures, we see some evidence of saturation, where the number of vortices per site becomes greater than 1/10, and continuing to follow (3.2) would eventually conflict with the limit N V /L 2 < 1. Provisionally then, we think of (3.2) as a dilute gas approximation.  • The prefactor N 0 , which we think of as related to the exponential of the fugacity, is conspicuously different for some trace classes from 1 × 1 plaquettes to 2 × 2 to 3 × 3. When N 0 rises drastically with plaquette size, as it does for trace 340.7 vortices from 1 × 1 to 2 × 2 plaquettes, we presume that a contributing effect is that the characteristic size of the vortex is too big for the smaller plaquette size to have a significant chance of capturing the full vorticity. Monotonic increase in N 0 is to be expected as plaquette size increases, because if the vortex is small, a large plaquette will enclose it completely at several different positions, all of which contribute to N 0 .
We suspect that a low-temperature expansion in terms of a dilute gas of vortices of many types can be used to account for much of the dynamics up to T c . It would be interesting to pursue this further because it could be a low-dimensional analog to the hadron gas treatment of the low-temperature phase of QCD, inspired originally by Hagedorn's statistical bootstrap approach [6]. Optimistically, one might try to estimate T c in terms of the growth in the number of different vortices with energy. The number of different vortex types increases exponentially with the length on H 2 between an initial point s and its image O loop s. (This is essentially the statement that the area enclosed by a circle in H 2 increases exponentially with its radius.) So it is not implausible that the number of different vortex types also increases exponentially with energy, facilitating a Hagedorn-style argument where a "vortex gas" eventually reaches a maximum possible temperature, which is T c . However, with our present numerical results, we cannot go very far in developing such a vortex gas model, for two main reasons. as obtained from a fit of N V to the Boltzmann form (3.2) for vortices in a given trace class. Moreover, the standard Kosterlitz-Thouless treatment of the XY model discourages us from thinking that the energy of a single vortex is well-defined in isolation. Second, the coefficients N 0 are not really accessible from numerics. Ideally they should correlate with the number of vortices of a given energy. We have suggested a Hagedorn treatment of vortices because the proliferation of vortex types with the length of the corresponding cycle on H 2 is reminiscent of the exponential growth of hadron species with mass. This viewpoint leads naturally to the parametrization (3.2), which approximately captures the population of different types of vortices over a significant temperature range. However, a simple alternative point of view is that vortices with small winding always dominate, and that T c occurs naturally when the dilute gas approximation for the smallest vortices breaks down, and the larger vortices have at most a modest effect on the thermodynamics. Referring to the 1 × 1 plaquette results in figure 5, the dilute gas approximation for the smallest vortices does break down for T only modestly below T c .
Finally, it would be interesting to explore the relationship of our work with [12] who study the Kosterlitz-Thouless transition in the context of a two-dimensional dilute bose gase.

Discussion and future directions
We have considered one of an infinite number of tilings of the hyperbolic plane H 2 . The {8, 8} is simplest tiling which yields a smooth target manifold but it is straightforward to generalize to hyperbolic Riemann surfaces of higher genera. It will also be interesting to consider quotients of hyperbolic space H n with n > 2. In addition, it would be interesting to consider deformations in the moduli space of the double torus which correspond to changing the geodesic lengths of the various cycles. Under such conditions one can expect a splitting of the phase transition into different critical temperatures corresponding to each cycle.

JHEP06(2016)145
A natural question to ask is whether the sigma models studied in this paper have welldefined UV completions. A general nonlinear sigma model is a theory of maps R 2 −→ T with Euclidean action functional given by The beta function for this theory can be computed perturbatively in the coupling parameter α 0 . Assuming that T is an n-dimensional maximally symmetric space for simplicity, Here µ and ν are curved indices for the target manifold T and α is a flat index for the worldsheet R 2 . Taking T = H n , the one-loop beta function is given by β(α ) = (n − 1)α 2 , which naively suggests a Landau pole in the UV. It is conceivable, however, that the higher derivative terms can balance the one-loop contribution and drive the theory to a conformal fixed point [7,8]. There also exist arguments [9] that models defined on compact targets with negative curvature H n /Γ should possess conformal fixed points which arise from the competition between infrared freedom at weak coupling and the discrete spectrum of the Laplacian at strong coupling. It is therefore an important open problem to verify if the transition seen in simulations is of the second-order type, and to ascertain if critical fluctuations are described by a conformal field theory. Answering this delicate question may require moving beyond the Villain approximation, by finding an exact embedding of hyperbolic double torus.
In this paper we have focused on the models with constant negative curvature and compact target H n /Γ but it is also interesting to consider models on the non-compact space H n . The non-compactness of the target makes it difficult to define correlation functions and a non-standard basis of observables may be required [10].
Finally, the modern conformal bootstrap techniques have so far only been applied to theories with compact symmetry groups. It would be very interesting to extend them to the non-compact groups.

JHEP06(2016)145
A Coordinate systems The following coordinate system covers the upper sheet of H 2 in R 1,2 , X 0 = cosh ρ , (A.1) The induced metric and area element are given by To sample points uniformly from the Poincaré disk we note that the form of the area element implies that we should choose φ uniformly in [0, 2π] and choose the cumulative distribution function v to be uniformly distributed in In practice, to sample points uniformly from the fundamental octagon we draw points uniformly from the Poincaré disk by picking φ and v uniformly and then rejecting points that do not lie within the octagon. Defining x = sinh ρ cos φ and y = sinh ρ sin φ we obtain the alternative parametrization The Poincaré disc model is obtained by defining r = tanh(ρ/2). Then (A. 10) and the relationship with the rectangular coordinates is

B Parametrization of hyperbolic polygons
The solutions to the geodesic equation on the Poincaré disc are given by [11] x(s) = cos φ cosh s − R sin φ sinh s √ 1 + R 2 cosh s + R ,

JHEP06(2016)145
where s ∈ (−∞, +∞). These functions define arcs on the Poincaré disc with the radius R and the centers at the point x 0 = √ 1 + R 2 cos φ, y 0 = √ 1 + R 2 sin φ, lying beyond the unit disc. In order to draw a the fundamental polygon corresponding to the {p, q} tessellation, we need to fix the angle α between any two neighboring lines from the center of the disk to the vertices to be α = 2π/p and the distance a of all the vertices from the center of the disk to be such that the interior angles are 2π/q. We will assume the corners of the polygon are at the points a exp (ikπ/2), a exp i(α + kπ/2). The lines connecting these points are given by the solutions to the geodesic solutions specified above. The p geodesics corresponding to the p exterior circles are completely specified by radii R ± and angles φ ± + kπ/2, given by [11] R ± = 1 2a T 2 ± + (1 − a 2 ) 2 , φ ± = arctan T ± 1 + a 2

C Simulation details
Our simulation uses a Metropolis-Hastings algorithm to produce equilibrium configurations of spins on an L × L periodic lattice with the energy given by the Villain approximation in eq. (2.9). Our simulation proceeds in steps wherein at each step a spin, s x at site x, is selected uniformly from the L 2 spins on the lattice. A new spin, s x , is then drawn independently and at random from the {8, 8} tessellation using the procedure established in section A. The spin s x is replaced by s x with probability, P s→s = min 1, exp − 1 T <x,y> (h(s x , s y ) − h(s x , s y )) (C. 1) where h is defined by eq. (2.9). The results presented in this paper are for simulations of L = 50 and L = 100 lattices. To ensure that the systems have reached equilibrium we run the simulation for t eq = 2 × 10 9 steps before beginning to take data. For the L = 50 lattices this amounts to t eq /L 2 = 8 × 10 5 steps per lattice site while for the L = 100 system this is t eq /L 2 = 2 × 10 5 steps per lattice site. Since t eq /L 2 establishes the relevant timescale for the system, the fact that we see no finite size effects in either the energy or the heat capacity gives us confidence that the systems are well equilibrated. Note, that at temperatures lower than those presented here we do indeed see such finite size effects.
After having equilibrated the system we output the state of the system every ∆t = 10 6 steps which amounts to ∆t/L 2 = 400 steps per site for the L = 50 system. We assume that ∆t is sufficiently long such that the system is allowed to uniformly sample equilibrium configurations. Again, the lack of finite size effects between the L = 50 and L = 100 systems supports this claim. We then compute the energy by averaging over 200 such states and the heat capacity by calculating C V = T −2 δE 2 where the latter quantity is the variance of the energy.