Emergence of the Circle in a Statistical Model of Random Cubic Graphs

We shall consider a discrete statistical model $(\Omega,\mathcal{A})$ defined schematically by the partition function

where $\Omega$ is some configuration space consisting of graphs, $\mathcal{A}:\Omega\rightarrow\mathbb{R}$ is an action functional and $\beta=\beta(\tilde{\beta},N)$ an a priori $N$-dependent function which acts as a scale-dependent parameter for the model. We call $\beta$ the inverse temperature of the system. Typically we will consider $\Omega$ consisting of random regular graphs at fixed $N$ subject to an additional constraint discussed subsequently, and investigate different values of $N$. Note that in principle when we say graph we mean abstract graph, but in practice abstract graphs are somewhat difficult to work with and simulations will use labelled graphs. The number of labellings of an abstract graph $\omega$ with $N$ vertices is $N!/|\text{Aut}(\omega)|$ where $\text{Aut}(\omega)$ is the automorphism group of $\omega$ and so at fixed $N$ we over-count each configuration in the partition function a fixed number of times unless global symmetries are present. Since a typical graph has no nontrivial automorphisms we have neglected this latter consideration. Note that since the Gromov-Hausdorff limit is characterised invariantly, there is in fact no strict need to consider abstract graphs. The action $\mathcal{A}$ is a discrete Einstein-Hilbert action, defined:

where $N_{\omega}(u)$ denotes the set of neighbours of $u$ in $\omega$ and $\kappa_{\omega}(e)$ is the Ollivier curvature of the edge $e\in E(\omega)$. A more complete presentation of the Ollivier curvature is given in appendix A, but for present purposes it is sufficient to recognise that the Ollivier curvature is a coarse version of the Ricci curvature in the following sense: consider two points $x$ and $y$ in a manifold $\mathcal{M}$ that are separated by a sufficiently small distance $\ell$, as well as the (unique) vector field $X$ parallel to the geodesic connecting $x$ and $y$. Then we have:

In this way the Ollivier curvature represents a discretisation of the manifold Ricci curvature. The action 2 thus corresponds to a discretisation of the (Euclidean) Einstein-Hilbert action as long as the edges incident to a vertex span the tangent space.

The Ollivier curvature is discrete in that it takes values in the rational numbers $\mathbb{Q}$ and is also local in the following sense: let $\omega$ be a graph; for each edge $uv\in E(\omega)$, we have

where $C(uv)\subseteq\omega$ is a subgraph of $\omega$ called a core neighbourhood of $uv$. For our purposes it is sufficient to assume that $C(uv)$ is the induced subgraph of $\omega$ with the vertex set

where $\pentagon(uv)$ is the set of non-neighbours of $u$ and $v$ that lie on a pentagon supported by the edge $uv$. Discreteness and locality are of course naively desirable properties for a quantised gravitational coupling.

It turns out that for certain classes of graph the Ollivier curvature may be evaluated exactly. We shall need a little notation.

• •

  $\triangle_{uv}$ denotes the number of triangles supported on the edge $uv$.

• •

  Consider the induced subgraph on the set $V(C(uv))/\set{uv}$. This will have $K$ connected components—each roughly corresponding to a cycle—labelled with the lower case letter $k$. We shall call these connnected components the components of the core neighbourhood.

• •

  $\square^{k}_{w}$ denotes the number of vertices in the $k$th component of the core neighbourhood that neighbour $w\in\set{u,v}$ such that the shortest cycle support by $uv$ that they lie on is a square.

• •

  $\pentagon^{k}_{w}$ denotes the same as $\square^{k}_{w}$ for $w\in\set{u,v}$ except the shortest cycle is a pentagon instead of a square.

Using this notation we have the following expression for the Ollivier curvature in cubic graphs [56]:

for each edge $uv\in E(\omega)$ where the sum over $k$ runs over the components of the core neighbourhood, $a\land b\coloneqq\inf\set{a,b}$ and $[a]_{+}\coloneqq\max(a,0)$ for any $a\in\mathbb{R}$. Note that the main property of $3$-regular graphs that allows this expression to be derived is the severe restriction on possible core neighbourhoods imposed by regularity of low degree.

We will not, in general, study the full configuration space of $3$-regular graphs, instead imposing additional kinematic constraints on $\Omega$. One particularly attractive constraint, derived in [57] is the so-called independent short cycle condition: we say that an edge $uv$ has independent short cycles iff any two short (length less than $5$) cycles supported on the edge share no other edges. Graphs satisfying this condition admit an exact expression, independently of any additional constraints on the core neighbourhoods due to regularity. Furthermore the quantities appearing in the exact expression have an unambiguous interpretation in terms of numbers of short cycles. In particular the constraint ensures that

for all $k$. Thus

where $\square_{uv}$ and $\pentagon_{uv}$ are the number of squares and pentagons supported on the edge $uv$ respectively. For $3$-regular graphs, we may thus express the curvature of an edge as

This expression for the curvature allows us to rewrite the action. Defining

we have:

where $\triangle_{\omega}$, $\square_{\omega}$ and $\pentagon_{\omega}$ denote the total numbers of triangles, squares and pentagons in the graph $\omega$ respectively. $\mathcal{A}_{MF}$ is determined by global quantities and thus represents a kind of mean field contribution to the action.

The precise significance of the independent short cycle condition is not entirely clear. In [57] it was a kind of ‘integrability’ constraint insofar as it allows one to write the action rather explicitly. It also functions similarly to standard hard core constraints in statistical mechanics and prevents short cycles from ‘condensing’ on an edge, though to a certain extent this is also guaranteed by regularity. For $3$-regular graphs, the main utility of the hard core condition appears to be the dynamical suppression of triangles it entails which we will see is an important requirement for the emergence of geometric structure in the model, though it appears that a somewhat weaker constraint is sufficient for this purpose.

Having expressions of the form 6 and 9 for the Ollivier curvature permits the efficient running of simulations studying statistical models $(\Omega(N),\mathcal{A})$, where $\mathcal{A}$ is the Einstein-Hilbert action above and $\Omega(N)$ is the class of cubic graphs on $N$ vertices. (Note that since each $\omega\in\Omega(N)$ is regular with odd degree, $N$ must be even.) We use elementary Monte Carlo techniques [77], evolving the graphs at each step via edge switches: given the random edges $uv$ and $xy$ in a graph $\omega$ such that $u\nsim x$ and $v\nsim y$, we construct a new graph $\tilde{\omega}$ by breaking the edges $uv$ and $xy$ and forming the edges $ux$ and $vy$, subject to any additional constraints on the configuration space that we may wish to impose.

As we have tried to make clear in the introduction, the main purpose of this paper is not to study the problem of (Euclidean) quantum gravity proper but the generation of flat geometries in a model of random graphs. We believe the analogy between the action 2 and the Einstein-Hilbert action are sufficient to call this model combinatorial quantum gravity, as we have done previously. At the same time it would be false to claim that we have no pretensions to addressing Euclidean quantum gravity and for the purpose of clarity it seems desirable to try and explain the present position of our approach in relation to ordinary Euclidean quantum gravity.

Clearly our hope is that there is some Ollivier curvature based action—which we expect to look rather like the action in equation 2—that will act as a discrete regularisation of the Einstein-Hilbert action on graphs $\omega$ which are sufficiently close to a manifold $\mathcal{M}$ in the sense of Gromov-Hausdorff. As mentioned in the introduction, precise convergence results for the Ollivier curvature are not yet known, and until they are we are somewhat restricted in the exact analysis of quantum gravity in general. The exception is the Ricci-flat sector where the convergence problem is absent—since the discrete and continuous Einstein-Hilbert actions vanish trivially. Thus on this sector Gromov-Hausdorff convergence alone guarantees agreement of our model with ordinary Euclidean quantum gravity.

What are the prospects of a precise convergence result? There are essentially two required steps: first we need to show that the Ollivier curvature is respected by Gromov-Hausorff limits, and secondly we need to specify a discrete action in terms of the Ollivier curvature that converges to the Einstein-Hilbert action. Recently there has been major progress with regards to the first step. In his original work, Ollivier [81] demonstrated the stability of the Ollivier curvature under Gromov-Hausdorff limits, in the sense that given a sequence of metric-measure spaces $X_{n}\rightarrow X$ and pairs of points $(x_{n},y_{n})\rightarrow(x,y)$ with $(x_{n},y_{n})\in X_{n}\times X_{n}$ we have $\kappa(x_{n},y_{n})\rightarrow\kappa(x,y)$. The convergence $X_{n}\rightarrow X$ is, however, Gromov-Hausdorff convergence augmented by additional assumptions controlling the Wasserstein distance between the push-forward measures of the points under isometric imbeddings. The difficulty, of course, is that it is these auxiliiary assumptions which represent the foremost challenge in showing the convergence of Ollivier curvature in general. Much of this challenge has been addressed in the recent paper [51] which shows that there is pointwise convergence of the Ollivier curvature—suitably rescaled—in random geometric graphs in arbitrary Riemannian manifolds. This is the first rigorous demonstration of the convergence of a notion of network curvature to its Riemannian counterpart known to the authors. One interesting feature is the need for two distinct length-scales that that are both sent to $0$ in the continuum limit: one is the microscopic ‘edge’ scale defined by the threshold for connecting points obtained by a Poisson point process and the other is an effective curvature scale governing the rescaling of the Ollivier curvature and the radius of the unit balls used for comparison in the random geometric graph. From a gravitational perspective, the kind of point-wise convergence described in [51] must be augmented by showing that the limit is in fact generally covariant; on the other hand it is perhaps a stronger requirement than is physically necessary since one expects only physically meaningful quantities to converge in general. Nevertheless this result considerably expands the potential of discrete-curvature quantum gravity models.

Assuming that there is indeed a precise convergence result, we may turn to some general problems of discrete approaches to quantum gravity. One basic question is whether the partition function 1 gives well-defined dynamics in the $N\rightarrow\infty$ limit. This of course requires energy-entropy balance, i.e. that $\beta(N)\mathcal{A}$ and $S$ have the same $N$-dependence where $S$ is the entropy of the model. In practice we use the energy-entropy balance condition to fix the $N$-dependence of $\beta$; $\mathcal{A}$ grows as $N$ so we need to know the $N$-dependence of the entropy. This depends rather strongly on our choice of $\Omega$; in the present paper we consider $\Omega$ consisting of random regular graphs, and the number of such graphs on $N$-vertices is known. Specifically, we have [110] that the number of $d$-regular graphs on $N$-vertices is

Using the Stirling approximation, taking the logarithm and using its properties we thus get the following naive estimate for the entropy:

We do not, however, simply consider random regular graphs, but instead random regular graphs satisfying an additional hard-core constraint. This constraint will have the effect of reducing the number of configurations and may lead to corrections in the entropy, though these—if they exist—are hard to compute. Indeed, below we find numerically that for $N$ large, $\beta$ is in fact constant with $N$ suggesting that the independent short cycle condition leads to a logarithmic correction to the entropy, at least in the case of $3$-regular graphs. We do not believe that this constant growth of $\beta$ is a generic feature of the model; instead we see it as a consequence of considering cubic graphs which, as we argue below, correspond to one-dimensional geometries and consequently a severely restricted set of action minimising configurations.

Conceptually the issue of the $N$-dependence of $\beta$ is an expression of a well-known issue with local actions in discrete quantum gravity [18] and can perhaps be interpreted as an expression of the nonlocality of the dimensionless action $\beta\mathcal{A}$. From a network theoretic perspective it is an expression of the infinite dimensionality of network phase transitions. To see how these perspectives relate, consider the standard square-lattice Ising model in $D$-dimensions. Such a model is finite-dimensional because the number of possible local interactions experienced by any bulk spin is fixed regardless of the system size. The average behaviour of these spins then only depends on the value of $\beta$, with the same average effect arising from the same value of $\beta$ for a bulk system whatever the system size. Since it is these local interactions which determine the presence or absence of long-range order in the system we may control the phase of the system with a $\beta$-independent of the system size. In a graph, where local interactions are modelled by edges, the number and type of possible interactions depends on the system size and so the parameter $\beta$ can only have the same average effect on a vertex as the system grows if $\beta$ also grows with the system size to compensate for the additional possible interactions.

Of course, like any discrete approach to quantum gravity, we must decide whether discreteness is fundamental as in causal set theory, or simply a regularisation technique that may be removed as a cut-off is removed in line with the asymptotic safety scenario. Both points of view require the recovery of a more or less geometric scaling limit, which as argued in the introduction we interpret in terms of Gromov-Hausdorff convergence. We adhere to the—more conservative—latter attitude which further demands a continuous phase transition; the main result of this paper is that both of these aims can be achieved by our model in the Ricci flat sector, i.e. precisely where our model potentially agrees with quantum gravity.

Finally, let us briefly comment on some defects of the Euclidean approach to quantum gravity. One major well-known problem is that the continuum Euclidean action is not positive definite for $D>2$ since one may choose a metric with conformal mode undergoing arbitrarily fast variations. Such a problem is of course immediately removed upon discretisation but may reappear in the $N\rightarrow\infty$ limit. In the model discussed here—which recall is not quantum gravity—this problem does not arise firstly because the Ollivier curvature is bounded between $-2$ and $1$ for unweighted graphs and secondly because the independent short cycle condition effectively excludes positive curvature (negative action) geometries. More generally, as described in sections 1.8 and 1.9 of [13], it is possible that the partition function is concentrated on configurations with bounded action near a non-Gaussian UV fixed point since the effective Euclidean action contains an entropic term coming from the number of configurations which share the same value of the action. Since our claims in this paper essentially amount to the existence of a UV fixed point, it seems quite possible that our approach may permit a similar escape from the problem of unboundedness.

The biggest problem with our approach from the perspective of quantum gravity proper is the fundamentally Euclidean nature of the approach; in particular this refers to the absence of causal structure and difficulties related to giving sense to some notion of ‘Wick rotation’. A related issue is the unitarity of the resulting quantum theory. We have little concrete to say on these matters at present.

What are the classical configurations of the model? By the specification 2, action minimising configurations maximise the total curvature $\sum_{e\in E(\omega)}\kappa_{\omega}(e)$. In [33], Cushing and collaborators have given a classification of positively curved cubic graphs which goes a long way towards a classification of $\Omega_{0}$. (Note that we say a graph is positively curved iff all its edges have curvature at least zero.) Essentially we find that a positively curved graph is either a discrete cylinder or a discrete Möbius strip. More precisely:

• •

  Consider a chain of $n$ squares, as in figure 1. A prism graph of length $n$, denoted $P_{n}$, is the graph obtained by identifying $u\cong x$ and $v\cong y$.

• •

  The Möbius ladder of length $n$, denoted $M_{n}$, is obtained by gluing $u\cong y$ and $v\cong x$ in figure 1.

With this terminology the classification of Cushing et al. [33] may be summarised as follows:

• •

  If a $3$-regular graph $\omega$ has positive Ricci curvature for all vertices then it is either a prism graph $P_{m}$ for some $m\geq 3$ or a Möbius ladder $M_{n}$ for some $n\geq 2$.

• •

  If $m=3$ the prism graph $P_{m}$ has edges $e\in E(P_{m})$ with $\kappa_{P_{m}}(e)>0$. Otherwise $P_{m}$ is Ollivier-Ricci flat, i.e. $\kappa_{P_{m}}(e)=0$ for all $e\in E(P_{m})$.

• •

  Similarly if $n=2$ the Möbius ladder $M_{n}$ (which is also the complete graph on $4$-vertices $K_{4}$) has strictly positive curvature for each edge $e\in E(M_{2})$. Otherwise $M_{n}$ is Ollivier-Ricci flat.

The key point to note is that for large $N$, a positively curved graph is Ollivier-Ricci flat and hence has total curvature zero. Moreover since the graphs in question have such an obvious geometric interpretation, it is tempting to see this as a model with a geometric classical phase. The issue, of course, is that one can imagine the situation where a graph has both positive and negative curvature edges with positive total curvature. Indeed such configurations may be constructed quite simply, and some examples are given in figure 2. Inspection of the expression 6 immediately shows that triangles necessarily appear in such configurations, and if we are to obtain a model with a classical phase $\Omega_{0}(N)=\Omega_{0}\cap\Omega_{N}=\set{P_{N},M_{N}}$ where $\Omega_{N}$ denotes the class of $3$-regular graphs on $N$ vertices, it is sufficient to ensure that triangles are suppressed. This can of course be achieved by fiat—by restricting to configurations of girth greater than $3$ or bipartite graphs, for instance—but following [57] we know that the suppression of triangles is a dynamical consequence of the model given the independent short cycle condition.

To see this first note that $\mathcal{A}_{MF}$ is an extensive quantity, i.e. $\mathcal{A}_{MF}\sim N$. The scaling of $\mathcal{A}_{P}$ and $\mathcal{A}_{Q}$ depends on the behaviour of $|P|$ and $|Q|$ as $N\rightarrow\infty$. In [57] it was argued that $\mathcal{A}_{P}$ and $\mathcal{A}_{Q}$ do play an important role as $\beta\rightarrow\infty$, but in the random phase $\beta\rightarrow 0$, short cycles are sparse and we may analyse the low $\beta$ dynamics by considering $\mathcal{A}_{MF}$ alone.

Naively we expect that $\mathcal{A}_{MF}$ is minimised by increasing the number of short cycles with the greatest effect coming from an increase in the number of triangles and the least effect coming from an increase in the number of pentagons. Indeed, an edge switch that converts a pentagon to a square leads to a reduction in the action of $1$, a pentagon to a triangle a reduction of $4/3$ and a switch from a square to a triangle a reduction of $1/3$. On these grounds we expect the number of triangles to increase in the random phase. However there is an alternative way of looking at the independent short cycle constraint which leads to different conclusions. In particular the independent short cycle condition can be viewed as an excluded subgraph condition with the abstract graphs in figure 3 excluded. The key point to note is that the subgraphs 3(a) and 3(c) ensure that neither can an edge support two triangles nor can it support a triangle and a square. As such any triangle must share each of its edges with a pentagon and the gain in the action of $1/3$ that arises from the conversion of a triangle to a square is more than made up for by the potential loss of $3$ in the action arising from the conversion of each of the three pentagons to a square. Hence, due to ‘kinematic’ constraints on the configuration space, triangles are strongly unfavoured in the local dynamics of the model in the random phase, and as such we expect them to be absent by the time the geometric phase is reached. On another note we expect pentagons to be similarly suppressed. Figure 4 indicates that these expectations are corroborated.

Figure 5 shows examples of classical configurations observed following simulations run in a configuration space consisting of graphs satisfying the independent short cycle condition. In the simulations we looked at an exponential sweep of different values of $\beta$ from $-10\leq\log(\beta)\leq 10$ and allowed the graph to ‘thermalise’ for $3N/2=|E(\omega)|$ sweeps at each value of $\beta$, where each sweep consists of $3N/2$ attempted Monte Carlo updates, i.e. one attempted update per edge on average. As desired one obtains both prism graphs and Möbius ladders indicating that we have an effective model in which the classical configurations are $\Omega_{0}(N)=\set{P_{N},M_{N}}$ with high probability. In fact, for smaller configurations (each possible value of $N$, from $N=20$ to $N=50$) we have studied the appearance of prism graphs and Möbius ladders systematically and have not found a single exception to one of these configurations appearing in 100 runs for each graph size. Larger graph sizes have not been studied systematically, but again we have not observed a single configuration $\omega$ at $\beta=\exp(-10)$ that is neither a prism graph nor a Möbius ladder, over the course of several runs of each graph size from $N=100$ to $N=500$ at intervals of $50$ and then from $N=500$ to $N=1000$ at intervals of $100$.

In the preceding section we provided strong numerical evidence that given a statistical model with configuration space $\Omega$ consisting of $3$-regular cubic graphs with independent short cycles, the classical phase $\Omega_{0}$ would consist overwhelmingly of prism graphs and Möbius ladders. This was to be expected: it can be proven that the classical phase consists of these graphs if we restrict to graphs without triangles while a strong heuristic argument and numerical evidence both point to this constraint effectively arising due to the low $\beta$ dynamics. We now consider the thermodynamic limit $N\rightarrow\infty$ of the classical configurations, given that $\omega_{N}\in\Omega_{0}(N)=\set{P_{n},M_{m}}$, $N=2n$. We show that in the thermodynamic limit, the classical configurations are effectively one-dimensional and argue that the correct limiting geometry should in fact be $S^{1}$. It turns out that thus intuition can be rigorously formulated in terms of so-called Gromov-Hausdorff limits as we prove in appendix B.

The first point to note is that dimensional evidence in the model already considered of $3$-regular graphs with independent short cycles suggests that the limiting geometry as $N\rightarrow\infty$ is one-dimensional. Dimensional data specifically refers to numerical data on the spectral and Hausdorff dimensions $d_{s}(\omega)$ and $d_{H}(\omega)$ of the graphs in question [28, 25, 37]; heuristically speaking, these quantities respectively define the dimension experienced by a particle diffusing in a space and the dimension governing the scaling of volume with distance in the space. Moreover both spectral and Hausdorff dimensions have become important quantities in quantum gravity: c.f. e.g. [50, 6, 11, 8, 53, 38, 25, 28, 94, 74]. For our purposes, their present significance lies in the fact that we find $d_{s}(\omega)\sim d_{H}(\omega)\sim 1$ for classical configurations $\omega\in\Omega_{0}$ as $N\rightarrow\infty$. Figure 7 shows the relevant plots. We shall briefly describe how the spectral and Hausdorff dimensions are defined here.

The spectral dimension is defined for spaces equipped with a Laplacian $L$ and controls the scaling properties of the eigenvalues of the Laplacian $\lambda\in\sigma(L)$; note that $\sigma(L)$ is the spectrum of the Laplacian which is assumed to be a finite-dimensional operator. Indeed, defining the heat kernel trace as the function

for $t$ in some open subset of $\mathbb{R}$, we may define the spectral dimension function

The spectral dimension of the space on which the Laplacian is defined is then defined as $d_{s}(t)$ in some limit of $t$ or for some range of values of $t$. For Riemannian manifolds $\mathcal{M}$, for instance, we are interested in the small $t$ asymptotics where it can be shown that

where $D=\dim(\mathcal{M})$. Thus if we define

we find $d_{s}(\mathcal{M})=D$.

For our purposes, we shall calculate the heat kernel trace $K(t)$ taking $L$ as the standard graph Laplacian:

where $D$ is an $N\times N$-dimensional diagonal matrix with the diagonal given by the graph degree sequence and $A$ is the adjacency matrix. $L$ is a finite symmetric matrix; its spectrum $\sigma(L)$ thus consists of the $N$-diagonal components of the diagonalisation of $L$. This ensures that $K(t)\rightarrow N$ as $t\rightarrow 0$ and $K(t)\rightarrow 0$ as $t\rightarrow\infty$; defining the spectral dimension as either the small or large $t$ limit of $d_{s}(t)$ thus gives $0$ trivially and for graphs such a prescription fails to specify an informative quantity of the space. Nonetheless it is possible to find a reasonable interpretation of the spectral dimension for discrete spaces by finding regions of the value $t$ for which the spectral dimension function $d_{s}(t)$ plateaus. The case of the torus is an instructive example: figure 6 shows the characteristic behaviour of the spectral dimension function for our purposes viz. vanishing asymptotics, a peak due to discreteness effects and a plateau at the expected integer dimension. Similar behaviour is displayed by high $\beta$ configurations in figure 7(a) indicating that $d_{s}(\omega)\approx 1$ for classical configurations $\omega\in\Omega_{0}$.

Alternatively, the Hausdorff dimension of a metric space equipped with some background measure effectively governs the volume growth of open balls in the space. The assumption is that for large $r$ the volume scales as

The Hausdorff dimension is then typically defined:

In a graph equipped with the shortest path metric, we take the background measure to be the counting measure and the volume of a ball of radius $r$ centred at a vertex $u$ is simply the number of points within a distance $r$ of $u$. For an unweighted graph this is simply all the points that can be reached from $u$ by a path of length at most $\lfloor r\rfloor$, the largest positive integer strictly less than $r$. If the graph is connected the diameter is finite and $B_{r}(x)=B_{R}(x)$ for all $r\geq R$ where $R$ is the diameter of the graph. Thus if we take the above definition $d_{H}(\omega)=0$ trivially for connected graphs and like the asymptotic definitions of the spectral dimension, the Hausdorff dimension so defined fails to be a useful measure of dimension for discrete spaces. As such we take the modified definition:

where $u\in V(\omega)$ and we assume that the diameter is large; the point is that if $|B_{r}(u)|=\alpha r^{d_{H}}$, then

and $\log(\alpha)/\log(r)$ becomes negligible as $r\rightarrow R$ if the diameter $R$ is sufficiently large. In practice, we will estimate $d_{H}(\omega)$ according to equation 29, that is by taking the gradient of a plot of $\log|B_{r}(u)|$ against $\log r$, and then take the mean of this estimate over all vertices. This estimate is valid as long as $|B_{r}(u)|\sim r^{d_{H}}$ for most values of $r$ in the relevant range, i.e. as long as the plot of $\log|B_{r}(u)|$ against $\log r$ is approximately linear; in fact, for such graphs this approach extends the definition 28 since it remains valid even if the diameter is small with respect to $\alpha$. Figure 8 is a typical example of a plot of $\log|B_{r}(u)|$ against $\log r$, indicating that our estimation procedure for $d_{H}$ is indeed valid. Calculated values of this estimate for classical configurations of various sizes are given in 7(b) and show that the Hausdorff dimension of these configurations is approximately $1$ for large $N$.

As mentioned above we have additional—somewhat circumstantial—evidence for the one-dimensionality of the limiting geometry coming from the apparent instability of the second cellular homology group of our classical geometries in the thermodynamic limit. The key point is that the classical configurations are decidedly not choosing their topology at random. The situation is somewhat puzzling: Ollivier curvature is a local quantity and should not be able to distinguish between the Möbius ladder and the prism graph of the same size; as such there seems to be nothing in the action to distinguish the two possible classical configurations and from this perspective we perhaps expect Möbius ladders and prism graphs to appear in equal number. The numerical evidence, however, belies this expectation. We have the following anecdotal evidence:

• •

  Every run of graphs of sizes $N=2n\in\set{100,200,300,400,500,600,700,800,900,1000}$ resulted in a prism graph.

• •

  Every run of graphs of sizes $N=2n\in\set{150,250,350,450}$ resulted in a Möbius ladder.

Precise numbers of runs in the above vary from a few to 15 runs for $N=100$; the data is not systematic but already makes the possibility of random choice between topologies highly implausible. We also see both prism graphs and Möbius ladders for fairly large graphs (about $400$ points) and the possibility of topology choice being a finite size effect seems slight. With these points in mind it is worth considering a more systematic study of the orientability of classical configurations obtained in simulations of small graphs: figure 9 shows the average value of the quantity

where $\omega$ is an observed classical configuration over 100 runs of the code for graphs of size $N=2n\in\set{20,22,...,48,50}$. As can be seen the classical configuration appears to alternate between orientable and nonorientable results as two vertices—or one square—are added to the graph. In particular noting that $n=\square_{\omega}$ we see that the quantity

appears to contain crucial topological information about the tree-level configurations that actually appear in the simulations: $\omega$ is orientable if $\mathfrak{o}(\omega)=0$ and nonorientable otherwise. It would be curious to see if this can be explained. It is tempting to see some relation to the fact that the diagram

commutes, where $\Omega_{0}$ is the set of classical configurations, $\text{Vect}_{1}(S^{1})$ the set of line bundles over $S^{1}$ and $w_{1}$ the first Stieffel-Whitney class; it is well-known that, in particular, $w_{1}$ is a bijection and $\text{Vect}_{1}(S^{1})$ consists of a cylinder and a Möbius band (the possible classical geometries observed in this model) so $\pi$ is to be interpreted as a mapping that sends prism graphs to cylinders and Möbius ladders to Möbius strips.

Given, then, that our classical configurations are converging to a one-dimensional space, what is the space and in what sense are they converging to that space? The limiting geometry is presumably (hopefully) a connected manifold $\mathcal{M}$. Again, as is well-known this means we have (up to homeomorphism) four possibilities [102]:

1. (i)

   if $\mathcal{M}$ is noncompact and without boundary then $\mathcal{M}=\mathbb{R}$.

2. (ii)

   if $\mathcal{M}$ is noncompact and has a nonvoid boundary then $\mathcal{M}=\mathbb{R}_{+}$.

3. (iii)

   if $\mathcal{M}$ is closed, i.e. compact and without boundary then $\mathcal{M}=S^{1}$.

4. (iv)

   if $\mathcal{M}$ is compact with boundary then $\mathcal{M}=[0,1]$.

Note that $1$-manifolds with boundary have a privileged set of points corresponding to the boundaries: $\set{0,1}=\partial[0,1]$ and $\set{0}=\partial\mathbb{R}_{+}$. If the limiting geometry were to be one of these manifolds with boundary, then, the classical configurations would also presumably have to be endowed with a privileged set of vertices that converge to boundary points. This cannot be the case insofar as we regard our classical configurations as abstract graphs and the limiting geometries are without boundary i.e. either $S^{1}$ or $\mathbb{R}$ respectively.

Naively, since the geometry in question is the limit of a discrete space it is rather more natural to assume that the limit is compact, but in fact it seems clear that the precise geometry obtained in the limit depends on the way that the limit is taken. Suppose we have a cylinder $[0,\ell]\times S^{1}_{r}$, where $\ell>0$ and $S^{1}_{r}$ is the circle of radius $r>0$. We may suppose that the cylinder is discretised in terms of a prism graph $P_{n}$ with the radius

Then we have effectively weighted each edge of the prism graph with an edge length $\ell$. A priori, we have two natural choices: keep $\ell$ fixed and let $r$ diverge as $n\rightarrow\infty$, or alternatively keep $r$ fixed and let $\ell\rightarrow 0$ as $n$ increases. The former corresponds to the non-compact limit $[0,1]\times\mathbb{R}$ which is two-dimensional and ruled out by the arguments above. The latter, on the other hand, corresponds to a one-dimensional compact limit $S^{1}_{r}$; it is also heuristically more in line with the idea that $P_{n}$ is a lattice regularisation of the underlying geometry $[0,\ell]\times S^{1}_{r}$. In particular, insofar as $\ell$ is a lattice cutoff we wish to take $\ell\rightarrow 0$ as $n\rightarrow\infty$. Roughly speaking, the limiting geometry is simply defined as the limit of cylinders $[0,\ell]\times S^{1}_{r}$ as the width $\ell\rightarrow 0$, which is of course simply the circle $S^{1}_{r}$. Similar considerations hold for Möbius strips.

These heuristic considerations on the limiting geometry have a rigorous formulation in terms of Gromov-Hausdorff limits. The fundamental idea of ‘metric sociology’ as Gromov termed it, is that there is a notion of distance between compact metric spaces, that gives us a precise notion of what it means for one compact metric space to converge to another. This is a particularly attractive framework for studying the problem of emergent geometry in a Euclidean context because every compact Riemannian manifold can be obtained as the Gromov-Hausdorff limit of some sequence of graphs [24]. We introduce the Gromov-Hausdorff distance in appendix B and show that the classical configurations and geometries discussed here each converge to $S^{1}_{r}$ as $N\rightarrow\infty$.

We shall study the phase structure of the theory via an examination of the quantities

where $C$ is called the specific heat of the system. Note that $\beta\mathcal{F}=-\log\mathcal{Z}$. In a phase transition we expect to find some nonanalyticity in the free energy $\beta\mathcal{F}$ in the thermodynamic limit and in the Ehrenfest classification the order of the derivative in which a discontinuity is observed gives the order of the phase transition. That is to say in a first-order transition we expect to observe a discontinuity in $\beta\mathcal{A}$ and a delta-function peak in $C$, while in a second-order transition we expect $\beta\mathcal{A}$ to be continuous and the specific heat to contain the discontinuity. Since such symptoms of nonanalyticity only arise in the thermodynamic limit [45] the difficulty lies in the assessment of the phase transition in finite systems.

It has already been argued that two distinct phases exist: the random phase at low $\beta$ and the classical phase as $\beta\rightarrow\infty$. Indeed figure 10 shows the behaviour of some observables including the action for an exponential range of values of the parameter $\beta$ and for a large range of graph sizes. The appearance of the two expected regimes as $\beta$ is varied is quite clear. From the same figure it appears that the phase transition occurs roughly in the region $0.05\leq\beta\leq 20$, and in fact we find that if we allow for longer relaxation times the upper value can be somewhat reduced. Figure 11 shows a plot of the quantity $\beta\braket{\mathcal{A}}_{\beta}$ in the critical region for a smaller range of graph sizes. While we should not draw too many conclusions on the basis of such plots, there does not appear to be a discontinuity anywhere in the plot and as such we already have an indication that the phase transition is not first-order.

A key step in the analysis of the phase transition is an estimation of the critical $\beta_{c}$. In general this procedure is a little subtle since the random graph models under consideration are expected to be infinite-dimensional statistical models. Indeed infinite dimensionality is typical in both random graph models [23, 52] and the statistical mechanics of networks [1, 83] and is also a feature of discrete quantum gravity models [57, 44] essentially as an expression of some expected nonlocality in the system. In practice this leads an $N$-dependent function $\beta_{c}$ with the precise $N$-dependence to be estimated via finite size scaling arguments as in [44]. If $\beta_{c}(N)$ appears to converge to some value $\beta_{c}$ as $N$ is increased we see that $N$-dependence of $\beta_{c}(N)$ detected is simply a finite-size effect and not a consequence of the infinite dimensionality of random graph models. Conversely if convergence is only observed after factoring out some power of $N$ then we see that infinite dimensionality of the model is at play in addition to finite size effects.

The most natural way to estimate $\beta_{c}(N)$ is to find the value of $\beta$ which maximises $C(N)$ where $C(N)$ denotes a numerical estimate of the specific heat $C$ for graphs of order $N$. In a first-order transition we expect $C$ to behave as a delta-function singularity at the transition point since, by definition, the action $\beta\braket{\mathcal{A}}_{\beta}$ has a discontinuity at the transition point. On the other hand, in a second-order transition we expect $C$ to diverge at $\beta_{c}$ as a natural consequence of the appearance of critical fluctuations and the divergence of the correlation length [45]. In both cases, in finite systems amenable to numerical analysis, the divergences are smeared out into finite peaks which become sharper as $N$ increases. Figure 12 shows the specific heat at various values graph sizes; we indeed observe sharp peaks forming as $N$ increases.

Figure 13 shows the value of $\beta$ which maximises the specific heat as a function of the graph size. This effectively shows the scaling of the observed critical point values with $N$. Fitting against a plot of $y=mN^{-1}+c$ suggests that the critical value $\beta_{c}(\infty)=1.84\pm 0.03$; more important than the precise value stemming from this quantitative estimate is the qualitative adequacy of the $N^{-1}$ fit, in line with earlier expectations, showing that the system is finite-dimensional and thus has a well-defined critical inverse temperature.

As already discussed above, a phase transition is first-order if there is a discontinuity in the quantity $\beta\braket{\mathcal{A}}_{\beta}$ and second-order if $\beta\braket{\mathcal{A}}_{\beta}$ is continuous while the specific heat $C$ is discontinuous. Such discontinuities are not easy to spot in finite systems while, as already argued, divergences in $C$ which do manifest in easily detectable traits of plots of finite systems cannot be used to distinguish between the order of the transition. Thus we need to consider other diagnostics if we are to identify the nature of the phase transition in question.

One characteristic feature of first-order transitions which is not apparent in second-order transitions is phase coexistence at the critical value $\beta_{c}$: the action at $\beta_{c}$ can take on any of at least two distinct values separated by the discontinuity that characterises the transition. Regarded as a random variable the action will thus be distributed as a superposition of the distributions in each of the distinct pure phases. Assuming two distinct pure phases with distributions $P_{0}$ and $P_{1}$ respectively, the distribution of $\mathcal{A}$ at the critical temperature is then given

for some $\alpha\in(0,1)$, where $E$ is some event. Assuming that we have enough observations we expect both $P_{0}$ and $P_{1}$ to be normally distributed around (or even entirely concentrated at) some point and the action distribution 39 takes the form of two superimposed Gaussians, becoming more pronounced as the size of $N$ increases. Assuming that the discontinuity is larger than the effective support of the two Gaussians combined, which will occur if we consider sufficiently large $n$, we thus expect to observe two distinct peaks in the frequency distribution of the action at the critical value. This can be observed directly by looking at a frequency histogram of the observations of $\mathcal{A}$, checked for in time-series data where we expect to observe sharp transitions between the two Gaussian centres, or observed in the following modified Binder cumulant:

Up to an overall factor, this is one minus the kurtosis of the distribution, and up a constant shift is the coefficient introduced in [30]. It can be shown that in a first-order transition, an observable taking the double Gaussian distribution 39 will have a non-zero minimum at the transition point and will take on the value $0$ elsewhere [30]. In a second-order transition, on the other hand, the action is distributed according to a single Gaussian with corresponding consequences for the observed frequency histogram and time-series data. In particular we expect $B$ to vanish identically everywhere in a second-order transition.

In finite systems, of course, phase coexistence is observed for both first and second-order transitions. However in the former case it is a fundamental feature of the transition while it is simply a finite-size effect in the latter; as such phase coexistence and the associated phenomena become more pronounced as we increase $N$ in a first-order transition and less pronounced in a second-order transition. More concretely, we expect transitions between phases to become less sharp in time-series data, distinct Gaussians to merge and $B$ to approach $0$ as $N$ increases in a second-order transition. Figures 14, 15 and 16 all corroborate these expectations indicating that we indeed have a second-order transition.

Finally let us consider the divergence of the normalised specific heat. This is characterised in terms of a critical exponent $\lambda$ defined by:

This may be estimated by looking at the collapse of the normalised specific heat $C/N$ in the critical regime. This is shown in figure 17 where we have also shown the normalised dimensionless action; we obtain an estimate of $\lambda=0.15$. Note that the maximum value of the specific heat is indeed increasing as may be seen from figure 18.