{'fraction_non_alphanumeric': 0.029417981961668546, 'fraction_numerical': 0.01252466178128523, 'mean_word_length': 4.940043946845244, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 9, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 2, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Quantizing Horava-Lifshitz Gravity via Causal Dynamical Triangulations

We extend the discrete Regge action of causal dynamical triangulations to include discrete versions of the curvature squared terms appearing in the continuum action of (2+1)-dimensional projectable Horava-Lifshitz gravity. Focusing on an ensemble of spacetimes whose spacelike hypersurfaces are 2-spheres, we employ Markov chain Monte Carlo simulations to study the path integral defined by this extended discrete action. We demonstrate the existence of known and novel macroscopic phases of spacetime geometry, and we present preliminary evidence for the consistency of these phases with solutions to the equations of motion of classical Horava-Lifshitz gravity. Apparently, the phase diagram contains a phase transition between a time-dependent de Sitter-like phase and a time-independent phase. We speculate that this phase transition may be understood in terms of deconfinement of the global gravitational Hamiltonian integrated over a spatial 2-sphere.


Connections
An intriguing body of evidence hinting at a deep connection between Hořava-Lifshitz gravity 1 and causal dynamical triangulations has recently accumulated in the literature. First, one of us demonstrated consistency of the spectral dimension computed in Hořava-Lifshitz gravity with the spectral dimension measured in causal dynamical triangulations in 3 + 1 dimensions [32]. Benedetti et al then verified within causal dynamical triangulations [19] the prediction of Hořava-Lifshitz gravity for the behavior of the spectral dimension in 2 + 1 dimensions [32]. Next, Ambjørn et al exhibited the remarkable resemblance of the phase diagram of causal dynamical triangulations [2] to the phase diagram of Lifshitz matter systems [29] as exemplified by the Lifshitz scalar field [34]. These authors and one of us also conjectured [2,34] that the apparent tricritical point of the former phase diagram could correspond to the tricritical limit of Hořava-Lifshitz gravity with dynamical critical exponent z equal to the dimension of space. Then, both Benedetti et al and Ambjørn et al noted the compatibility of certain solutions of Hořava-Lifshitz gravity with the minisuperspace model fit to the expectation value of the geometries emerging from causal dynamical triangulations [2,19]. Recently, Sotiriou et al successfully fit the behavior of the spectral dimension of causal dynamical triangulations at intermediate scales to a dispersion relation derived from Hořava-Lifshitz gravity [48]. Just days ago, Budd argued that the kinetic term in the semiclassical effective action of causal dynamical triangulations exhibits a Hořava-Lifshitz-like form [22].
This mounting evidence motivated us to extend the Regge action-or, more precisely, the discrete path integral measure-used in causal dynamical triangulations to include dependence on the broader class of terms appearing in the action of Hořava-Lifshitz gravity. Our purpose is threefold: firstly, to test the applicability of causal dynamical triangulations to modified classical theories of gravitation; secondly, to explore quantum Hořava-Lifshitz gravity with nonperturbative techniques; and, thirdly, to illuminate further the links between Hořava-Lifshitz gravity and causal dynamical triangulations. As an initial step toward these goals, we have begun to investigate an appropriate reduction of (2 + 1)-dimensional projectable Hořava-Lifshitz gravity using causal dynamical triangulations. This model, though removed from the phenomenologically interesting case of 3 + 1 dimensions, serves as a simplified yet nontrivial beginning for our research.
After briefly introducing the formalisms of Hořava-Lifshitz gravity and causal dynamical triangulations in section 2, we explain our adaptation of causal dynamical triangulations to Hořava-Lifshitz gravity in section 3. We present the results of our initial numerical studies-evidence for the existence of extended phases of geometry, the structure of these phases, and their relation to the classical solutions-in section 4. Finally, we summarize our conclusions and discuss ongoing and future work in section 5. To streamline our presentation, we relegate to appendix A the construction of the relevant classical solutions of Hořava-Lifshitz gravity and to appendix B certain geometric properties of causal dynamical triangulations.

Hořava-Lifshitz Gravity
Hořava-Lifshitz gravity is a field theory of the dynamical metric on spacetime manifolds that carry the additional structure of a preferred foliation by spacelike hypersurfaces. The relevant set of gauge symmetries is then the group Diff F (M) of diffeomorphisms of the spacetime manifold M that preserve the preferred foliation F . The foliation structure leads to one important novelty: the possibility of anisotropic scaling with a nontrivial dynamical scaling exponent z measuring the degree of anisotropy between space and time. In Minkowski spacetime the relativistic scaling is thus replaced by the anisotropic scaling for constant b > 0. Starting with z > 1 at short distances markedly improves the ultraviolet behavior of the theory, potentially rendering it power-counting renormalizable.
One might be tempted to rewrite this theory in a manifestly relativistic fashion by integrating in the remaining part of the group Diff(M) of full spacetime diffeomorphisms. This procedure leads to a reformulation of the model as a specific scalar-tensor theory with higher-derivative interactions in which the spacelike hypersurfaces of constant scalar field dynamically determine the leaves of the foliation F . Such a relativistic rewriting is only equivalent to the original nonrelativistic formulation at the classical level with, moreover, subtle regularity conditions on the scalar field's dynamics. At the quantum level the relativistic rewriting immediately besets the theory with the notorious problem of time, whereas the nonrelativistic nature of the original formulation of Hořava-Lifshitz gravity potentially renders the spacetime metric's dynamics more directly compatible with quantum mechanics. Additionally, as we shall show, the original formulation of Hořava-Lifshitz gravity is nicely suited to the framework of causal dynamical triangulations, which also utilizes a preferred foliation structure in the microscopic definition of the path integral for gravity.
In a local smooth coordinate chart (t, x) adapted to the preferred foliation F , Diff F (M) consists of all reparametrizations of the form for arbitrary functions f and ζ. Note that, although reparametrizations of the space coordinates may be time-dependent, reparametrizations of the time coordinate must be space-independent. Given its preferred foliation, the structure of Hořava-Lifshitz gravity is naturally discussed in the Arnowitt-Deser-Misner formalism [18]. In this formalism the spacetime metric tensor g(t, x) is decomposed into the metric tensor γ(t, x) on a spacelike hypersurface Σ of constant time coordinate t, the shift vector N(t, x), and the lapse function N (t, x) such that one can reassemble the standard line element as Note, however, that recombining the spatial metric tensor, the shift vector, and the lapse function into the spacetime line element (2.3) is somewhat misleading since, in the regime with z > 1, different terms contributing to ds 2 carry different scaling dimensions. The shift vector and the lapse function play the role of gauge fields associated with Diff F (M). Since the time coordinate reparametrizations (2.2a) are independent of the space coordinates, a natural choice is to restrict the corresponding gauge field N (t, x) also to be only a function of the time coordinate. We choose to make this restriction, yielding the so-called projectable version of the theory. (See, for instance, [31,20,46,52,53].) Of course, we have interest in studying the more general nonprojectable version in which the lapse function is a spacetime field and of which the projectable version is a dynamical limit. The gauge symmetries then permit new terms in the action involving spatial derivatives of the lapse function; finding the appropriate realization of such terms in the framework of causal dynamical triangulations is beyond the present work's scope. We instead focus on the projectable version of Hořava-Lifshitz gravity, which has a particularly clear translation into the language of causal dynamical triangulations.
Employing the Arnowitt-Deser-Misner decomposition of the metric tensor, we now construct the action of projectable Hořava-Lifshitz gravity. We aim to build a power-counting renormalizable unitary classical theory of gravitation. Even in relativistic gravity, including higher curvature terms can render the theory renormalizable [26,50]; typically, however, such terms come at the cost of sacrificing perturbative unitarity and propagating unphysical degrees of freedom [49]. These issues stem from the inclusion of higher temporal derivative terms, which necessarily accompany higher spatial derivative terms in a relativistic theory. By only permitting higher spatial derivatives, we can in principle avoid problems with perturbative unitarity and unphysical degrees of freedom. Working in the Arnowitt-Deser-Misner formalism allows for the rather straightforward inclusion of higher spatial derivatives and exclusion of higher temporal derivatives as we desire in constructing the action of Hořava-Lifshitz gravity.
We build the general action in d spatial dimensions as the sum of a kinetic term quadratic in temporal derivatives and a potential term of mass dimension 2d in spatial derivatives. The dynamical critical exponent z thus equals d. Accordingly, we form the kinetic term from invariants of ∂ t γ(t, x). This derivative alone is not covariant under Diff F (M), but the extrinsic curvature tensor K(t, x) of a spacelike hypersurface Σ satisfies this criterion. With K(t, x) having components for the covariant derivative ∇ associated with γ(t, x), the most general kinetic term invariant under Diff F (M) Here the parameter λ arises from the generalized DeWitt supermetric compatible with the theory's gauge symmetries, and K(t, x) is the trace of the extrinsic curvature tensor K(t, x) [30,31]. We form the potential term from invariants of γ(t, x) and its spatial derivatives. The most general potential term is x)] is a scalar functional of the spatial metric tensor and its spatial derivatives up to order 2d. Putting together (2.5) and (2.6), the action of projectable Hořava-Lifshitz gravity is Note that the coupling constant G is related but not equal to the Newton constant G N : the low energy Newton constant G N is determined by rescaling the time coordinate by an effective speed of light factor selected so that the term linear in the spatial Ricci scalar within V [γ(t, x)] is correctly normalized with respect to the kinetic term [30,31]. This rescaling also dictates how the constant term in V [γ(t, x)] relates to the low energy cosmological constant Λ. For d > 2, with the anticipated z = d scaling at short distances, there is a proliferation of marginal and relevant contributions to the potential term; for d = 2, however, the most general action containing only the marginal and relevant terms for z = 2 is quite compact: with R 2 (t, x) the Ricci scalar of the spatial metric tensor γ(t, x) and the coupling constant Λ related to the low energy cosmological constant Λ by the rescaling described above. The potential term's relative simplicity in (2 + 1)-dimensional projectable Hořava-Lifshitz gravity further motivates our initially studying this case. The equations of motion (A.2) derived from variation of the action (2.8) with respect to the spatial metric tensor are independent of the coupling constant β: for d = 2 the R 2 term is a total derivative, and the Gauss-Bonnet theorem determines its integral over a spacelike hypersurface Σ of the preferred foliation in terms of the Euler number of Σ. Similarly, the equations of motion (A.1) obtained from variation of the action (2.8) with respect to the shift vector are also insensitive to β. The coupling constant β does, however, appear in the equation of motion derived from variation of the action (2.8) with respect to the lapse function. Since the lapse function only depends on the time coordinate, its equation of motion takes the form of a spatially integrated Hamiltonian constraint for the global Hamiltonian At this stage there is an additional choice to make beyond that of restricting to projectable Hořava-Lifshitz gravity. Recall the distinction that arises between noncompact and compact spatial topology when imposing the local Hamiltonian constraint in relativistic theories of gravity. In the former case, after imposing appropriate asymptotic fall-off conditions and accounting for possible boundary contributions, the zero mode H ⊥ of the local Hamiltonian constraint does not vanish on physical states. Instead, this mode is one of the conserved asymptotic charges, namely the total energy. In the latter case the zero mode H ⊥ of the local Hamiltonian constraint vanishes on all physical states since the relativistic gauge symmetries mix it together with the other gauge symmetry generators. Now, in Hořava-Lifshitz gravity, as pointed out in [33,34,54], this situation is somewhat different: while the case of noncompact spatial topology remains essentially unchanged, a novelty arises in the case of compact spatial topology. Imposing the vanishing of H ⊥ on physical states is not mandatory: either one may treat H ⊥ as a gauge symmetry generator and thus impose the condition H ⊥ = 0, or one may treat H ⊥ as a conserved global charge measuring the energy levels of physical states in the Hilbert space. Of course, the second option requires us to check that this Hamiltonian's spectrum is bounded from below on physical states. This novel situation was explicitly encountered in [54] whose authors constructed from bosonic systems on a rigid lattice the renormalization group fixed points corresponding to the free field limit of Hořava-Lifshitz gravity with z = 2 and z = 3. There the specific choices that one could make at the level of the microscopic lattice Hamiltonian gave rise to the option of imposing or not imposing the condition H ⊥ = 0 for the effective low energy gravitons.
When we discretize the action (2.8) in section 3 for use as a functional on simplicial manifolds, we find that the second option-treating H ⊥ as a conserved global charge-appears more natural in the setting of causal dynamical triangulations. The geometric restrictions derived by employing simplices of fixed edge lengths suggest that the lapse function has been effectively set to a constant, though without leading in any discernible way to the imposition of the constraint (2.9). Also, the R 2 term with coupling constant β in the action (2.8) does not affect any observables for the topological reasons stated above. The value of H ⊥ does, however, depend on β, which suggests that enforcing the constraint (2.9) is not completely consistent. In section 4.3, where we compare the geometries emerging from our numerical simulations to the classical solutions for the constraint (2.9) not imposed, we find evidence supporting our selection of the second option. These comparisons are, moreover, not without weight: as briefly discussed in [34] and further demonstrated in [35], the classical phase diagram for the theory with variable lapse differs significantly from that of the present setting. The classical theory that we quantize in section 4 is thus defined by the action (2.8) for fixed lapse.

Causal Dynamical Triangulations
The causal dynamical triangulations approach aspires to define the continuum limit of a quantum theory of gravity by appealing solely to those nonperturbative tools applied with great success to the quantization of local gauge field theories. (See [38,39] for reviews.) In particular, the approach invokes a lattice regularization of the path integral for gravity and then utilizes the renormalization group to explore its continuum limit. In causal dynamical triangulations one thus attempts to define this path integral as the partition function of a statistical model of dynamical geometry. (See [3,16] for recent reviews.) As in several previous programs of lattice quantization of gravity, the statistical ensemble of geometries is comprised of simplicial manifolds triangulated by a fixed set of (d + 1)-simplices. In 2 + 1 dimensions this set consists of those 3-simplices or tetrahedra pictured in figure 1.
The novel feature of causal dynamical triangulations lies in the imposition of an additional restriction on the geometries permitted to enter the path integral: these simplicial manifolds must possess a global foliation by spacelike hypersurfaces of constant discrete time. The global foliation was originally introduced to enable a Wick rotation from Lorentzian to Riemannian signature. The causal structure of a Lorentzian triangulation could thus be faithfully retained in its Wick rotated version, only the latter triangulation being amenable to numerical analysis [7,8,10]. A simplicial manifold is endowed with the foliation structure as follows: every spacelike hypersurface, all of a chosen fixed topology, is triangulated with equilateral d-simplices, and then the vertices of adjacent spacelike hypersurfaces are connected by timelike edges so as to produce only the (d + 1)-simplices of the fixed set. Spacelike edges have length squared l 2 SL = a 2 , defining a lattice spacing for the triangulation, and timelike edges have length squared l 2 T L = −ηa 2 for η > 0 [7,8]. Note that our parameter η corresponds to the parameter α typically used in the causal dynamical triangulations literature.
More concretely, one aims to approximate the path integral for general relativity, by the path sum In (2.11) the path integral is taken over all physically distinct spacetime metric tensors g(t, x) interpolating between the initial and final boundary geometries specified by the spatial metric tensors γ(t i , x) and γ(t f , x), and S EH denotes the Einstein-Hilbert action In (2.12) the path sum is taken over all causal triangulations T interpolating between the initial and final boundary geometries specified by the triangulations Γ i and Γ f [7,8], and S R denotes the Regge action Here, h is a (d + 1 − 2)-dimensional hinge having area A h and deficit angle δ h , and V s is the spacetime volume of a (d + 1)-simplex [51]. The measure factor 1 C T is the inverse of the order of the automorphism group of T , included to account for discrete symmetries arising in the structure of T .
Before proceeding further with our discussion of causal dynamical triangulations, we make a few clarifying remarks on the implications of the above prescription for regularizing the continuum path integral (2.11). The restriction to simplicial manifolds admitting a preferred foliation effectively changes the path integral measure and integration domain. Even though we continue to employ the Regge action, these changes may result in our model belonging to a universality class different from that of Euclidean dynamical triangulations. In 1+1 dimensions, where the path sum (2.12) can be evaluated analytically, this is known to be the case [1]: by preventing the birth of baby universes, the foliation requirement alters the critical exponents, indicating that Euclidean and causal dynamical triangulations occupy different universality classes. In higher dimensions the same situation seems very likely to hold true; otherwise, we should expect our model to fall into one of the universality classes in which no smooth macroscopic limit exists, and the spacetime geometry exhibits a branched polymer or crumpled behavior.
The change in the path integral measure and integration domain may also translate into a change in our theory's effective action such that it no longer assumes the general relativistic form. Recall that, although we may locally treat the choice of a preferred foliation as a gauge choice, we cannot generally make such a choice globally. A gauge choice typically introduces a Jacobian into the path integral measure, and a mismatch between the path integral measure and the gauge choice can lead to a breaking of the gauge symmetry. Starting with the Regge action, the question of whether counterterms sensitive to the preferred foliation are generated thus remains open. If such counterterms appear, then the resulting effective action likely corresponds to Hořava-Lifshitz gravity at some particular values of its couplings. The Regge action's bare coupling constants G N and Λ are only indirectly related to the continuum renormalized values of the Newton constant and the cosmological constant at long distances, so we must not naively identify the former with the latter.
Continuing with our discussion of causal dynamical triangulations, the action (2.14) simplifies considerably for the simplicial manifolds contributing to our ensemble since only a small set of simplices with fixed geometries is used to construct them. Still, the path sum (2.12) has resisted all attempts at analytical computation for d > 1, so exploration of its properties has primarily employed numerical techniques. In particular, Markov chain Monte Carlo methods are used to simulate the path sum (2.12) [4,5,9,10,11,12,14,15,24,40]. To render the path sum (2.12) amenable to such analysis, we must first Wick rotate the real time action (2.14) to imaginary time, a well-defined process owing to the foliated structure of causal dynamical triangulations. We perform the Wick rotation by analytically continuing the parameter η in the lower half complex plane [8]. For our case of interest-2 + 1 dimensions with spacelike hypersurfaces having the topology of S 2 and periodic time coordinate having the topology of the 1-sphere S 1 -the Regge action becomes after Wick rotation and application of various topological relations for η = 1 [9]. Here, N 0 is the number of vertices and N 3 is the number 3-simplices in the triangulation T . The coupling constants k 0 and k 3 are related to the bare couplings G N and Λ as We have thus transformed the path sum (2.12) into the statistical partition function Computer simulations of the partition function (2.17) in both 2 + 1 and 3 + 1 dimensions have thus far provided considerable support for the existence of an extended phase of geometry possessing not only a semiclassical limit on large scales, but also a quantum regime on small scales [4,9,10,11,12,13,14,15,19,24,40]. In particular, the average observed geometry matches well that of (possibly deformed) Euclidean de Sitter spacetime at both the classical and semiclassical levels [4,9,10,11,12,13,14,15,19,24]. Furthermore, studies of the spectral dimension of this extended phase of geometry have revealed an apparent dimensional reduction to effective 2-dimensionality on small scales [13,19,40]. This phenomena of dynamical dimensional reduction-particularly, the extrapolated value of the minimal dimensionality-has elicited comparisons of causal dynamical triangulations to both the asymptotic safety approach and Hořava-Lifshitz gravity. In the former theory the effective change in the spectral dimension apparently results from large anomalous dimensions near the nontrivial fixed point, while in the latter theory the spectral dimension flows to the short distance value of 2 as a result of the anisotropic scaling near the Gaussian fixed point [32,42]. Further evidence for such dimensional reduction has also surfaced in other approaches to the quantization of general relativity [23].
As in Euclidean dynamical triangulations, additional phases of geometry also emerge [43]. At first these phases were viewed as unphysical, but current interpretations favor their role as further phases in the vicinity of a multicritical fixed point [2]. Specifically, in both 2 + 1 and 3 + 1 dimensions there exists a phase characterized by spacelike hypersurfaces that effectively decouple from one another [9,12,40], and in 3 + 1 dimensions a second additional phase, distinguished by its large Hausdorff dimension, appears [12,40]. Recently, Ambjørn et al have argued that the transition between this last phase and the physical phase is of second order [5]. This finding raises the possibility of rigorously defining a continuum limit of causal dynamical triangulations.

A Discrete Action for Hořava-Lifshitz Gravity
We now derive a discrete form of the action (2.8) suitable for analysis with the techniques of causal dynamical triangulations. As above we assume a topology of S 2 × S 1 , primarily motivated by the relative ease of numerically analyzing such compact spacetimes. As in the lattice regularizations of local quantum field theories [38,39], the precise details of the discretization do not matter since universality ensures that many of the details at the lattice spacing scale become irrelevant in the long distance limit. Our primary goal in constructing a discrete analogue of the action (2.8) is thus to build an action sufficiently specific so that the continuum limit lies in the universality class of Hořava-Lifshitz gravity yet sufficiently generic so that the continuum limit does not depend on all of the renormalized coupling constants in this universality class. Of course, we can only judge whether or not we have achieved these goals once we have thoroughly studied the quantum theory of the model defined below.
Working from the continuum action (2.8), we allow two technical criteria to guide our construction of its discrete version: first, the discrete action should manifestly reduce to the Regge action used in causal dynamical triangulations when the bare coupling constants λ and α assume their general relativistic values, and, second, the transfer matrix corresponding to the discrete action defined on the space of boundary geometries should yield a well-defined Hamiltonian. These two criteria apply to the classical discrete action constructed below. In the quantum theory defined via the path integral based on this action, we naturally anticipate a nontrivial relationship between the bare coupling constants and the renormalized coupling constants. Depending on this relationship, the quantum theory's long distance continuum limit may or may not be relativistic.
To implement the first criterion, we use the Gauss-Codazzi equation, to rewrite the action (2.8) as up to irrelevant boundary terms. We neglect the term in the action (2.8) with coupling constant β as it contributes only an additive constant for the spacetime manifolds under consideration. In this form the action (3.2) straightforwardly reduces to the Einstein-Hilbert action (2.13) when the coupling constants λ and α take on their general relativistic values of one and zero. We use the discrete action (2.15) for the Einstein-Hilbert portion of the action (3.2), leaving us the task of determining the discrete analogues of the K 2 and R 2 2 terms. We construct these terms below, postponing discussion of the second criterion.

Curvature Squared Terms in Regge Calculus
There exist well-established prescriptions for constructing the Ricci scalar and the trace of the extrinsic curvature tensor in Regge calculus. The former is defined in terms of deficit angles δ h about (d + 1 − 2)dimensional hinges h of a (d + 1)-dimensional simplicial manifold [51]: The latter is defined in terms of angles ψ h between the normal vectors to the two d-simplices intersecting at the (d + 1 − 2)-dimensional hinge h within the spacelike hypersurface Σ of a (d + 1)-dimensional simplicial manifold [21,28]: Technically, these curvatures have a distributional definition on the hinges, which complicates the construction of curvature squared terms. In particular, simply taking the square of (3.3) or (3.4) to define the discrete versions of the R 2 2 or K 2 terms leads to a mathematically ill-defined continuum limit [6]. Accordingly, we adhere to the philosophy of [6,27] when building discrete analogues of curvature squared terms, adopting the alternative scheme of volume sharing. As discussed in [6], we make the identification for the number N (v) of spacelike triangles containing v, the volume is and the share-volume is since all of the spacelike triangles are equilateral. Employing (3.5) and noting that the most natural discretization of the time integral is we make the identification up to multiplicative factors, where V τ (T ) denotes the set of vertices v belonging to the spacelike hypersurface labeled by the discrete time coordinate τ .

K 2 Term
For the K 2 term the choices of objects o and top-dimensional objects containing o are not as clear. The extrinsic curvature captures in part how the spacelike hypersurface is embedded in the spacetime manifold. The discrete analogue of the K 2 term must reflect this geometric information, requiring that we appropriately account for how tetrahedra connect to the spacelike hypersurface. This observation suggests that we take tetrahedra as the top-dimensional objects contributing to the share-volume. Furthermore, in the continuum K 2 scales as an inverse length squared, implying that the objects o be spacelike triangles. We thus need to assign a deficit angle δ to a spacelike triangle. We largely follow the treatment of [37]. Consider a spacelike hypersurface of the triangulation T . For each spacelike triangle we may define a future-directed normal vector at the center of the (3, 1)-tetrahedron of which it forms the base. A natural measure of extrinsic curvature at the common edge e between two adjacent such spacelike triangles is the deficit of the angle that the normal vector traces out as it is parallel transported from its own (3, 1)-tetrahedron to the adjacent (3, 1)-tetrahedron. In figure 2 we illustrate this construction for the particular case of three (2, 2)-tetrahedra stacked between the two adjacent spacelike triangles. In general, the deficit angle for the common edge is B. We thus assign to a spacelike triangle the deficit angle which, we note, is reminiscent of a trace. The volume of a spacelike triangle is in appendix B. The share volume assumes this value since a given spacelike triangle has four (3, 1)tetrahedra and N ↑ (2,2) ( ) (2, 2)-tetrahedra in its immediate future. According to the prescription (3.5), we make the identification up to multiplicative factors, where T SL τ (T ) denotes the set of spacelike triangles belonging to the spacelike hypersurface labelled by the discrete time coordinate τ .
The discretization (3.15) of the K 2 term does not, however, respect our second criterion. Following [8], to ensure the existence of a well-defined Hamiltonian on the space of boundary geometries, we must make (3.15) time-reversal invariant. A straightforward calculation shows that this invariance guarantees positivity of the squared transfer matrix, which, along with the transfer matrix's symmetry, yields a well-defined Hamiltonian [8]. To realize time-reversal invariance, we add an analogous term for past-directed (2, 2)-tetrahedra. The complete discrete analogue of the K 2 term is Note that we did not require such an adjustment for the R 2 2 term since it depends only on the intrinsic geometry of the spacelike hypersurface.

Imaginary Time Action
Putting together (3.10) and (3.16), our discrete action for Hořava-Lifshitz gravity becomes . (3.17) Wick rotating the action (3.17) to imaginary time, we find that We give the values of θ in appendix B. We use the action (3.18) for a = 1 and η = 1 in the path integral quantization analyzed below. The value of the length a has no a priori physical meaning, and any value of the parameter η > 1 2 is permitted.

HL [T ]
(4.1) employing Markov chain Monte Carlo methods. In our simulations so far, we have fixed the topology of the spacelike hypersurfaces to be that of S 2 and have fixed the total number T of spacelike hypersurfaces, introducing a discrete time coordinate τ that enumerates these hypersurfaces. In the simulations reported below, we have set T = 64. Additionally, we impose periodic boundary conditions on this time coordinate, endowing it with the topology of S 1 . Consequently, the initial triangulated spacelike hypersurface Γ i and the final triangulated spacelike hypersurface Γ f of each causal triangulations T entering the partition function (4.1) are identified. Next, we hold the number N 3 of tetrahedra in each triangulation T approximately fixed; otherwise, the weights e −S (E) HL [T ] appearing in the partition function (4.1) may grow without bound, eventually causing our computer code to crash. To implement the constraint of fixed number of tetrahedra, we add to the action (3.18) a term of the form |N 3 −N 3 |, where is a Lagrange multiplier andN 3 is the target number of tetrahedra. Essentially, this Lagrange multiplier term modifies the value of the coupling constant k 3 . In all of the simulations reported below, we have set = 0.02 andN 3 to approximately 10200. These parameter values typically yield a one percent variation in N 3 over any given ensemble. This number of tetrahedra is sufficiently large for physical effects clearly to outweigh finite size effects and sufficiently small for our limited computing resources to survey a reasonable portion of the coupling constant space.
Finally, we tune to a set of coupling constants {k 0 , λ, α, k c 3 (k 0 , λ, α)} on the critical surface of the coupling constant space defined by the partition function (4.1) for a fixed N 3 . As our notation suggests, we first select values for k 0 , λ, and α and then tune to the associated critical value k c 3 of k 3 . The critical value k c 3 is that for which N 3 remains approximately constant while a simulation runs. In this sense our model is only welldefined at the critical value: for any other value the number of tetrahedra either increases without bound or plummets to zero. With these conditions established, a simulation begins with the generation of an initial triangulation having the topology S 2 × S 1 composed ofN 3 tetrahedra. Using the Pachner moves adapted to causal dynamical triangulations, as described, for instance, in [8], we run a standard Metropolis algorithm to generate an ensemble of spacetimes representative of the weighting defined by the partition function (4.1), sampling only after a period of thermalization. We sample the representative spacetimes generated every one hundred sweeps, a single sweep comprisingN 3 attempted Pachner moves. Once collected, we estimate the expectation values of observables as averages over the ensemble.
Testing our code is a nontrivial matter: without a known nonperturbative quantization of (2 + 1)dimensional projectable Hořava-Lifshitz gravity, we possess no standard of comparison for our results. Of course, this situation also pertains to causal dynamical triangulations formulated with the Regge action. Our code is a modified version of that reported in [40], which has yielded independent corroboration of the results of [9,10,11,12,13,14,15,19]. We have run the modified code at the general relativistic values of the coupling constants λ and α to check that we correctly reproduce the results of these references. In figure 3 we present depictions of two representative spacetimes-one in the physical phase and one in the decoupled phase of causal dynamical triangulations-generated by our code. In figure 4 we plot the ensemble average spectral dimensions-discussed further below-for the two ensembles to which the representative spacetimes in figure  3 belong. Up to finite size effects currently under investigation, these measurements agree quantitatively with those of [19,24,40]. In figure 8(a) below we also show the ensemble average discrete 2-volume as a function of discrete time for the second of these two ensembles. This data, as well as the fit to it, are consistent with the findings of [9,24]. Beyond these checks, we rely on the plausibility of our new results as a test of our code's correctness. on a 3-dimensional subspace. This critical surface approximates the so-called infinite volume surface, which corresponds to the limit in which N 3 increases without bound while the lattice spacing a vanishes such that the product N 3 a 3 remains constant. Supposing that our model possesses a second order phase transition, at which we could define its continuum limit, this transition must be located at a phase boundary on the critical surface. Now, a 3-dimensional subspace of largely unknown extent represents a formidably expansive space to explore numerically. Accordingly, we have limited our initial investigations to the subspace of the coupling constant space consisting of the λ − α plane near the origin for fixed k 0 . Specifically, we set k 0 = 1.00, select values for λ and α, and then tune to the value of k 3 on the critical surface. With λ = 1.00 and α = 0.00 this corresponds to a point in the physical phase of causal dynamical triangulations for the Regge action.
Within this region we have generated forty-seven ensembles of representative spacetimes, each for a different set {k 0 , λ, α, k c 3 (k 0 , λ, α)} of the coupling constants. We display in figure 5 the critical surface and the associated phase structure as ascertained thus far. Our explorations indicate the existence of three phases: a phase contiguous with the physical phase of causal dynamical triangulations for the Regge action that we call phase C; a phase emerging for sufficiently large values of λ that we call phase D; and a phase emerging for sufficiently large values of α that we call phase E. We are not entirely certain that phases D and E are distinct: our measurements might instead indicate modulation of a single phase across the relevant region in the λ − α plane, a possibility under active investigation. On the other hand, we are quite certain that phases D and E are not artifacts: our simulations exhibit the characteristic lengthening of the autocorrelation time near the phase boundaries, and the geometric properties of phases D and E persist under an increase in the total number of spacelike hypersurfaces. The presence of these novel phases of extended geometry counts as the first explicit indication that the spatial curvature squared terms in the action (3.18) exert a significant effect as opposed to being renormalized to irrelevance. This nicely matches the scaling behavior of the curvature squared terms expected from the analytic approach to Hořava-Lifshitz gravity.
In figure 6 we depict spacetimes representative of each of the three phases. Specifically, we plot the number N SL 2 of spacelike triangles-a measure of the discrete 2-volume-as a function of discrete time τ . Representative spacetimes in phase C are characterized by a single correlated accumulation of tetrahedra spread across a substantial range of τ . Representative spacetimes in phase D are characterized by an intermittent series of accumulations of tetrahedra each spread across a small range of τ . Representative spacetimes in the phase E are characterized by a moderately uniform distribution of tetrahedra spread across the entire range of τ .

Physicality and Semiclassicality
We now present preliminary evidence suggesting that the ensemble average geometries in all three phases show signs of being both physical and semiclassical. This evidence stems from analyses of two geometric observables associated with our ensembles of representative spacetimes: the spectral dimension as a function of diffusion time and the discrete 2-volume as a function of discrete time. Based on these analyses, we draw certain comparisons to the relevant solutions of the classical theory now in imaginary time presented in

Spectral Dimension
The spectral dimension of a space measures its effective dimensionality as experienced by a random walker. We determine the ensemble average spectral dimension by the method described, for instance, in [19]. Specifically, we directly measure the return probability P r (σ) as a function of diffusion time σ for each representative spacetime in a given ensemble. We then compute the ensemble average spectral dimension as employing an appropriate discretization of the derivative. In figure 7 we display plots of the ensemble average spectral dimension as a function of diffusion time for each of the three phases. Before interpreting these plots, we should comment on the evident bifurcation for small values of diffusion time. This effect reflects the discrete nature of our triangulated spacetimes: random walks of even and of odd lengths yield diverging estimates for the spectral dimension on scales sufficiently short in comparison to the discretization scale. Increasing the number of tetrahedra comprising each spacetime pushes the bifurcation scale towards smaller diffusion times, revealing the physical nature of the spectral dimension on such scales. Now compare the plots of figure 7 to those of ordinary causal dynamical triangulations in figure 4. All three of the Hořava-Lifshitz spectral dimensions much more closely resemble that of figure 4(b)-for the physical phase of causal dynamical triangulations-than that of figure 4(a)-for the decoupled phase of causal dynamical triangulations. The resemblance between the plots of figures 7(a) and 4(b) is particularly close, as we might expect, since the two ensembles are both in phase C. In particular, the spectral dimensions shown in figure 7 all peak at values between 2 and 3 for small diffusion times and then gradually decay for large diffusion times. We expect this decay at large diffusion times: these spacetimes are compact and have curvature.
For (2 + 1)-dimensional causal dynamical triangulations the spectral dimension in the physical phase reaches 3 for ensembles characterized by larger values of N 3 [19,40]. Such measurements provide a key piece of evidence for the semiclassical nature of the phase's ensemble average geometry. This suggests that the maxima of the spectral dimension plots in figures 4(b) and 7(a) are depressed by finite size effects. Assuming that such effects are comparable for our ensembles in phases D and E, we are led to conclude that the maxima of the spectral dimension in these phases is slightly below the topological value of 3. In Hořava-Lifshitz gravity the spectral dimension is generally predicted to flow under renormalization, so potentially the plots of figure 7 reflect such behavior [32]. We are currently running simulations in phases D and E for larger values of N 3 in the hope of resolving these issues. Nevertheless, we maintain that the similarities in form between the spectral dimensions depicted in figure 7 and in figure 4(b) hint at the physicality and semiclassicality of phases C, D, and E.

Discrete 2-Volume
The foliated structure of causal triangulations allows for the measurement of certain quantities as functions of the discrete time coordinate. The number N SL 2 of spacelike triangles on a given spacelike hypersurface is one such quantity. Although the value of N SL 2 on a single spacelike hypersurface is not physically meaningful, the set {N SL 2 (Σ τ )} T τ =1 does contain physical information. For an ensemble in phase C, there exist straightforward techniques for performing a coherent ensemble [15].) These techniques rely on the characteristic feature of phase C spacetimes: the single accumulation of tetrahedra as depicted in figures 3(b) and 6(a). Intuitively, each method functions to align an appropriately defined center of this accumulation with the central value of the discrete time coordinate. Once accomplished for all of the representative spacetimes in an ensemble, the coherent ensemble average of {N SL 2 (Σ τ )} T τ =1 is defined by the discrete timewise average. We employ a method that we call equal discrete 2-volume splitting. The name refers to the algorithm for appropriately defining the center of the accumulation of tetrahedra. Formally, given a representative spacetime, form all divisions D i of that spacetime into two sets each of T 2 spacelike hypersurfaces, maintaining . . .
Then select the two particular divisionsD eq andD eq that most nearly equalize the discrete 2-volume summed over the T 2 spacelike hypersurfaces in each set of the division: the quantity is minimized for bothD eq andD eq . Note thatD eq andD eq only differ in the order of their two sets. Next relabel the discrete time coordinate over the set of values {− T 2 , − T 2 + 1, . . . , T 2 − 1, T 2 } so that the values {− T 2 , . . . , − 1 2 } label the first set of spacelike hypersurfaces and the values { 1 2 , . . . , T 2 } label the second set of spacelike hypersurfaces in bothD eq andD eq . Finally choose the division for which the values of N SL 2 (Σ − 1 2 ) and N SL 2 (Σ 1 2 ) are greatest. With all of the spacetimes in an ensemble so aligned, perform a discrete timewise average of N SL 2 over the representative spacetimes. In figure 8 we show the results of the equal discrete 2-volume splitting average for two ensembles in phase C. The data points indicate the equal discrete 2-volume splitting average value of {N SL 2 (Σ τ )} T τ =1 ; the light vertical bars indicate one standard deviation of error. The thin curve is a one parameter fit of the data points within the central accumulation of tetrahedra to the functional form which is a discretization of the 2-volume as a function of global time of Euclidean de Sitter spacetime. (See [24] and [15] for a derivation of this discretization in 2 + 1 and 3 + 1 dimensions, respectively.) Substantial evidence already demonstrates that, for the general relativistic values of λ and α, the ensemble average geometry in phase C closely approximates Euclidean de Sitter spacetime [4,9,10,11,12,13,14,15,19,24]. The plot in figure 8(b) provides the first evidence that the ensemble average geometry also has this property for the broader range of nonrelativistic values of both λ and α. Hořava-Lifshitz gravity admits Euclidean de Sitter spacetimes: for λ > 1 2 and α = 0, the solutions (A.10) and (A.11) coincide in their description of such spacetimes. Note, however, that for relatively small values of α, the solutions (A.10) and (A.11) do not deviate too markedly from Euclidean de Sitter spacetime. Given the presence of finite size effects and the inherent error in our measurements, one of these solutions for α = 0 may well fit better the ensemble average geometry in phase C. Our current data are not sufficiently detailed for us to make a conclusive statement. Furthermore, there exists some evidence [19] that the ensemble average geometry in (2 + 1)-dimensional causal dynamical triangulations is a deformed version of Euclidean de Sitter spacetime. Supposing that this is the case, the process of taking the continuum limit must generate some additional terms besides those in the Regge action.
The two plots in figure 8 quite closely resemble one another-in both spatial and temporal extents of the central accumulation of discrete 2-volume-even though that of figure 8(a) is for general relativistic values of λ and α while that of figure 8(b) is for nonrelativistic values of λ and α. To dispel the suspicion that the K 2 and R 2 2 terms are renormalized to irrelevance in phase C, we present in figure 9 four plots of the equal discrete 2-volume splitting average of {N SL 2 (Σ τ )} T τ =1 each for a different value of λ. Clearly, varying λ affects the ensemble average geometry that emerges.
For phases D and E we currently do not know how to coherently average {N SL 2 (Σ τ )} T τ =1 over an ensemble. In an effort to determine a method, we computed the ensemble average power in the discrete Fourier transform of {N SL 2 (Σ τ )} T τ =1 , the results of which we display in figure 10. As these plots show, there is no notable periodicity present in these ensemble's average geometry since virtually all of the power falls in the zero frequency mode. This lack of periodicity may in fact be indicative of the semiclassical nature of phase E. To achieve a more faithful comparison, we have considered the above two ratios instead of the deviations themselves. Although not definitive, the fact that the ratios for phase E fall between those for phase C lends credence to our interpretation of the semiclassical nature of phase E. Supposing that further analysis supports the conjecture that time-independent geometries dominate phase E, how could one then interpret the C-E phase transition? A possible answer involves the global gravitational Hamiltonian H ⊥ of (2.10). Recall that, in our implementation of Hořava-Lifshitz gravity into the framework of causal dynamical triangulations, we chose not to impose the condition (2.9) as a constraint, instead treating H ⊥ as a conserved global charge, essentially the total energy. In our statistical ensemble whether H ⊥ vanishes or has a nontrivial spectrum on physical states thus becomes a question of dynamics. If the dynamics enforce the vanishing of H ⊥ on physical states, then we might interpret this phenomenon as dynamical confinement of the gravitational charge. A similar phenomenon has recently been discovered in (2 + 1)-dimensional relativistic chiral gravity [44]. There not just one gravitational charge but an infinite hierarchy of conserved chiral charges are confined, that is, dynamically vanish on all physical states of finite energy.
Now consider the hypothetical phase in which H ⊥ is dynamically confined. On spacetime geometries of the Friedmann-Lemaître-Robertson-Walker type, the Hamiltonian constraint equation becomes the Friedmann equation for the scale factor. The Friedmann equation precludes the ground state geometry from being time-independent, forcing a cosmological, de Sitter-like evolution of the universe. In the context of causal dynamical triangulations, we observe this behavior in phase C.
In the hypothetical deconfined phase, on the other hand, the situation is more reminiscent of a typical condensed matter system: the total Hamiltonian simply measures the system's energy levels with the ground state identified as the lowest energy state, typically static. We apparently observe this behaviour in our phase E. Accordingly, we speculate that the C-E phase transition may be viewed as the deconfinement of the gravitational charge H ⊥ .

Outlook
Motivated by a suite of striking similarities between Hořava-Lifshitz gravity and causal dynamical triangulations, we constructed a discrete version of the action for (2 + 1)-dimensional Hořava-Lifshitz gravity adapted to the formalism of causal dynamical triangulations. Using this action in Markov chain Monte Carlo simulations of the corresponding path integral, we found significant evidence for the existence of extended phases of geometry approximating certain classical solutions of Hořava-Lifshitz gravity. Quantum Hořava-Lifshitz gravity thus appears amenable to and compatible with the techniques of causal dynamical triangulations.
Since we have only just initiated this study, the prospects for further research are expansive. In the near term we have three primary goals for ongoing research. First, we plan to map more extensively our model's coupling constant space to determine the extent of the phases so far discovered and to identify the locations of phases not yet discovered. Specifically, we hope to ascertain how the curvature squared terms affect the decoupled phase of causal dynamical triangulations. Second, we wish to further our analysis of the semiclassical natures of phases C, D, and E along the lines of [4,15,19,24]. Third, we want to establish the orders of our model's phase transitions. In this direction we have attempted to determine order parameters for the C-D and C-E phase transitions; unfortunately, none of our trial order parameters have yet distinguished between these adjacent phases. If our speculation about the confinement-deconfinement nature of the C-E phase transition is correct, then the appropriately defined ground state energy may act as an order parameter.
In the long term the richness of the literature on Hořava-Lifshitz gravity provides a host of directions for continuing research. First of all, recall that the simplest versions of (2 + 1)-dimensional Hořava-Lifshitz gravity, both projectable and nonprojectable, possess a single propagating scalar degree of freedom. This is in stark contrast to the case of not only (2 + 1)-dimensional general relativity, but also the nonrelativistic generally covariant version of Hořava-Lifshitz gravity constructed in [33]. This scalar mode has generated considerable controversy with regard to its phenomenological viability [46,52], but its full dynamics, especially around the stable ground state, remains poorly understood. The (2 + 1)-dimensional theory should provide a clear window into the scalar mode's dynamics since the complications of propagating tensor modes are absent. By studying the nonperturbative dynamics of the scalar mode using causal dynamical triangulations, we hope to shed light on this issue. Our first challenge is the identification of an observable in causal dynamical triangulations that measures the number of local propagating degrees of freedom or that at least distinguishes between the absence and presence of local propagating degrees of freedom. Given the relative ease of extracting correlation functions of observables associated with the spacelike hypersurfaces of causal triangulations, we are currently working to understand the behavior of our model's conformal mode in relation to the scalar mode. We hope that such an investigation might illuminate the nature of the scalar dynamics.
Relatedly, perturbations of projectable Hořava-Lifshitz gravity about (3 + 1)-dimensional Minkowski spacetime generate instabilities [41]. This finding is of course not surprising: Minkowski spacetime is not the ground state of the commonly analyzed models. The situation is much improved on the background of de Sitter spacetime [25,53]. As the coupling constant λ approaches unity, however, the higher derivative terms become relevant, leading to a breakdown of the linearized analysis. The authors of [45,53] have suggested that a nonperturbative Vainshtein mechanism might take effect, rendering this limit continuous to the general relativistic values of the coupling constants. Potentially, we could uncover this transition behavior in Markov chain Monte Carlo simulations. Indeed, the critical surface depicted in figure 5(a) provides a small but intriguing piece of evidence for that possibility: the apparent geometric feature along the λ = 1 direction.
To assess the phenomenological viability of at least the simple version of Hořava-Lifshitz gravity on which we have focused our study, we need to study the renormalization group flow of the coupling constants for the purpose of comparing the long distance behavior to that of general relativity. In [36] Henson proposes a coarse graining procedure for causal dynamical triangulations on which one could try to base a renormalization group procedure. If this scheme proves apt -a question that we are currently exploring -then we plan to employ the procedure to study the renormalization group flows of our model. Ultimately, we would like to explore the general conjecture formulated by several groups [2,34] suggesting that Hořava-Lifshitz gravity and causal dynamical triangulations belong to the same universality class. This circumstance would neatly explain the remarkable resemblances between these two approaches to the quantization of gravity.
for constant of integration C. Letting which has the solution for a second constant of integration δ. In terms of a(t), the solution (A.9) is for positive values of α and for negative values of α. We assume that Λ > 0 as in the previous studies of causal dynamical triangulations. Note that the solution (A.11) has a finite time extent dictated by the zeros of a 2 − (t) whereas the solution (A.10) has no such restrictions. Also note that for vanishing α the solutions (A.10) and (A.11) both reduce to that of Euclidean de Sitter spacetime.

B On the Geometry of Causal Dynamical Triangulations
For the tetrahedra depicted in figure 1, the Lorentzian dihedral angles between spacelike and timelike faces are θ As previously mentioned, Wick rotation consists in analytically continuing η in the lower half complex plane [8]. If the argument of a square root becomes negative as a result of the Wick rotation, then we replace it by the negative of the argument multiplied by −i. The Lorentzian dihedral angles (B.1) are thus continued