Topology induced first-order phase transitions in lattice quantum gravity

Causal Dynamical Triangulations (CDT) is a lattice formulation of quantum gravity, suitable for Monte-Carlo simulations which have been used to study the phase diagram of the model. It has four phases characterized by different dominant geometries, denoted phase A, B, C and Cb. In this article we analyse the A − B and the B − C phase transitions in the case where the topology of space is that of the three-torus. This completes the phase diagram of CDT for such a spatial topology. We observe that the order of a phase transition of spacetime geometries can depend on the topology of spacetime.


Causal Dynamical Triangulations
Causal Dynamical Triangulations (CDT) is a lattice formulation of quantum gravity [1,2] based on the path integral approach of Feynman, where certain piecewise linear geometries provide the lattice structure and the UV cut-off.The action used is the Einstein-Hilbert action provided for such geometries by Regge [3].It is the hope that it can be as successful as Lattice QCD has been in regard to testing non-perturbative aspects of the theory via numerical simulations.The difference between CDT and Lattice QCD is that in the latter case the lattice structure is fixed and is just representing the underlying flat spacetime, while in CDT the geometry is encoded in the lattice connectivity.Since we want to integrate over geometries in the path integral, we have to use lattices with different connectivity structure.More precisely, within the framework of CDT, the formal path integral of quantum gravity is defined as a sum over geometries constructed by gluing together simplicial building blocks satisfying certain topological constraints.The most important constraint is that the quantum geometries admit a global time-foliation into spatial slices of fixed topology which is conserved in the time evolution.In all cases studied before, the topology of the spatial slices was chosen to be either spherical (S 3 ) or toroidal (T 3 ).Even though the formulation seems to be discrete in nature, it is important to stress that the lengths of links (i.e., the edges of the simplicial building blocks) play a role of the UV cut-off, which is supposed to be removed in the continuum limit, if it exists.In 3+1 dimensional spacetime the simplicial manifolds are constructed from two types of building blocks: the {4, 1} and the {3, 2} simplices, where the first number denotes the number of vertices at the (discrete) lattice time t and the second one is the number of vertices at time t ± 1. Within each four-simplex, all vertices with the same time coordinate are connected by space-like links and vertices at the neighboring time layers are connected by time-like links.The lengths of all spatial links are fixed to be a s and of all time-like links to be a t = −α a 2 s , where the asymmetry parameter α > 0 if the signature of spacetime is Lorentzian.The simplices are then glued together along their three-dimensional faces (tetrahedra) in such a way that no topological defects are introduced and the assumed spatial topology and time-foliation structure of the simplicial manifolds are preserved.The imposed global time-foliation allows for a well-defined analytic continuation between Lorentzian and Euclidean signatures.The quantum amplitude can be calculated as a weighted sum over all such discretized fourdimensional geometries, later called triangulations or configurations, joining two geometric states.In practice, one usually identifies the initial and final states and thus assumes time-periodic boundary conditions, leading to the global topology of the CDT manifolds being either T 1 × S 3 or T 4 , in what we call the spherical or the toroidal CDT, respectively, depending on the (fixed) spatial topology choice.In the continuum formulation, the path integral is formally defined as: [g]  (1.1) where [g] denotes the equivalence class of metrics with respect to diffeomorphisms and D M [g] is the integration measure over nonequivalent classes of metrics.This expression becomes meaningful only if one introduces a regularization.In CDT we make a particular choice by introducing the piecewise linear regularization described above, and when performing the Wick rotation from Lorentzian to Euclidean signature, the path integral becomes a partition function: with the Regge action [4]: where N 0 is the number of vertices in a triangulation T , N 41 and N 32 are the numbers of {4, 1} and {3, 2} simplices, respectively.The action is parametrized by a set of three dimensionless bare coupling constants, κ 0 , κ 4 and ∆ corresponding respectively to the inverse gravitational bare coupling constant, the dimensionless cosmological constant and a function of the asymmetry parameter α.The model is analytically solvable in 1+1 dimensions [5] and the 2+1 dimensional model can be associated with matrix models [6][7][8].The 3+1 dimensional model discussed here can be analyzed via Monte Carlo simulations.

Phase transition methodology
The phase structure is such that for given κ 0 and ∆ there exists a κ crit 4 such that for κ 4 > κ crit 4 the average value of N 4 = N 41 + N 32 will be finite, while for κ 4 < κ crit 4 the partition function is not well defined.We are interested in the limit where N 4 → ∞.In order to obtain that in a controlled way, in the numerical simulations where N 4 will fluctuate, we add to the action (1.3) a quadratic volume fixing term: (N 41 − N41 ) 2 .Then, to achieve that N 41 actually fluctuates around N41 in the large-volume limit, one has to add the volume fixing term and fine-tune the bare cosmological constant κ 4 to its pseudo-critical value κ crit 4 ( N41 ).Only for N41 → ∞ will κ crit 4 ( N41 ) approach the true critical value κ crit 4 .The nontrivial behavior of κ 4 near κ crit 4 is transferred to the volume dependence of κ 4 on N41 .It means that for each N41 we effectively fix the pseudo-critical value of the coupling constant κ crit 4 ( N41 ) and we study the properties of the model parametrized by the two remaining coupling constants: κ 0 and ∆.Observing the geometric structures of individual triangulations and their extracted averaged properties, one may find regions where these properties differ from one another.In this way the two-dimensional phase diagram of the model can be determined, and it reveals a rich structure with many phase transition lines, see Fig. 1.The phase transitions in CDT have so far been analyzed using standard techniques.It has been shown that the peak in the variance of an order parameter typically designates, with high accuracy, the location of the phase transitions, and also the Binder cumulant can be used to predict the value of the infinite volume limit of a given coupling constant using finite volume scaling measurements [9][10][11][12][13].Typically, one selects a value of a coupling constant (κ 0 or ∆) and creates a set of configurations with fixed values of one coupling (e.g., κ 0 ) and changed values of the other coupling (e.g., ∆) with step-size δκ 0 or δ∆.After thermalizing these configurations one collects statistics of various observables for different geometries.Collecting large enough statistics is necessary to calculate averages of these observables with high accuracy.Each line (or grid) of measurements gives information about the location of the phase transition and the pseudo-critical value of the observed quantities for fixed N41 . 1 Taking the limit N41 → ∞ gives the continuum behavior of the chosen observable, and the critical exponents can be extracted.The finite size dependence of the observables (which are functions of ∆ and κ 0 ) is non-trivial.In this paper we try to use this dependence to determine the properties of the CDT phase transitions in the infinite volume limit, i.e., when N41 → ∞ and κ 4 → κ crit 4 .One should also note that in the phase transition studies described herein, we have tested and used a new method of finding locations of the phase transition points based on machine learning techniques, namely the so-called logistic regression method which turned out to be very efficient.More details can be found in Appendix 1.

The phase diagram of CDT
The geometry of a CDT triangulation can be investigated by analyzing various observables. 2hey can be general quantities characterizing global properties of a triangulation, e.g., the ratios: O 1 ≡ N 0 /N 41 or O 2 ≡ N 32 /N 41 .They can be as well more local quantities, e.g., the distribution of the spatial three-volume as a function of the (lattice) time coordinate, which can be measured both in the spherical and the toroidal version of CDT, as it is presented in Fig. 2. Other examples of local observables are the recently introduced spacetime coordinates [14], which map a triangulation onto four scalar fields, enabling one to measure and visualize four-volume density distribution in all spacetime directions within the toroidal version of CDT, see Fig. 3. Geometric properties of generic triangulations dominating the path integral (1.2), and thus the typical behavior of the above-mentioned observables, will vary depending on the choice of the bare coupling constants: κ 0 and ∆. 3 Even though the bare action of CDT (1.3) is very simple, the 3+1-dimensional model reveals a surprisingly rich phase structure, see Fig. 1, which seems to be universal, independent of the spatial topology choice [15].So far, four phases with distinct geometric properties were discovered -called A, B, C and C b , respectively.
For large values of the bare inverse gravitational constant κ 0 one encounters phase A. This phase is characterized by the spatial volumes at different (lattice) time layers being uncorrelated [16], which is seen clearly in the upper left panel in Fig. 2. Approaching phase A from the phase C, where one has an effective minisuperspace action [17] (see below), the coefficient in front of the kinetic term, which couples spatial volume in the neighboring time layers, vanishes reflecting the lack of correlation between spatial volumes.On the A-side the system is described by an effective action without kinetic term.On the other hand, the generic geometry of the spatial slices themselves seems to be isotropic and homogeneous, see Fig. 3, left upper panel.
For small or negative values of the asymmetry parameter ∆ one observes phase B,  where the typical spatial volume distribution is collapsed into a single time layer with the three-volume distribution pinched to the cut-off size in the remaining part of a generic triangulation, see Fig. 2. Thus, the four-dimensional geometry is effectively reduced to a three-dimensional geometry, in the sense that all four-volume is contained in a thin slab surrounding a single time layer.Moreover, in such a generic configuration, the isotropy and homogeneity of three-dimensional space is itself maximally broken, in the sense that the whole volume is concentrated around a few vertices.This is illustrated in Fig. 3, top right panel.
For not too large κ 0 , increasing ∆ one moves from the B phase to the so-called bifurcation phase, C b .In this phase the number of constant time layers where the spatial volume is significantly different from the minimal value starts to increase and forms a blob, see Fig. 2.However, when one increases the four-volume N 41 sufficiently, the extension of this blob ceases to increase as N 1/4 41 .This implies that in the limit where N 41 → ∞, the stalk, where the spatial volumes are of cut-off scale, will dominate if we allow the time direction to be of order at least N 1/4 41 .In addition, if we look at the spatial geometry inside time layers belonging to the blob, it shows the same characteristics as the spatial geometry in the B phase: isotropy and homogeneity are broken for a generic configuration by the presence of a few high-order vertices.These high-order vertices are correlated from slice to slice within the blob.As we continue to increase ∆ or κ 0 this effect becomes less and less pronounced [18]: the order of the high-order vertices decreases and the blob broadens but it does not disappear.A stalk will remain.The lower left panel in Fig. 3 shows that not only in the time direction, but also in the spacial directions, the geometry is more extended than in the B phase, but still inhomogeneous.Increasing ∆ or κ 0 even further bring us into phase C.
In phase C, which is also denoted the de Sitter phase or the semiclassical phase, the spatial volume profiles seem to depend on the choice of spatial topology, as illustrated in the right panels in Fig. 2. In the case where the spatial topology is S 3 the volume profiles form a blob, like in the bifurcation phase.The crucial difference is that the extension of the blob scales as N 1/4 41 and a typical spatial volume in the blob scales as N 3/4 41 .Thus, one might obtain a semiclassical geometry of the blob in the limit N 41 → ∞.Similarly, the generic spatial slices seem to be isotropic and homogeneous without the high-order vertices of phase B and phase C b .In the case where the spatial topology is T 3 we do not observe a blob, but rather a constant volume profile with superimposed (small) fluctuations.In both cases, the spatial geometries (except for the geometry in the stalk in the case where the spatial topology is S 3 ) of generic configurations seem to be isotropic and homogeneous on large scales while on shorter scales one can observe volume density fluctuations forming voids-and-filaments [14,19], surprisingly similar to structures formed by matter content of the real Universe, as illustrated in the lower right panel in Fig. 3.In accordance with this, the volume profiles are well described by minisuperspace actions similar to the Hartle-Hawking minisuperspace action [20][21][22][23][24], and the kinetic term in the minisuperspace actions has the same coefficient in the two cases.
The A−C transition has been determined to be a first-order transition [9,12], while the the B − C b transition has been determined to be a second-order transition [9,13] in both S 3 and T 3 spatial topology choices.The C −C b transition was found to be a secondorder transition in the case where the spatial topology was S 3 [10,11] and (seemingly) a first-order transition when the spatial topology was T 3 [15].The B−C transition and the A−B transition could be studied by MC simulations in the case where the spatial topology is T 3 .Preliminary results for the B−C transition were reported in [25], but herein we will improve the data for the B −C transition and complete the phase diagram by measuring the yet unexplored A − B transition.As a result, we will see an interesting pattern emerge, to be discussed in Sec. 5.

The A−B phase transition
The phase transition between phases A and B was not analyzed earlier in detail.This transition is important because there is a point in the phase diagram which is the common endpoint of the A−C, A−B and B −C phase transitions (the A−B −C triple point), see Fig. 1, which is a potential candidate of the UV fixed point of CDT.Thus, understanding the nature of all surrounding phase transitions is crucial.We perform a finite volume scaling analysis and determine the order of the A−B phase transition for three different fixed values of κ 0 = 4.8, 4.6 and 4.5 (and the corresponding observables will be denoted with their corresponding κ 0 coupling as a lower right index).At this phase transition, the time-reduced collapsed configurations on the B side of the transition compete with time-uncorrelated configurations on the A side (see Fig. 2).For consistency, the scaling exponent is measured for the three fixed values of κ 0 independently by varying the values of ∆ to find the (pseudo-)critical values ∆ crit ( N41 ) together with the scaling exponents γ.We performed simulations for various volumes (ranging between N41 = 20k and 720k) to find the scaling exponents for each fixed κ 0 measurement set.Our assumption about the critical behavior of ∆ crit ( N41 ) in the limit N41 → ∞ is a (standard) power function, the same as used in the previous phase transition studies cited before: where C is a constant, ∆ ∞ is the value of the asymmetry parameter corresponding to the infinite volume limit at fixed κ 0 and γ is the scaling exponent.The best fits for the scaling exponents γ were measured to be γ 4.8 = 1.088 ± 0.101, γ 4.6 = 1.029 ± 0.178, and γ 4.5 = 1.151 ± 0.379 for the three κ 0 values considered.As the best fitted exponent values are all consistent with the first-order behavior, i.e., γ = 1, we used fixed γ = 1 to fit relation (3.1) to our data, as presented in Fig. 4. By extrapolating the fits to the infinite volume limits ∆ crit were determined to be ∆ ∞ 4.8 = −0.110± 0.001, ∆ ∞ 4.6 = −0.077± 0.001 and ∆ ∞ 4.5 = −0.059± 0.001 for the three cases respectively.
In the following analysis, we will assume that near the phase transition point not only the critical coupling(s) ∆ crit ( N41 ), but also other relevant order parameters scale similarly to eq. (3.1) with the critical exponent fixed at γ = 1.This is indeed supported by our data.In Fig. 5 we show the results for the O 2 ≡ N 32 /N 41 observable measured for each side of the phase transition independently.The fits shown in Fig. 5 have a close to zero value of O ∞ 2 for all cases, but on the A side of the transition the approach to O ∞ 2 ≈ 0 becomes slower with decreasing κ 0 , i.e., slower when we approach the A−B −C triple point.
A final observation corroborating the first-order nature of the A−B transition is the  fact that the values κ crit 4 ( N41 ) observed in the A phase and the B phase data closest to the transition point do not converge to a common value for N41 → ∞.This is shown in Fig. 6 (where we have used data from Fig. 4 to express N41 in terms of ∆ crit ).It is seen that when we increase N41 , which corresponds to higher ∆ crit in the plot, the gap between κ crit 4 (A-side) and κ crit 4 (B-side) increases.This shows that the MC simulations collect very different configurations on the A and the B phase transition side, a situation typical for a first-order transition.
Summarizing, in this section we have analyzed the behavior of the A−B phase transition for three different fixed values of κ 0 via the scaling of the asymmetry parameter ∆ crit ( N41 ) and O 2 on the two sides of the phase transition.For all measured values of κ 0 , the fitted values of the scaling exponent γ are consistent with the first-order behavior (γ = 1).Thus, we conclude that the A−B phase transition is a first-order phase transition.The order of the A−B phase transition does not exclude the possibility that exactly at the endpoint (A−B −C triple point) the phase transition could be of higher-order.

Revisiting the B −C phase transition
Recently, using manifolds with toroidal spatial topology, we were able, for the first time, to measure the properties of the B−C phase transition [25] at fixed κ 0 = 4.0.This result was important since this region of the phase diagram is between the two common endpoints of CDT phase transition lines (the A − B − C and the C b − B − C triple points), see Fig. 1.We determined the critical behavior based on measurements for a sequence of values of ∆ approaching the phase transition from B and C sides of the phase diagram.Using the approach presented in [25] we established a range of coupling constants, within which, for a given finite value of N41 , the system jumps between the two phases.We determined the limits of this hysteresis region and measured the (volume dependent) pseudo-critical values of the coupling constants ∆ and κ 4 as close to the phase transition as we could get with respect to the hysteresis.We observed that, similarly to the A − B transition case, the position of the phase transition line moves towards larger ∆ and κ 0 values when N41 is increased.Again, the true phase transition is observed in the infinite volume limit N41 → ∞ which corresponds to κ 4 → κ crit 4 (κ 0 , ∆).In [25] we assumed the standard finite volume scaling for ∆ crit ( N41 ), κ crit 0 ( N41 ) and κ crit 4 ( N41 ), analogous to eq. (3.1): We showed that in the limit N41 → ∞, the hysteresis shrinks to zero.Furthermore, we observed that on both sides of the phase transition ∆ crit ( N41 ) and κ crit 4 ( N41 ) scale in the same way.The scaling analysis of this behavior could be characterized by a common critical exponent γ ≈ 1.62 ± 0.25, potentially signaling a higher-order phase transition.As we will show below, an equally good fit to the data can be obtained using γ = 1 and a higherorder correction given by a universal small constant shift of N 41 for data measured on both sides of the transition line.This indicates that the transition might be first-order, which is supported by the analysis of the order parameters O 1 ≡ N 0 /N 41 and O 2 ≡ N 32 /N 41 showing a strong discontinuity between the two phases even in the infinite volume limit, as already noted in [25].
Herein we have repeated the analysis for κ 0 = 4.0 with much larger statistics and we have also performed new measurements for κ 0 = 4.2, which is closer to the A − B − C triple point.The lattice volumes used ranged from N41 = 40k and 1600k.We denote the measurements performed at fixed κ 0 = 4.0 and 4.2 as vertical since in the (κ 0 , ∆) coupling-constant plane only ∆ varies.Similarly, we have performed measurements with fixed ∆ = 0 and −0.02 varying κ 0 values.We refer to these as horizontal measurements.In the horizontal measurements, the lattice volumes ranged from N41 = 40k to 800k.n the following, we will present results in pairs related to the two values of the fixed coupling constant in the horizontal and the vertical measurement, respectively.Even though we now have much better statistics than in [25] and have measured the B −C phase transition at various bare coupling constant locations, approaching the line in different ways (vertically and horizontally), a fit to eq. ( 4.1) still does not allow a very good determination of the scaling exponents γ, γ and γ , indicating that for the considered range of lattice volumes the assumed functional form (4.1) might not be optimal.Typically, we could only determine γ to be between 1.6 and 2.6, but on the other hand γ = 1 does not seem to be a good fit.This leads to the conclusion that higher-order corrections should be included as they can be important for smaller volumes.Below we propose a very simple finite size correction of the kind considered in [26], namely we perform a shift N41 → N41 − const.Thus, instead of using eq.( 4.1), we will use a slightly modified scaling relation, which can be viewed as corresponding to the first-order transition (fixed γ = 1), but with a specific form of a finite size correction: Both scaling relations of eq. ( 4.1) and eq.( 4.2) have the same number of free parameters, i.e., they require a three-parameter fit.We can now compare the two different classes of fits.The parameters for the best fits of the critical scaling of ∆ crit (vertical measurements) and κ crit (horizontal measurements) using equations (4.1) and (4.2) were found via the least squares method.In Table 1 we show values of χ 2 corresponding to the fits presented in Fig. 7.
In Fig. 7 the vertical measurements (for fixed κ 0 ) are in the left plot and the horizontal ones (for fixed ∆) are in the right plot.In both cases, the dashed lines show the fits of  (right panel).The dashed and solid curves represent fits to eq. ( 4.1) and eq.(4.2), respectively.In the left panel (vertical measurements), green data-points are for fixed κ 0 = 4.0 and blue are for κ 0 = 4.2.In the right panel (horizontal measurements), green data-points are for fixed ∆ = −0.02and blue are for ∆ = 0. eq.(4.1) and solid lines are the fits of eq.(4.2).Table 1 and Fig. 7 clearly show that the two classes of fits are really close to each other in all cases; thus, using this information alone, it is hard to distinguish the order of the B −C transition.Additional hint comes from analysis of the κ crit data on both sides of the hysteresis region as close to the transition as one can get.This is presented in Fig. 8 where we compare fits of eq.(4.1) and eq.( 4.2) to these data.The vertical measurements are shown in the left panel of Fig. 8 and the corresponding χ 2 values in Table 2.The horizontal ones are in the right panel of Fig. 8 and the χ 2 values can be found in Table 3. Again, the difference between the two classes of fits is not large, however, our data favor slightly the first-order fits.Summarizing, we have shown that the first-order type of critical scaling relation (with the proposed finite size correction) fits our data at least as well or even better than the (standard) higher-order fits.Therefore, taking into account the behavior of order parameters which, as already shown in [25] and discussed above, on both sides of the transition behave discontinuously even when extrapolated to the infinite volume limit, we conclude that the data favor a first-order B−C transition, but with the admittedly strange property that in the limit N41 → ∞ the region of hysteresis shrinks to zero.

Discussion and conjecture
In this paper, we have analyzed the previously not measured A − B phase transition at three different locations of fixed κ 0 = 4.8, 4.6, and 4.5.The finite size scaling analysis of ∆ crit ( N41 ) revealed a scaling exponent γ close to value 1, which points to a first-order phase transition.Next, we analyzed the B −C phase transition, improving the statistics and extending the previous research to another fixed value of the κ 0 = 4.2.We also analyzed the finite size scaling when approaching the critical line along the perpendicular direction in the coupling constant space where ∆ was fixed at 0 and −0.02.We concluded that the B−C phase transition was also a first-order transition, although we had to use a slightly modified fitting function to obtain acceptable agreement with data.Let us stress that the phase transitions considered are non-standard from the point of view of ordinary phase transitions, since they involve changes in spacetime, rather than changes in field configurations on a fixed spacetime.From the above and our earlier measurements, an interesting pattern emerges, relating the order of the CDT phase transition to changes in the topology of spacetime and it leads to the following conjecture: phase transitions which involve a change in topology will be first-order transitions.
First, let us clarify that by topology we mean effective topology.This is exemplified by the configurations we observe in the C phase if the topology of space is S 3 and we have timeperiodic boundary conditions.Thus, from the outset the topology of the triangulations is that of T 1 × S 3 .This is encoded in the initial configuration when we start our MC simulations.By construction, any MC update of the configuration does not change the imposed topology.However, after many MC updates, a typical generic configuration will effectively look like an S 4 -triangulation where the north-pole and the south-pole of this S 4 are connected by a thin stalk of cut-off scale.While technically the triangulation still has the topology T 1 × S 3 , it is clear that if allowed by the updating algorithm, the preferred configurations would have the topology of S 4 .Secondly, recall from the discussion of the nature of generic configurations in Sec. 2, that we observed a second-order transition between phase B and phase Cb, where there is no change in the effective topology.Likewise, in the case when the spatial topology is S 3 there is no real change in topology when we move from phase C b to phase C and again we observe a second-order transition.All other phase transitions observed are first-order transitions and they are related to a change of the effective topology.In the case of toroidal spatial topology, the C b − C transition involves a change from S 4 in phase C b to T 4 in phase C. Similarly, the B − C transition involves a change from S 4 to T 4 .The effective topology of phase A may be best characterized as the topology of a disjoint union of spatial geometries of various extensions (see upper left panel in Fig. 2).Thus, it is different from the topologies encountered in phase C and phase B. Again, the observed phase transitions between phase C and phase A, as well as the transition between phase B and phase A, are first-order transitions.
With hindsight, it is not surprising that a change in the topology of spacetime at the phase transition might result in a first-order transition, since often one would need some major rearrangements of the configurations to implement the change.In this way, a barrier for such a change can be created and lead to pronounced hysteresis.Such phenomena are of course well known in field theories, where the topology of field configurations might change.The new thing here is that we are discussing the topology of spacetime itself.We were alerted to the relation between topology change of spacetime and phase transitions when we studied scalar fields coupled to geometry [19,27], but the phenomena seen there are in fact already present in the pure theory of geometry as discussed here.
In the CDT theory, our main interest is to find a second-order phase transition where we might be able to define a continuum limit and maybe even a UV fixed point.Since only phase C seems to offer acceptable infrared configurations, it is natural to look for phase transition lines in or at the border of phase C. From the results reported here, we cannot use toroidal spatial topology in such studies, since the transitions will be first-order.One should note that there is still a chance that even though the phase transition lines bordering the C phase are first-order, the triple points where the lines end and meet may be of higher-order.This may be indeed the case of the C b − B − C triple point as the higher-order B − C b transition line ends in this point as well.Nevertheless, this scenario is not very natural and it is hard to be proven numerically.However, if the spatial topology is S 3 , the C b − C transition is a second-order transition.In addition it is then possible that the B−C transition will also be second-order and it makes the C b −B−C triple point an interesting candidate for a UV-fixed point.At least it would then be a point where three second-order transition lines meet.In addition, the C b − C transition has in that case the simple interpretation as the phase transition related to the breaking of isotropy and homogeneity.One could then imagine that the shadow of such a breaking could be important for inhomogeneities in our universe.
to learn (fit) the logistic regression model.Then we use the learned model to classify other data, measured closer to the transition region, as belonging to one of the two phases.The phase transition point is signalled by a sudden change of the probability (A-1) of the measured data to belong to a given phase from ∼ 0 to ∼ 1, see Fig. 9.We have decided to look at data measured for each chosen value of the fixed coupling constant and each lattice volume independently, i.e., for each such a choice we create an independent training set, fit our model and use it to find the critical point.Take, e.g., the A−B transition analysis done for some fixed value of κ 0 and fixed N41 volume based on data measured for different values of the ∆ coupling : ∆ min , ∆ 1 , ...., ∆ k , ∆ max , see Section 3. The procedure to find ∆ crit ( N41 ) is the following: 71.take all data measured for the highest value of ∆ = ∆ max (deep in phase A) and label it as P 1 (or "phase A"), 2. take all data measured for the lowest value of ∆ = ∆ min (deep in phase B) and label it as P 0 (or "phase B"), 3. take a subset of data from points 1 and 2 as a (labeled) training set, 4. learn (fit) the logistic regression model using the training set8 ; this step also includes optimizing metaparameters of the model (λ 1 and λ 2 regularization coefficients), 5. check accuracy of the fitted model using all data from points 1 and 2 (in all cases analyzed the accuracy was 100%, i.e., the model was able to distinguish between the P 1 and P 2 data perfectly), 6. use the fitted model to classify other, yet unlabeled, data measured for ∆: ∆ 1 , ...., ∆ k , i.e., closer to a phase transition region than data from points 1 and 2, 7. compute the probability of all classified data points to be in phase P 1 : Pr(p j ∈ P 1 ), see eq. (A-1), 8. compute the mean value Pr(p j ∈ P 1 ) ∆ (and the error: standard deviation) of the probability of all data points measured for a given ∆, 9. find the transition point ∆ crit ( N41 ), where Pr(p j ∈ P 1 ) ∆ jumps from ∼ 0 to ∼ 1, see Fig. 9.
Then we repeat the analysis (model learning and classification) for a different lattice volume N41 , and then also for each different κ 0 independently.
It is important to stress that in the data analysis described above we used only geometric information about triangulations (configurations) generated in MC simulations, i.e., the set of features used in the ML model fitting and classification does not include any information about the values of the CDT coupling constants or other MC simulation parameters.Therefore we can say, that the method learned to distinguish the phases purely on (some) geometric features of the triangulations.We used a set of 30 geometric observables, including both global parameters of the triangulations: N 0 , N 1 , N 2 , N 4 , N 41 , M O, i.e., the total number of vertices, links, triangles, four-simplices, {4, 1}-simplices and maximal coordination number of vertices, respectively, as well as "per-time-slab" data: N 41 (t), N 32 (t), N 23 (t), N 14 (t), N 0 (t), M O(t), i.e., the number of {4, 1}, {3, 2}, {2, 3} and {1, 4}-simplices in a given slab between (lattice) time coordinate t and t + 1, the number and maximal coordination number of all vertices with a time coordinate t, respectively.In all cases t = 1, 2, 3, 4 (with periodic boundary conditions).In order to increase statistics of our data set and also to encode information about a time shift symmetry of the CDT model, we have quadrupled the data by performing a time shift of all "per-time-slab" features by (periodically) changing their time coordinates t = (1, 2, 3, 4) → (4, 1, 2, 3) → (3, 4, 1, 2) → (2, 3, 4, 1).The values of the global features kept unchanged.The above-mentioned observables also enable one to compute traditional order parameters, which were earlier used to distinguish between the CDT phases, e.g., O 1 ≡ N 0 /N 41 or O 2 ≡ N 32 /N 41 and thus check if the automatic ML classification algorithm is consistent with our former results, see Fig. 9.We implemented the logistic regression model in Wolfram Mathematica 12 using the built-in function: Classify[..., Method -> "LogisticRegression"] with standard parameters.In all cases analyzed the model classified both A−B and C−B phase transition data correctly, with a 100% accuracy, see Fig. 9 as an example.

Figure 1 .
Figure1.Phase-diagram of CDT with the four phases.The plot shows the position of phase transitions measured for fixed lattice volume N41 = 160k in the toroidal CDT.In the spherical CDT case the phase structure is similar, although the phase transition between C and B phases (denoted by the green solid line) was inaccessible by numerical simulations.

Figure 2 .
Figure 2. Spatial volume profiles of generic CDT configurations in different phases.Top: Spherical CDT: A, B, C b , C; Bottom: Toroidal CDT: A, B, C b , C, respectively.

Figure 3 .
Figure 3.The structure of the geometry of phase A (top left), B (top right), C b (bottom left), and C (bottom right), respectively, using the scalar field mapping method defined in [14].The plots show four-volume density distribution of generic CDT configurations in each phase projected on two spatial directions.

Figure 6 .
Figure 6.Values of the critical κ crit 4 ( N41 ) coupling measured as a function of ∆ crit ( N41 ) in phase B (red) and phase A (blue) closest to the phase transition points.The darker the color the lower the corresponding κ 0 coupling, i.e., the closer to the A−B −C triple point.Data suggest the lack of a common value of κ ∞ 4 for the A and B phases at the infinite volume limit.

Figure 7 .
Figure 7. Finite volume scaling of the coupling constants ∆ crit (left panel) and κ crit 0

Figure 8 .
Figure 8. Finite volume scaling of the κ crit 4 in vertical measurements with fixed κ 0 = 4.2 and 4.0 (left panel) and horizontal measurements with fixed ∆ = 0 and −0.02 (right panel).The dashed and solid curves represent fits to eq. (4.1) and eq.(4.2), respectively.In both plots blue color corresponds to data measured in phase B and red in phase C -closest to the phase transition (hysteresis) region.Darker colors correspond to lower values of the fixed coupling constant.

Figure 9 .
Figure 9. Machine learning analysis of the A−B transition for fixed κ 0 = 4.8 and N41 = 100k.Red points are the mean probabilities Pr(p j ∈ P 1 ) ∆ that all MC measurements performed for a given value of ∆ belong to phase A. The probabilities were computed using the logistic regression model (A-1) learned (fitted) using only part of data measured for the lowest and the highest ∆ (denoted by empty dots).For comparison we also plot the traditional order parameters O 1 (rescaled by 5×) and O 2 (rescaled by 35×).

Table 2 .
The χ 2 values related to the individual fits of eq.(4.1) (free γ) and eq.(4.2) (γ = 1) to finite volume scaling of the κ crit 4 in vertical measurements (κ 0 fixed).The C or B indices indicate which side of the phase transition the fitted data belong to.