The transfer matrix in four-dimensional CDT

The Causal Dynamical Triangulation model of quantum gravity (CDT) has a transfer matrix, relating spatial geometries at adjacent (discrete lattice) times. The transfer matrix uniquely determines the theory. We show that the measurements of the scale factor of the (CDT) universe are well described by an effective transfer matrix where the matrix elements are labeled only by the scale factor. Using computer simulations we determine the effective transfer matrix elements and show how they relate to an effective minisuperspace action at all scales.


Introduction
Minisuperspace models of the universe provide us with simple quantum mechanical models of fluctuations of the scale factor of the universe. In the simplest models one assumes spatial homogeneity and isotropy. Classically this implies that we can write ds 2 = −dt 2 + a 2 (t)dΩ 2 (1) where dΩ 2 is the line element of a homogeneous and isotropic three-dimensional space. Assuming space is compact it is S 3 . Under these assumptions the dynamical variable is the scale factor a(t) of the universe and the quantum field theory of the universe is reduced to quantum mechanics of a single variable a(t).
The 4D Causal Dynamical Triangulation model (CDT) of quantum gravity is by construction a (regularized) quantum field theory model which is compatible with spatial homogeneity and isotropy (for reviews see [1,2]). It uses the path integral formulation and assumes there exists a foliation in (proper) time. When the average over all geometries of this kind is performed one indeed finds that the average geometry, i.e. what one naively would think is closest to a "classical geometry", can be described by a line element of type (1). Even more, it turns out that to a good approximation the scale factor a(t), measured as the average of the third root of the spatial volume at (proper) time t is well described by the simplest minisuperspace action. Contrary to the standard use of minisuperspace, where one postulates the reduction to geometries described by a metric like (1) and then only quantizes the single degree of freedom a(t), the CDT discussion of a(t) is exact (to the extent that CDT describes quantum gravity). The a(t) entering in the CDT discussion is obtained by including all geometries in the path integral. It is natural to ask the following question: how well does the simplest minisuperspace action describe the CDT data generated by Monte Carlo simulations. The tool we will use when trying to answer this question is the transfer matrix.
The time foliation present in CDT provides us with a transfer matrix. The geometries considered in the (regularized) path integral are piecewise linear geometries constructed in such a way that at discretized times t n the spatial slices are triangulations of S 3 . The transfer matrix relates a given spatial (piecewise linear) geometry at time t n to a given spatial geometry at time t n+1 . In ordinary (Euclidean) lattice field theory reflection positiveness of the transfer matrix ensures a unitary time evolution. The CDT transfer matrix has a similar property [3,2] 1 . In this article we will analyze an "effective" transfer matrix of CDT to 1 The idea of a time foliation and the fact that there exists a related unitary time evolution are features which CDT shares with Hořava-Lifshitz gravity [4]. However, no spatial higher derivative terms are explicitly added to the action, like in Hořava-Lifshitz gravity, and it is possible that fixed points of the lattice theory can be identified with the non-trivial UV fixed points conjectured in the asymptotic safety scenario suggested by Weinberg [5] and investigated in [6,7,8].
be defined below.

CDT and the CDT transfer matrix
The use of piecewise linear geometries was introduced in the context of general relativity by Regge [9] as a natural tool to work with a discretized version of the Hilbert-Einstein action, but without the use of coordinates. The idea was to approximate a given smooth geometry by a continuous, piecewise linear geometry. This piecewise linear geometry is uniquely defined by a triangulation where the length of the links are given. The observation by Regge was that the standard Einstein-Hilbert action in D dimensions, for such piecewise linear geometries has a geometric interpretation as a sum over deficit angles of the D − 2-dimensional sub-simplices in the triangulation. While the idea of Regge was to approach a given classically smooth geometry by a sequence of suitable piecewise linear triangulations, a different use of piecewise linear geometries was made in the formalism of Dynamical Triangulations (DT) [10]. Although we have no mathematical rigorous definition of the path integral over geometries, it is natural, in analogy with the path integral in ordinary quantum mechanics, to assume that the summation over geometries will involve not only smooth geometries but all continuous geometries. A subclass of these geometries is the piecewise linear geometries and a further subclass is the piecewise linear geometries defined by triangulations obtained by gluing together equilateral D-simplices such that they form a manifold of fixed topology. In DT the assumption is that this set of geometries is in a suitable sense dense in the set of continuous geometries when we take the link length a to zero, and in this way the link length a will act as an ultraviolet cutoff, just like in ordinary lattice field theory. Further, the natural signature of space-time in the DT formulation is Euclidean. We will assume this is the case in the rest of this article. That this approach works in principle, i.e. that it is able to reproduce a continuum quantum field theory which is diffeomorphism invariant, is well documented for D = 2. Two-dimensional Euclidean quantum gravity coupled to a conformal field theories with c ≤ 1 can be solved analytically both in the continuum [11,12] and using the formalism of DT [13] and agreement is found. In higher dimensions DT was studied numerically, using Monte Carlo simulations both in three dimensions [14] and four dimensions [15]. However, no convincing continuum limit has been obtained so far in higher dimensions [16], and this was one of the motivations for changing the class of triangulations used in the path integral in CDT 2 . In the CDT formalism one sums over geometries with a (proper) time foliation. In principle one starts out with space-times with a Lorentzian signature (contrary to the situation in DT) and the foliation is in proper (Lorentzian) time. However, each piecewise linear geometry used in the CDT path integral allows a rotation to Euclidean proper time. The set of Euclidean geometries we obtain in this way is a subset of the DT Euclidean geometries and this restriction seemingly cures some of the higher dimensional DT diseases, while in two dimensions the relation between the restricted theory and the full DT theory has been worked out in detail: one obtains the CDT theory from the DT theory by integrating out all baby universes (which results in a non-analytic mapping between the coupling constants of the two theories), and (somewhat surprisingly) one can restore the DT theory from the CDT theory by the inverse mapping [17,18]. Using four-simplices (which is the case having our attention in this article) as building blocks one can, for suitable choices of bare coupling constants, observe a four-dimensional (Euclidean) universe [19]. For these choices of coupling constants the shape of the universe is consistent with an interpretation as an (Euclidean) de Sitter space, at least as long as one looks at the scale factor [20]. This is the region of coupling constants which will have our interest (for other choices of the coupling constants one obtains more degenerate configurations [21]).
An interesting feature of the CDT model, rotated to Euclidean signature, is that it possesses a transfer matrix [3,2]. At (discrete) time t n we have a spatial hypersurface with a spatial geometry characterized by a three-dimensional triangulation T 3 (t n ). At time t n+1 we have spatial geometry defined by another three-dimensional triangulation T 3 (t n+1 ). We assume for simplicity that the topology 3 of the spatial triangulations is that of S 3 . In CDT we sum over all four-dimensional triangulations of the "slab" between t n and t n+1 compatible with the topology S 3 × [0, 1] and such that each four-simplex which "fills" the slab has subsimplices which are (sub)simplices of T 3 (t n ) as well as T 3 (t n+1 ). This leads to 4 types of four-simplices in the slab: type (4,1), (3,2), (2,3) and (1,4), where the numbers denote the number of vertices in T 3 (t n ) and T 3 (t n+1 ) respectively. The number N (4,1) (t n ) of (4, 1) simplices in the slab is equal to the number N 3 (t n ) of three-simplices in T 3 (t n ), and similarly the number of (1,4) simplices in the slab is equal to the number N 3 (t n+1 ) of three-simplices in T 3 (t n+1 ). The transfer matrix M, i.e. the amplitude between the T 3 (t n ) and T 3 (t n+1 ), is now given as the sum over all such triangulation, where the summation is over all four-dimensional triangulations of a slab, with boundary triangulations T 3 (t n ) and T 3 (t n+1 ), C T 4 is the order of the automorphism group of the triangulation T 4 , and where S[T 4 ] is the Regge action of the fourdimensional triangulation of the slab. The transfer matrix is defined on the vector space spanned by the set T 3 of three-dimensional triangulations. This space is infinite dimensional and has the natural scalar product where C T is the order of the automorphism group of the triangulation T . The transition amplitude for a three-dimensional triangulation T to develop into a three-dimensional triangulation T after t tot + 1 (integer) time steps is We will define the partition function corresponding to t tot time steps as path integral with periodic boundary conditions after t tot time-steps: This is the partition function we have used in our computer simulations. The measurements performed so far, using Monte Carlo simulations, have been concentrated on the measurement of the scale factor, or more conveniently the threevolume n t i ≡ N 3 (t i ) at the spatial slice at time t i , as well as the correlation between the three-volumes at time t i and time t j . These "observables" can be expressed using the transfer matrix M. The probability of measuring the spatial volume n t i at time t i is given by In (7) T 3 (n t i ) denotes the subset of three-dimensional triangulations where the number of three-simplices is n t i andρ(n t i ) the projection operator on the subspace spanned by these triangulations: Similarly the correlator between n t 1 and n t 2 , separated by ∆t = t 2 −t 1 is given by and this expression can clearly be generalized to multi-correlators. As mentioned in the Introduction it was possible to explain the observed distributions P ttot (n t 1 ) and P ttot (n t 1 , n t 2 ) using a simple minisuperspace model. How does the concept of a minisuperspace labeled by states n t i , i = 1, . . . , t tot relate to the transfer matrix M which is defined on the much larger space spanned by the vectors T ? Let us define the "effective" transfer matrix M by n|M |m represents the average of the matrix elements T n |M|T m like where N n denotes the cardinality of T n . N n grows exponentially with n.
In (10) and (11) it is misleading to think of the "state" |n as the (suitably) normalized sum of the N n vectors |T n (although by doing so one would of course obtain the correct expectation value M n,m using such vectors). Such a vector would again be a single vector located in the N n -dimensional space spanned by the |T n 's. It is more appropriate to think of the "state" associated with n as arising from a uniform probability distribution of states |T n and in this way to think ofρ(n) as the associated density operator. However, once we have reduced our consideration to the matrix n|M |m we are of course free to find eigenvectors for this matrix and expand them in the abstract basis |n , and we will indeed do that.
The statement that we can use the matrix n|M |m as an effective transfer matrix is the statement that the standard deviation of the N n N m numbers T n |M|T m is sufficiently small. In fact the difference between tr M 2 and tr M 2 can exactly be expressed as a sum over deviations squared for each n, m: In the following we will assume that we can work with an effective transfer matrix n|M |m . Eq. (10) is an attempt to define this effective transfer matrix from first principles and in principle one can check by computer simulations if it is a good approximation. We will here take the pragmatic attitude to assume there exists such an effective transfer matrix and use it to analyze the computer generated data. The consistency of this analysis is indirectly evidence that an object like n|M |m provides a good approximation of our data. Thus we will use the "effective" version of (7)- (9): where ρ(n) should be distinguished fromρ(n).
In particular we can measure the matrix elements n|M |m up to a normalization by considering t tot = 2. We have: This method requires a major change in our general computer program which assumes t tot ≥ 3. We updated it but we were not completely convinced that the new version is stable although it gave exactly the same results as the method we finally used. For t tot = 3, 4 we get: From the measurements of P (3) (n 1 , n 2 ) and P (4) (n 1 , n 3 ) we can determine the matrix elements n|M |m up to a normalization: There is nothing magic about the above choice. One could have chosen t tot = 4 and t tot = 6 and formed the combinations from which one can again extract n|M |m like in (19). We have indeed checked that measurements of P (4) (n 1 , n 2 ) and P (6) (n 1 , n 4 ) lead to the same M matrix as extracted from measurements of P (3) (n 1 , n 2 ) and P (4) (n 1 , n 3 ), up to a normalization. In earlier work we have shown [20,22] that the following minisuperspace action 4 : 4 In fact we used slightly different form of the potential terms: µ nt+nt+1 This parametrization was more convenient to extract the parameters of the action from the measured covariance matrix of volume fluctuations. In this article we implement a modified form (22) which better fits our data. describes well the measured n t and the fluctuations n t n t − n t n t in the bulk where n t is large. The effective action (22) suggests that the effective transfer matrix is a good approximation in the bulk. We will in the following try to determine the transfer matrix from the data, also in the range where n t is not necessarily large and we will try to improve the expression (23).

How to perform the computer simulations
The simplest version of the discretized CDT theory has three parameters, two related to the cosmological constant and the gravitational constant, and an additional parameter which controls the asymmetry between the edge lengths in the spatial and time directions. This latter parameter seems not to be a genuine coupling constant since it just labels the different length assignment of spatial and time-like links. The action used is still the Einstein-Hilbert action (as formulated by Regge for piecewise linear geometries), adjusted for this asymmetry. However, because we study the theory in a truly non-perturbative region of coupling constant space the effective action is determined by a competition between the classical action used and a contribution coming from the measure term. The contribution from the measure term is "entropic" in nature: it counts the number of configurations with the same action and is thus independent of the other parameters. Effectively this promotes the asymmetry parameter to a genuine coupling constant (see [2] for a detailed discussion).
In the numerical simulations the topology of the manifold is assumed to be S 3 × S 1 with periodic boundary conditions in the (Euclidean) time, as mentioned above. The four-simplices used to construct the simplicial manifolds of CDT are characterized by their position in spatial and time directions. As also mentioned above we have four types of four-simplices: (4, 1)-simplices, with four vertices at time t and one vertex at t + 1, (3, 2)-simplices with three simplices at t and two at t + 1 and the "time-reversed" (1, 4)-simplices and (2, 3)-simplices. All simplices of a particular type are identical.
The discretized (Regge) Einstein-Hilbert action becomes extremely simple because we are essentially only using the two kinds of building blocks to construct the four-dimensional triangulation T [3,2]: where N 0 is the total number of vertices in the triangulation, N (4,1) the total number of type (4, 1) plus (1,4) simplices and N (3,2) the total number of simplices of type (3,2) plus (2,3). κ 0 , κ 4 and ∆ are the (bare) dimensionless coupling constants obtained by the discretization of the continuous action (2). κ 0 is proportional to the inverse bare gravitational constant, κ 4 related to the cosmological constant while ∆ is related to the asymmetry between the spatial and time-like links. ∆ = 0 corresponds to spatial and time-like links having the same length. An additional geometric parameter is the length t tot of the periodic time axis. The partition function has a critical value κ crit 4 (κ 0 , ∆), depending on κ 0 and ∆, such that Z is divergent for κ 4 < κ crit 4 . The existence of this critical value is reflecting the fact that the number of triangulations with a fixed number of four-simplices N 4 grows exponentially with N 4 . In principle we want to fine tune κ 4 to this critical value since we really want a limit where N 4 → ∞. In practice the simulations have so far been carried out by keeping N 4 (or N (4,1) ) fixed. In this way we have been trading the coupling constant κ 4 with N 4 and the partition function Z(N 4 ) is related to Z(κ 4 ) by a Laplace transformation: The phase diagram now depends on κ 0 and ∆ and we refer to [2,21] for a detailed discussion. Here we will be working in the interesting, so-called de Sitter phase where we, for a given (large) N 4 and sufficient large t tot , observe a (Euclidean) de Sitter universe, i.e. a four-sphere where the temporal extension is proportional to N 1/4 4 while the rest of the time extension (assuming t tot is large enough compared to N 1/4 4 ) is a stalk of almost no spatial extension. Presumably this stalk only exists because our computer algorithm does not allow the spatial extension to shrink to zero. In Fig. 1 we have shown a typical situation with t tot = 80, N (4,1) = 160000 and we observe a bulk region (the "blob", approximately from t = 20 to t = 60) where n t , the three-volume, i.e. the number of tetrahedra at the time-slice t, is large, and the rest is the stalk region where n t is very small. Fig. 1 shows both the average over many configurations and a typical configuration which appears in the path integral. When taking the average over many configuration we align the center of mass of the blobs (see [20] for a detailed discussion).
We can also measure the probability distribution P ttot (n t ) of n t in the blob for a given t in Fig. 1. It is shown in Fig. 2 (left figure). It is well approximated by a Gaussian distribution around the mean value n t . This is in contrast to the situation in the stalk where the probability distribution splits in three families [23], as shown on the right part of Fig. 2.
Figs. 1 and 2 are based on computer simulations of the type mentioned above: N (4,1) is kept fixed. Technically this has been done by adding a term ε(N (4,1) − N (4,1) ) 2 to the action, ε being a suitably small parameter: This term ensures that N (4,1) is going to fluctuate not too far fromN (4,1) . The precise value N (4,1) depends on the choice of κ 4 . We now want to study the transfer matrix. However, the structure of the transfer matrix is incompatible with a global constraint of this type, so we have to change the updating procedure. We have done this in two different ways. The first way is to drop the constraint term ε(N (4,1) −N (4,1) ) 2 and only use the discretized Einstein-Hilbert action (24). The way to obtain an average N (4,1) is to fine tune κ 4 to κ crit 4 . The closer κ 4 is to the critical value the larger N (4,1) . In practice this fine tuning can be difficult and the larger the system the more difficult the fine tuning. Thus we can and will only use it for small systems. For larger systems we apply a different strategy which also constrains the value of N (4,1) , but which is compatible with the transfer matrix structure: we change the global constraint imposed on N (4,1) in (27) to a local constraint in t: Of course this constraint will drastically change the profile n t , since n t will now fluctuate around n vol . Thus we will have different transfer matricesM and M , and different probability distributionsP (n 1 , n 2 , . . .). The new probability distribution for n t is shown for various n vol in Fig. 3. However, we can reconstruct some of the probability distributions associated with the action S R if we know it (i.e. measure it) for the actionS R . The proba- bility for measuring (n 1 , n 2 , . . . , n ttot ) is given bỹ P ttot (n 1 , n 2 , . . . , n ttot ) = n 1 |M |n 2 n 2 |M |n 3 · · · n ttot |M |n 1 trM tt tot , and is directly related to the distribution without the volume fixing term, P ttot (n 1 , n 2 , . . . , n ttot ) ∝ P ttot (n 1 , n 2 , . . . , n ttot )e − (n 1 −n vol ) 2 · · · e − (nt tot −n vol ) 2 . (30) We calculate the transfer matrixM in the same way as M : n|M |m =P (3) (n 1 = n, n 2 = m) P (4) (n 1 = n, n 3 = m) .
To calculate the original transfer matrix M we have to cancel the volume fixing term, which is easily done: From equations (28), (29) and (30) For each choice of n vol we observe n t with some approximate Gaussian distribution centered around n vol , where the width depends on our choice of ε, and we use the associated probabilities to construct n|M |m , as described above. To reconstruct the matrix in a larger region of the n-space we have to merge data from different n vol regions. Since the matrix is determined only up to a normalization, the way to do this is to make sure there are regions of overlap between the n t distributions and in these regions choose a suitable calibration procedure such that we can merge the data. We will later describe how this is explicitly done. Needless to say the procedure we are employing here is a kind of multi-canonical Monte Carlo method (see [24] for a review).

The effective action at large three-volumes
We can measure the transfer matrix n|M |m for large n, m as described above. As already noted it is well approximated by the matrix M (th) : where the effective Lagrangian is We now ask how well? We make a best fit of the parameters Γ, µ, λ and N (fixing n 0 = 0). The measured M , M (th) from (33) as well as their difference, are shown in Fig. 4 for n vol = 1400 and for the n t range 1200 < n t < 1600. The values of parameters Γ, µ and λ for different values of n vol , obtained from the best fits of the matrix M (th) (33) to the measured matrix M are presented in Table 1. Again n 0 is chosen to be zero.

The kinetic term
To get a better estimation of the parameters associated with the effective action (34), we first try to fit only to the parameters of the kinetic term which is by far  Table 1: The values of Γ, µ and λ for different n vol , obtained from best fits of M (th) to the measured M . the dominating term from a numerical point of view. We do that by keeping the sum of the entries, i.e. n + m, fixed such that the potential term is not changing. In this way we can try to determine Γ and even n 0 which we had put to zero in the fits mentioned above in order not to have too many fit-parameters. The matrix elements for constant n + m = c show the expected Gaussian dependence on n (see Fig. 5): where the terms in the effective action which only depend on c are included in the normalization. We expect the denominator of the kinetic term k(c) to behave like k(n + m) = Γ · (n + m − 2n 0 ). As shown on Fig. 6 this is indeed true and the parameter Γ is common for all ranges. Fig. 6 presents measured coefficients k(c) for various c's and ranges of n t denoted by distinct colors together with a linear fit. The best

The potential term
The potential part of the effective Lagrangian may be extracted from the diagonal elements of the transfer matrix L ef f (n, n) = − log n|M |n + c(n vol ) = 1 Γ µn 1/3 − λn .
However, because of different normalizations of the transfer matrices for different ranges (hence the dependence of the constant c(n vol ) on n vol ), the fit of L ef f (n, n) to the transfer matrix data cannot be performed in a straightforward way. The transfer matrices have first to be merged properly via a scaling procedure, i.e. by adjusting the c constant in (36). This is done in the following way. For example, for n vol = 1400 the range of spatial volumes for which we measured the transfer matrix is n t = 1180 . . . 1630, while for n vol = 1800 the range is n t = 1580 . . . 2040. Thus, there is a non-vanishing intersection n t = 1580 . . . 1630 for which elements of both matrices were measured. We scale the second matrix, so that the mean value of the diagonal elements on the intersecting region is equal for both matrices. After applying this procedure for successive ranges, we finally get scaled transfer matrices which can be merged. The result of such merging is shown on Fig. 7,

A global effective action fit
Summing up, the effective action determined via the transfer matrices is strikingly well described by eq. (34). As in the last subsection we can merge the scaled matrices for all the choices of n vol (see Fig. 8), and fit expression (34) to the aggregated data. The best fit gives Γ = 26.17, n 0 = 7, µ = 15.0 and λ = 0.046.
26.17 ± 0.01 7 ± 1 15.0 ± 0.1 0.046 ± 0.001 Previous method* 23 ± 1 − 13.9 ± 0.7 0.027 ± 0.003 Table 2: The values of Γ, n 0 , µ and λ fitted in different ways. * We also present the parameters of the effective action extracted from the covariance matrix of volume fluctuations in our earlier work [23]. We summarize the results of fitting the parameters of L ef f in different ways in Table 2. As long as we are concentrating on the large (bulk) values of n t the various data clearly do not allow us to improve the expression (34). For comparison we also present the parameters of the effective action measured indirectly from the covariance matrix of volume fluctuations. This method is based on the analysis of the effective propagator around the semi-classical solution [23]. The parameters of the action from the previous method agree quite well with those measured directly from the transfer matrix. The small difference may result from slightly different parametrization of the potential term in the effective action (see footnote 4).

Miscellaneous
The spectral decomposition of matrices presents us with an interesting way to compare the measured transfer matrix n|M |m with the "theoretical" matrix n|M (th) |m given by eq. (33) with parameters obtained from a best fit, as described above (in this case the fit from the n vol = 1400 data). We obtain very good agreement between eigenvalues and eigenvectors for the two matrices. The first 6 eigenvalues are presented in Fig. 9 and the first 4 eigenvectors are shown in Fig. 10.
As a final check of the consistency of the effective transfer matrix with data we We useM instead M , because the power M ttot involves summation over all possible volumes n which are not accessible until we suppress them with the term e − (n−n vol ) 2 (which is present inM ). We compare the theoretical expectation (37) with the measuredP ttot (n) for t tot = 3 and t tot = 4. The comparison is shown on Fig. 11 and the error is of order 0.02%.

The transfer matrix for small three-volumes
We have seen that the effective transfer matrix is very well described by the simplest effective Lagrangian (34). However, it is natural to expect that the corresponding effective action is only a first approximation. Appealing to isotropy and homogeneity one would expect that a potential R k (t) term, where R(t) refers to the scalar curvature of the three-dimensional space at time t, translates into a term nt+n t+1 2 1−2k/3 . The leading term R(t) is already part of the effective action where it appears as the term nt+n t+1 2 1/3 . For the data coming from large n t we have seen that the difference between our measured M and M (th) coming from the effective action with only the nt+n t+1 2 1/3 term is already at the noise level.
From these large n t data we have no chance with the present statistics to study higher k corrections. Thus we now turn to the small n t region. The difficulty of analyzing the small n t region is that we might encounter discretization effects as is apparent from Fig. 2. The study of this region started in [23], but the present approach offers the great advantage that we can perform high statistics study of small systems, while in the earlier studies the interesting region was a small part of a larger system and thus the statistics becomes less good. We want to measure the effective transfer matrix by measuring P ttot (n i , n j ) for t tot = 3, 4 as we have already done, but now for small n t , i.e. for much smaller systems. We also want to check as well as possible that the concept of an "effective" transfer matrix actually works. While it seemed to work well for large n t , it is not obvious that this will remain true for small n t : physics might be different and discretization effects might also spoil such a picture. Our systems will be so small that we do not have to introduce the auxiliary transfer matrix M and the volume fixing parameter n vol . Thus we only use the Regge Einstein-Hilbert action S R , eq. (24), to generate the probability distributions P (3) (n 1 , n 2 ) and P (4) (n 1 , n 3 ) and construct n|M |m from eq. (19).

Eigenvectors analysis.
Given the transfer matrix M we can perform a spectral decomposition in terms of eigenvalues λ i and (orthonormal) eigenvectors |α i : Since the measured M is only determined up to a normalization, we will assume λ 1 = 1 and |λ 1 | ≥ |λ 2 | ≥ .... If the gaps are significant between the first, second and third eigenvalues , it is clear that the large t tot limit of Z(t tot ) and the large ∆t limit of the simplest correlation functions P ttot (n t , n t+∆t ) will be completely dominated by the first two eigenstates: Thus, in the limit where t tot ∆t 1 we have (recalling that λ 1 is normalized to 1) ≈ n|α 1 2 α 1 |m 2 + λ ∆t 2 n|α 2 α 2 |m n|α 1 α 1 |m  Figure 11: The probability distributionP ttot (n) (dots) and n|M t tot |n trM t tot (line) for t tot = 3 (left) and t tot = 4 (right).
The average can be written as: For the correlator we have: and the long distance behavior is an exponential fall off e −µ∆t , µ = − log λ 2 /λ 1 (where we have reintroduced λ 1 for clarity). ≈ 0.3222. Note that according to our normalization the biggest eigenvalue is set to one.
How well are these approximate relations satisfied? First we observe that there is indeed a clear gap between the first eigenvalues, as illustrated in Fig. 12.
The closer κ 4 is fine-tuned to κ crit 4 the smaller the gap, but even very close to the critical value we observe a clear gap. That we have gaps even at κ crit 4 just illustrates that we indeed consider small systems.
Next we ask how large t tot has to be in order that the approximations made in eq. (40)-(42) are valid. That can of course be read off from the eigenvalues and already for t tot = 4 the approximation is very good for n ttot . For the correlator one has of course to consider larger t tot . In Fig. 13 we have shown for t tot = 12 the expected exponential decay with exponent log(λ 2 /λ 1 ) compared to the actually measured correlator. The agreement is very good even for small ∆t where it is not obvious that ignoring the eigenvectors with eigenvalues smaller than λ 2 is a valid approximation. ∆t Empirical, T = 12 Empirical full CDT Theoretical Figure 13: The "theoretical" correlator (42) of spatial volumes between timeslices separated by ∆t and calculated including the first two eigenvectors of the transfer matrix M calculated for κ 4 = 0.3223 (green line). The correlator decays exponentially with ∆t (log scale). There is a very good agreement with the correlator measured directly in simulations for T=12 (red) and in "full CDT" stalk range (blue). The bars indicate measurement errors.

The "Full-CDT" approximation.
We have also compared our approximate large t tot probability distributions P (n) (eq. 40) with the data taken from the stalk range of full CDT simulations (including the blob and the stalk range). The distributions approach very well "full CDT" measurements as κ 4 tends to critical value (Fig. 14). The validity of the transfer matrix model is further confirmed by the behaviour of the correlator n t m t+∆t − n t m t+∆t measured directly in the stalk range of "full CDT". As illustrated in Fig. 13 the measured correlator falls off as e −µ∆t . The parameter µ is well explained by the ratio of the first two eigenvalues of the transfer matrix M calculated for κ 4 closest to the critical value.

The effective action for small three-volumes
In principle the matrix M presented above would allow us to determine an effective Lagrangian for small n t , i.e. even in the stalk range of the CDT configurations, precisely as we did for large n t : However, we are confronted with the existence of three families of states, as is apparent in Fig. 2. We can however define a reduced matrixM performing a summation over the three families, i.e.
where the rectangular matrix U has a form: The elements of the matrixM behave much more smoothly and can be analyzed, using the effective action idea. We normalizeM by choosing its largest eigenvalue to be one.
Let us assume that L ef f has the form: where the functions k(·) and v(·) are to be determined. We follow the same strategy as when n t was large: first we analyze matrix elements for constant n + m = c in order to keep the potential term v(n + m) constant. One observes a Gaussian dependence on n (Fig. 15): Fitting (46) for different c's one can easily check that the kinetic coefficient k(c) is linear (Fig. 16). Therefore one can write: where: v() = Γv(). The best fit of Γ and n 0 is presented in Table 3. In the above we recognize the familiar kinetic term present in the effective action for the blob range. Let us also test the assumption that the potential part is similar by analyzing the diagonal elements ofM where the kinetic term is zero. From (43) and (47)  Let us assume v(2n) = −λn + µn 1/3 + δn −ρ , inspired by the bulk expressions already derived and where we have included a new term δn −ρ , inspired by earlier remarks about effective powers of the threedimensional curvature R(t). We present the fit of the logarithm of the diagonal elements n|M |n in Fig. 17 and the parameters of the fit in Table 3. Table 3: Fitted parameters of the effective action for the stalk (50). For comparison we also present estimates of parameters of the effective action for the blob calculated from the large n t simulations (see Table 1 and 2).
It is quite surprising that the same effective action is still present in the stalk, despite the volume behavior seems, at the first sight, quite different from that in the blob range. It is even more surprising that the parameters of the fit agree quite well with the effective Newton constant Γ measured in the blob range. The parameter µ is slightly bigger but of the same order of magnitude as the potential coefficient from the blob range. Only λ which is related to the size of the dynamically created universe is quite different, but that should be no surprise since λ semiclassically is related to the size of the universe. Finally the value of ρ is difficult to explain from the point of view of higher powers of R(t) which should give ρ = 2k/3 − 1, k = 2, . . ., but also it should be mentioned that it is not very well determined from the fits.
Summing up: S stalk ef f = t 1 Γ (n t − n t+1 ) 2 n t +n t+1 −2n 0 +µ n t +n t+1 2 (50) Having determined by a best fit the parameters of the effective action we can calculate the "theoretical" transfer matrixM (th) using (43) and (50). To appreciate the quality of this approximation we present a plot of six lowest eigenvalues of the measuredM and the "theoretical" matrixM (th) as well as the comparison of their six lowest eigenvectors. In each case the continuous line corresponds to theM (th) (Fig. 18 -19).
7 Discussion and conclusions.
CDT comes with a transfer matrix T |M|T . The way CDT is defined allows us to measure certain distributions, say P ttot (n t ) and P ttot (n t , m t+∆t ), of threevolumes n t . These distributions have an exact definition in terms of the transfer matrix M and the density matricesρ(n) which are projections onto the subspace  Figure 19: The first six eigenvectors of the measuredM (blue dots) and "theoretical" matrixM (th) (red line).M (th) was calculated using the effective action for the stalk range (50).
of three-dimensional combinatorial triangulations T (3) of S 3 spanned by the triangulations with n tetrahedra. (see eqs. (7)- (9)). While the transfer matrix is defined on the large vector space spanned by the elements in T (3), the actual data coming from Monte Carlo simulations seem to allow for a much simpler description in terms of an "effective" transfer matrix M , only labeled by abstract vectors |n referring only to the three-volume. Not only that: basically over the whole range of n t the data are described by a transfer matrix which can be represented as n|M |m = N e −L ef f (n,m) , where L ef f (n, m) is given by and with a corresponding effective action S ef f = t L ef f (n t , n t+1 ).
The last term in (52) unfortunately is not very well determined. For large n t we can not really observe it. The first terms seem to fit that data perfectly with the present statistics. For small n t one can detect a term like δ n+m 2 −ρ , but as mentioned, still ρ is not well determined, and in addition the value we obtain depends on the specific merging of the three different distributions one observes for small n t . Thus we cannot really claim that we have a result which is discretization independent. We are caught in an unfortunate dilemma: we want to go to small n t in order to observe this term, which indicates corrections to the simplest minisuperspace action. However, taking n t small also brings us into the region where discretization effects are likely to be important. It is possible that one can find a window where n t is small enough for the term to be observed via high statistic measurements, but where n t is large enough to avoid discretization effects, but we have not yet pursued this in a systematic way. This discussion highlights an important advantage of the present method: since t tot is small, we are effectively simulating much smaller systems than in the traditionally "full" CDT computer simulations. In this way we can actually obtain measurements of high statistics with relatively moderate computer resources and in a finite amount of time.
The amazing accuracy with which the effective transfer matrix seems to be described by eq. (51) indicates that one obtains a good approximation to the partition function by writing: where S ef f [n t i ] is the effective action (52)-(53). We have strictly speaking only shown that this expression is a good approximation for some special values of the bare coupling constant of the Einstein-Hilbert action S R , given by eq. (24). However, without much doubt any choice of the bare coupling constants which places us well inside the so-called de Sitter phase will allow for a description in terms of an S ef f [n t i ], just with different Γ, µ and λ. Let us assume that this is also true in the two other phases which have been observed in the "full" CDT theory. If this is the case one can actually use S ef f [n t i ] and the expression (54) to study the phase structure of CDT. This is of course much easier than using the full system. This is precisely what has been done in a recent paper [25]. Seemingly one obtains a good qualitative description of the CDT phase diagram (and also new interesting phase structures), corroborating the conjecture that the functional form (52)-(53) might be sufficient to describe CDT for all choices of the bare coupling constants of the Einstein-Hilbert action (24). Checking this is an obvious task for the future. However, an even more interesting application of the multi-canonical Monte Carlo simulation method developed here is that it might allow us to investigate the CDT phase transitions in more detail. A possible UV scaling limit of the CDT theory has to be associated with these phase transitions.