Towards an UV fixed point in CDT gravity

CDT is an attempt to formulate a non-perturbative lattice theory of quantum gravity. We describe the phase diagram and analyse the phase transition between phase B and phase C (which is the analogue of the de Sitter phase observed for the spherical spatial topology). This transition is accessible to ordinary Monte Carlo simulations when the topology of space is toroidal. We find that the transition is most likely first order, but with unusual properties. The end points of the transition line are candidates for second order phase transition points where an UV continuum limit might exist.


Introduction
Since the middle of last century physicists have been pursuing the idea of unifying the four fundamental interactions, the strong, the weak, the electromagnetic and the gravitational interactions. The framework of Quantum Field Theory (QFT) unified the first three of them in the so-called Standard Model. Including gravity remains an unsolved problem in a QFT context 1 . Difficulties appear when one tries to formulate a quantum version of Einstein's theory of General Relativity. The naive quantization leads to a perturbatively non-renormalizable theory which cannot be simply included in the unified model of all interactions. The idea of asymptotic safety introduced by Weinberg [1] is an attempt to formulate a non-perturbative QFT of gravity. It assumes that the renormalization group flow in the bare coupling constant space leads to a non-trivial finite-dimensional ultraviolet fixed point around which a new perturbative expansion can be constructed which leads to a predictive quantum theory of gravity. The so-called Exact Renormalization Group program [2][3][4][5][6] has tried to establish the existence of such a fixed point with a fair amount of success, but relies in the end, despite the name, on truncation of the renormalization group equations. Thus it would be reassuring if other non-perturbative QFT approaches could confirm the exact renormalization group results.
Lattice QFT is such a non-perturbative framework and it is well suited to deal precisely with the situation where one identifies fixed points, since these are where one wants to reach continuum physics by scaling the lattice spacing to zero in a way which keeps physics fixed. It has been very successful providing us with results for QCD which are not accessible via perturbation theory. There exists a number of lattice QFT of gravity. One of them, the so-called Dynamical Triangulation (DT) formalism [7][8][9][10][11][12] has provided us with a "proof of concept", in the sense that it has shown us, in the case of two-dimensional quantum gravity [13][14][15][16], that the continuum limit of the lattice theory of gravity coupled to conformal field theories agree with the corresponding continuum theories. Of course there are no propagating gravitational degrees of freedom in two dimensions, but the main issue with the lattice regularization is whether or not diffeomorphism invariance is recovered when the lattice spacing goes to zero. That is the case in the DT formalism, and for the conformal field theories living on the lattice one obtains precisely the non-trivial critical scaling dimensions obtained also in the continuum, i.e. scaling dimensions which are different from the ones in flat spacetime (the so-called KPZ scaling [17][18][19]). The DT formalism was extended to higher dimensional gravity [20][21][22][23][24][25][26][27], but there it was less successful [28,29]. It is not ruled out that the theory can provide us with a successful version of quantum gravity, but if so the formulation has to be more elaborate than the first models (see [30][31][32][33] for recent attempts). However, there is one modification of DT which seems to work in the sense that lattice theory might have a non-trivial continuum limit, the so-called Causal Dynamical Triangulations model (CDT). The model is more constraint than the DT models because one assumes global hyperbolicity, i.e. the existence of a global time foliation.
The CDT model of four-dimensional quantum gravity is realized by considering piecewise linear simplicial discretizations of space-time. The simplicial building blocks can be glued together, satisfying the basic topological constraints of global hyperbolicity (as mentioned) and a simplicial manifold structure. The quantum model is now defined using the Feynman path integral formalism, summing over all such geometries with a suitable action to be defined below. The spatial Universe with a fixed topology evolves in proper time. Geometric states at a fixed value of the (discrete) time are triangulated, using regular three-dimensional simplices (tetrahedra) glued along triangular faces in all possible ways, consistent with topology. The common length of the edges of spatial links is assumed to be a s . Tetrahedra are the bases of four-dimensional {4, 1} and {1, 4} simplices with four vertices at time t connected by time links to a vertex at t ± 1. All time edges are assumed to have a universal length a t . To construct a four-dimensional manifold one needs two additional types of four-simplices: {3, 2} and {2, 3} (having three vertices at time t and two vertices at t ± 1). The structure described above permits for every configuration the analytic continuation between imaginary a t (Lorentzian signature) and real a t (Euclidean signature). Even after Wick rotation the orientation of the time axis is remembered. The spatial and time links may have a different length, and are related by αa 2 s = a 2 t . The quantum amplitude between the initial and final geometric states separated by the integer time T is a weighted sum over all simplicial manifolds connecting the two states. In the Lorentzian formulation the weight is assumed to be given by a discretized version of the Hilbert-Einstein action. (1.1) where [g] denotes an equivalent class of metrics and D M [g] is the integration measure over nonequivalent classes of metrics. A piecewise linear manifold where we have specified the length of links defines a geometry without the need to introduce coordinates. In the CDT approach the integration over equivalent classes of metrics is thus replaced by a summation over all triangulations T satisfying the constraints. After a Wick rotation the amplitude becomes a partition function where S R is a suitable form of the Einstein-Hilbert action on piecewise linear geometries. There exists such an action, which even has a nice geometric interpretation, the so-called Regge action S R for piecewise linear geometries [34]. In our case it becomes very simple because we have only two kinds of four-simplices which we glue together to form our piecewise linear four-manifold: where N 0 is the number of vertices in a triangulation T , N 41 and N 32 are the numbers of {4, 1} plus {1, 4} and {3, 2} plus {2, 3} simplices, respectively. The action is parametrized by a set of three dimensionless bare coupling constants, K 0 , related to the inverse gravitational constant, K 4 -the dimensionless cosmological constant and ∆ -a function of the parameter α, the ratio of the spatial and time edge lengths (for a detailed discussion we refer to [35] and to the most recent review [36] and for the original literature to [37,38]). The amplitude is defined for K 4 > K crit 4 and the limit K 4 → K crit 4 corresponds to a (discrete) infinite volume limit. In this limit, the properties of the model depend on values of the two remaining coupling constants. The model was extensively studied in the case, where the spatial topology was assumed to be spherical (S 3 ) [39][40][41][42][43][44][45]. The model could not be solved analytically and the information about its properties was obtained using Monte Carlo simulations. It was found that the model has a surprisingly rich phase structure, with four different phases. The most interesting among the four phases is phase C, where the model dynamically develops a semiclassical background geometry which in some respect is like (Euclidean) de Sitter geometry, i.e. like the geometry of S 4 . Both the semiclassical volume distribution and fluctuations around this distribution can be interpreted in terms of a minisuperspace model [46][47][48][49]. For increasing K 0 phase C is bounded by a first-order phase transition to phase A, where the time correlation between the consecutive slices is absent. For smaller ∆ phase C has a phase transition to a so-called bifurcation phase, where one observes the appearance of local condensations of geometry around some vertices of the triangulation [50][51][52][53]. The phase transition is in this case of second or higher order. For still lower ∆ the bifurcation phase is linked with the fourth phase, the so-called B phase, where one observes a spontaneous compactification of volume in the time direction, such that effectively all volume condenses in one time slice. The phase transition between the bifurcation phase and the B phase is also of second or higher order [44]. The behavior of the model near continuous phase transitions is crucial if one wants to define a physical large-volume limit (a careful discussion of this can be found in [54]). In this respect phase C stands out, the reason being that only in this phase the large scale structure of the average geometry is "observed" (via the Monte Carlo simulations) to be four-dimensional, isotropic and homogeneous, and one can define an infrared semiclassical limit with a correct scaling of the physical volume [42,46]. Via phase C we thus want a renormalization group flow in the bare coupling constant space towards an UV fixed point (the asymptotic safety fixed point), while keeping physical observables fixed. The natural endpoint of such a flow would be a point in the phase diagram where several phases meet. In the early studies it was speculated that there could be a quadruple point, where all four phases meet. Unfortunately the numerical algorithm used was not efficient in this most physically interesting range in the coupling constant space. As a consequence it was not possible to analyze the model in this range.
The present article discusses a new formulation of the model, where the spatial topology is assumed to be that of a three-torus (T 3 ) [55][56][57], rather than that of a three-sphere, which was the topology used in all the former studies. It was found that the four phases in this case are the same as in the spherical model, with the position of phase boundaries shifted a little 2 . The additional, important bonus in this new formulation comes from the fact that the physically interesting region in the bare coupling constant space mentioned above becomes numerically accessible with the standard algorithm used in the earlier studies. We could then observe that the speculative quadruple point, maybe not surprisingly, separates into two triple points, connected by a phase transition line between phase C and the B phase, and not separated by the bifurcation phase (see Fig. 1). An important point is that we now have access to these triple points directly from phase C and it is thus possible to have a renormalization group flow from the infrared to the potential UV fixed point entirely in the "physical" C phase. The phases of the model were identified for a system withN 41 = 160k, analyzing the structure of geometry at the grid of points in the coupling constant plane shown in Fig. 1, the different phases represented by dots with different colors. In the presented phase diagram the precise position of phase transitions was not determined. This requires a careful study of the infinite volume limit and scaling of the position of phase transition lines with the lattice volume. The most interesting region is the one separating phase C and B where we may observe two triple points. The present paper is the first step in the analysis of this most physically interesting region. We will perform a detailed analysis of the behavior of the model at K 0 = 4.0 in the neighborhood of the phase transition line. We will try to determine the order of the phase transition at this point. We will show that the transition seems to be a first order transition. The results presented in this article show that the most interesting region in the bare parameter space can successfully be analyzed using the standard Monte Carlo algorithm used in the earlier simulations.
2 The phase structure of CDT As mentioned, the phase diagram of the CDT model with a toroidal spatial topology permits us to investigate the properties of the model in an important range of the bare coupling constants, previously inaccessible to numerical measurements. For systems with a spherical spatial topology a detailed analysis of the phase diagram was performed following two lines in the bare coupling constant space. These were the vertical line with varying ∆ at K 0 = 2.2 and the horizontal line at ∆ = 0.6. In the first case it was possible to analyze the phase transition between C and bifurcation phases and between the bifurcation and B phases. In the second case a transition between the C and A phases was studied (see [58] for recent results). The belief coming from the analysis of the spherical case was that if we decrease the value of ∆ for a fixed value of K 0 we necessarily move from C phase to the bifurcation phase and only, for still lower ∆, to the B phase. However, changing to toroidal spatial topology we discovered that this is not the case, probably also in the spherical topology. There exists a range of bare coupling constants where C and B phases are directly neighboring. This happens close to the ∆ = 0 line in the range of K 0 between, approximately, 3.5 and 4.5. One may expect the existence of two triple points (instead of the previously conjectured quadruple point): one triple point where C, A and B phases meet, and a second triple point where C, bifurcation and B phases meet. Finding the precise location of the triple points may be numerically more difficult than analyzing the generic transition between phase C and B. As a first step in the detailed analysis we have chosen to determine the position and the order of the phase transition between C and B phases along a vertical line at K 0 = 4.0. This is approximately in the middle between the position of the two triple points. Since the characteristic behavior in the two phases corresponds to different symmetries of the configurations (we have translational symmetry in time in the C phase and a spontaneously breaking of this symmetry in the B phase) we expect a relatively large hysteresis when we cross the phase boundary. We want to find methods which make the hysteresis effect as small as possible. We also expect relatively large finite size effects. An important point in the analysis will be to check how the hysteresis behaves when the system size goes to infinity.
The analysis presented in the paper is based on a study of systems with a fixed time period T = 4 and different (almost) fixed volumes N 41 . In the earlier studies, it was shown that reducing the period T does not produce significant finite-size effects [58]. On the other hand, in particular in the C phase, the average volume per time slice for a fixed total volume gets relatively large, which is very important. In the Monte Carlo simulations we enforce the lattice volume N 41 to fluctuate around a chosen valueN 41 , so that the measured N 41 =N 41 . This is realized by adding to the Regge action (1.3) a volume-fixing term (2.1) In the thermalization process it is essential to fine-tune the value of K 4 in such a way that one gets stability of the system volume. This is realized by letting the value of K 4 dynamically change by small steps, until the required stable situation is realized. If a value of K 4 is too high, we observe that system volume stabilizes below the target valueN 41 . Similarly, if we take it too small, the volume will be too large. Only for K 4 ≈ K crit 4 (N 41 ) fluctuations of volume are centered aroundN 41 with the width controlled by . During the thermalization part of the Monte Carlo simulations the algorithm tries to find the optimal value of K 4 for a given fixed set of parameters K 0 , ∆ andN 41 . The whole process of measurements is organized in the following way: • We start a sequence of thermalization runs at a set of ∆ values in the neighborhood of the expected position of the phase transition. The initial configuration of the system is taken to be the small hyper-cubic configuration discussed in reference [55]. We choose the target volumeN 41 and let the system size grow towardsN 41 and adapt the K 4 value from the guessed initial value. The initial K 4 can be chosen either a little below or a little above the guessed critical value.
• We find that on the grid of ∆ values we can determine ranges corresponding to the appearance of two different phases, with a relatively sudden jump between the phases. In general the jump is observed between two neighboring values on the grid of ∆.
The corresponding values of K 4 are markedly different in the two phases. Typically the value is smaller for the C phase than for the B phase. We can determine the phase of the system by the measured values of the order parameters (see later for definitions), which are very different in the different phases.
• The value of ∆ where the phase transition is observed depends on the initial value of K 4 used in the thermalization process. As a consequence, we observe in general two values ∆ crit low (N 41 ) and ∆ crit high (N 41 ). Both values are determined with the accuracy depending on the grid of ∆.
• We repeat the analysis on a finer grid, which covers the range where we observed phase transitions. We found the most effective procedure is to restart the Monte Carlo evolution from the same small initial configuration as before, but using as the initial values of K 4 the ones determined for the C or the B phase from earlier runs in the neighborhood of the transitions, corresponding to ∆ crit low (N 41 ) or ∆ crit high (N 41 ) respectively.
• A finer grid permits to determine the two positions of the phase transition with better accuracy. The different position of jumps between the two phases (low or high) can be interpreted as the hysteresis effect in a process where we slowly increase the value of the ∆ parameter or slowly decrease its value. We observe that the size of the hysteresis for a particular choice ofN 41 does not decrease within reasonable thermalization times. By taking a finer grid in ∆ we can only determine the end points of a hysteresis curve with a better accuracy. We illustrate the situation in Fig. 2. The lines shown were obtained from the measured values of ∆ and K 4 for N 41 = 160k.
• In the range of ∆ values between ∆ crit low (N 41 ) and ∆ crit high (N 41 ), depending on the initial value of K 4 a system ends either in the B or C phase. This can be interpreted as a range of parameters, where the two phases may coexist. The distribution of the values of the order parameters (to be defined below), characteristic for the two phases, is very narrow. As a consequence, a tunnelling between the two phases is never observed after we have reached a "stable" ensemble of configurations in the thermalization stage.  The thermalization path chosen above means in practice, that in the beginning, the system grows in a relatively random way from the initially small configuration to the desired target volumeN 41 and then the geometry evolves to a stable range in the configuration space. The first step can be interpreted as a step in the direction typical for the phase A, where correlations between the spatial configurations in the consecutive time slices are small or absent. Only afterwards we reach the domains corresponding to the two phases we study. As a consequence, we expect that the described method will be very well suited to the future analysis of the triple point involving the A phase.
The behavior of the pseudo-critical values K crit 4 (N 41 ) is very similar to that of ∆ crit (N 41 ). This can be seen in Fig. 3 Figure 3. The pseudo-critical value K crit 4 (N 41 ) as a function of ∆ crit (N 41 ). The data points measured for increasing lattice volumeN 41 are going from left to right. Center of the black ellipse corresponds to the estimated position of (∆ crit (∞) , K crit 4 (∞)) and its radii corresponds to the estimated uncertainties. Colors of the fits follow the convention used in Fig. 2.
Although the size of the hysteresis shrinks with volumeN 41 , the plots indicate that the shrinking process is relatively slow and thus in order to get rid of the hysteresis one should use extremely large lattice volumes, not tractable numerically. The dependence of ∆ crit on the lattice volume, ranging betweenN 41 = 40k andN 41 = 1600k is presented in Fig. 4. As it was explained above, the plot contains four sets of data corresponding to the four different points describing the hysteresis (see Fig. 2). The data points can be fitted with the curve The best fit for the combined sets of data (with fixed ∆ crit (∞) = 0.022 determined above) was obtained for γ = 1.64 ± 0.18. An alternative fit with γ = 1 (and the same value of  where the measured value of γ = 1.62 ± 0.25 agrees well with the result obtained for ∆ crit . The fits are represented by curves with different colors, which again follow the convention used in Fig. 2. On the scale used in this plot the green and blue curves practically overlap.

Order parameters
To identify the phases of CDT with toroidal spatial topology we follow methods used in the previous studies. These are based on the analysis of order parameters which have a different behavior in the different phases. We use order parameters which characterize both global and local properties of the simplicial manifolds. The global order parameters were called O 1 and O 2 , where In each phase the distributions of N 0 and N 32 are very narrow, and practically Gaussian.
Phases B and C are characterized by very different average values for the two distributions.  The dependence of the order parameters O 1 and O 2 onN 41 at the endpoints of the hysteresis is presented in Fig. 6. The colors follow the convention used in Fig. 2.
The data presented on the plots correspond for eachN 41 to the four values of the ∆ crit (N 41 ) points, following again the notation of Fig. 2. It is seen that although both pseudo-critical values K crit 4 (N 41 ) and ∆ crit (N 41 ) become very close for increasingN 41 , this is not the case for the order parameters, which in fact behave in a way similar to that characterizing the first order transition. It means that a transition between the B and C phases becomes very rapid. On the other hand, due to the observed hysteresis, the method used in this analysis chooses a position of measured values for the order parameters slightly away from the true transition point (located inside the hysteresis region) and thus in fact we were not able to perform stable simulations exactly at K crit 4 (N 41 ) and ∆ crit (N 41 ) corresponding to such a transition point 3 .
A similar behavior is observed for the set of local order parameters O 3 and O 4 defined by Here n t is the number of tetrahedra shared by {4, 1} and {1, 4} four-simplices with bases at time t and t n t = t  This is accompanied by the appearance of two singular vertices located at times t ± 1 and linked to a macroscopic number of four-simplices in a triangulation. should approach zero. The behavior of these two order parameters is presented in Fig. 7.

Conclusion and Discussion
In the present article we made a detailed study of the phase transition observed between the phase C and the phase B at the value of the dimensionless gravitational coupling constant K 0 = 4.0. The transition appears to be located close to ∆ = 0. The identification of this region, and the possibility that one can move all the way to the triple points of the phase diagram, staying entirely inside the "physical" C phase, is a good news for the renormalization group program started in [54] (and temporarily put on hold by the discovery of the bifurcation phase). The renormalization group analysis is probably the cleanest way to connect CDT lattice gravity approach to asymptotic safety. The analysis of the relevant coupling constant region was made possible by switching from spherical spatial topology to toroidal spatial topology. In this first study of the interesting region we positioned ourselves in the middle of the B-C phase transition line, between the two triple endpoints and from the analysis of the Monte Carlo data we conclude that the transition is most likely of first order. Since endpoints of phase transition lines often are of higher order, the triple points might well be of second order and one of them could then serve as a UV fixed point for a quantum theory of gravity. We are actively pursuing this line of research. Let us end by some remarks about our quantum gravity model, viewed as a statistical system of four-dimensional geometries. Despite the almost trivial action (1.3), the model has an amazingly rich phase structure, with four different phases, each characterized by very different dominating geometries. In addition, some of the phase transitions have quite unusual characteristics. The transition between phase B and the bifurcation phase is a second order transition [44], but superficially, for a finite volume, it looked like a first order transition. However, analyzing the behavior as a function of the increasing lattice volume the first order nature faded away. Moving towards larger values of K 0 , i.e. towards the region we have been investigating in this article, the transition became more and more like a first order transition. With the spherical spatial topology used in [44] one could not get to the region investigated in the present article, but it is natural to conjecture that passing the triple point moving from the bifurcation-B line to the C-B line, the transition changes from second order to first order. However, this first order transition is still somewhat unusual. Firstly, it has kept the characteristics of the second order bifurcation-B transition that the finite size behavior of the pseudo-critical points, given by eqs. (2.2) and (2.3) have non-trivial exponents γ. Secondly, the hysteresis gap goes to zero with increasing volume, which is a non-standard behavior in the case of a first order transition. However, the jumps of the order parameters seem volume independent and that is the main reason that we classify the transition as being a first order transition. The large finite size effects we observe might be related to the global changes of dominant configurations which take place between phase C and phase B, and these global rearrangements might, for finite volumes, have a different "phase-space" in the case of spherical and toroidal topologies. That might explain why our Monte Carlo algorithm can access the B-C transition only in the case of toroidal topology. The statistical theory of geometries is a fascinating area which is almost unexplored for spacetime dimensions larger than two.

Acknowledgements
DN would like to thank Renate Roll the fruitful discussions and hospitality during his stay at Radboud University in Nijmegen. JGS acknowledges support from the grant UMO-2016/23/ST2/00289 from the National Science centre, Poland. JA acknowledges support from the Danish Research Council grant Quantum Geometry, grant 7014-00066B.