Criticality and conformality in the random dimer model
Sergio Caracciolo, Riccardo Fabbricatore, Marco Gherardi, Raffaele Marino, Giorgio Parisi, Gabriele Sicuro
CCriticality and conformality in the random dimer model
S. Caracciolo, R. Fabbricatore, M. Gherardi, R. Marino, G. Parisi, and G. Sicuro ∗ Dipartimento di Fisica dell’Università di Milano, and INFN, sez. di Milano,Via Celoria 16, 20100 Milano, Italy Laboratoire de Théorie des Communications, EPFL, 1015, Lausanne, Switzerland Dipartimento di Fisica, INFN – Sezione di Roma1, CNR-IPCF UOS Roma Kerberos,Sapienza Università di Roma, P.le A. Moro 2, I-00185, Rome, Italy IdePHICS Laboratory, EPFL, 1015, Lausanne, Switzerland
In critical systems, the eect of a localized perturbation aects points that are arbitrarily far away fromthe perturbation location. In this paper, we study the eect of localized perturbations on the solution of therandom dimer problem in 2 𝐷 . By means of an accurate numerical analysis, we show that a local perturbationof the optimal covering induces an excitation whose size is extensive with nite probability. We compute thefractal dimension of the excitations and scaling exponents. In particular, excitations in random dimer problemson non-bipartite lattices have the same statistical properties of domain walls in the 2 𝐷 spin glass. Excitationsproduced in bipartite lattices, instead, are compatible with a loop-erased self-avoiding random walk process.In both cases, we nd evidences of conformal invariance of the excitations, that are compatible with SLE 𝜅 withparameter 𝜅 depending on the bipartiteness of the underlying lattice only. Let us suppose that we are given a graph (cid:71) . A dimer on (cid:71) is an edge of (cid:71) with its endpoints, that are said to be “cov-ered” by the dimer. A dimer covering of (cid:71) is a subset ofits edge set such that each vertex is covered by one dimeronly. In other words, a dimer covering is a perfect matchingon the graph (cid:71) [1]. Despite this apparently abstract deni-tion, dimer coverings appear in dierent models and theoriesin statistical physics. Deposition of diatomics molecules, de-fects in crystals, or simple magnetic systems, are all examplesof problems that can be studied as dimer models [2].For this reason, the properties of dimer coverings of a givengraph have received considerable attention. In the 1960s,Kasteleyn [2, 3] and, independently, Temperley and Fisher [4]computed the asymptotic expressions of the number of dimercoverings on innite lattices in 1 𝐷 and 2 𝐷 . Their results areat the basis of the Fisher–Kasteleyn–Temperley algorithm forcounting perfect matchings in planar graphs. Remarkably,they also made clear the correspondence between the solu-tion of a dimer covering problem and the solution of the 2 𝐷 Ising model on a lattice [2, 5].In the present paper, we study a disordered version of thedimer covering problem, namely the random dimer model(RDM). In the RDM, covering an edge by a dimer costs a(xed) edge-dependent random price, so that a total cost ofthe covering is the sum of the costs of the covered edges. Byconsequence, there is an optimal way of performing a cover-ing. The non-trivial “optimal conguration” can be thoughtas a “ground state” of a disordered system. In this paper,we will focus in particular on dimer coverings of 2 𝐷 latticesand we will show that there is a correspondence between theRDM and glassy systems in 2 𝐷 . Unlike the 2 𝐷 Ising model,that is critical at nite temperature and has trivial groundstate, the RDM has a nontrivial ground state and a criticalbehaviour exactly at “zero temperature”, in analogy with thephysics of 2 𝐷 spin glasses. Moreover, we will give evidencesof conformal invariance of the optimal solution, a not obvi- ∗ [email protected] ous property that we will relate to observed conformality of2 𝐷 spin glasses at zero temperature.More specically, we assume that a graph (cid:71) ( (cid:86) , (cid:69) ) is given,with vertex set (cid:86) , | (cid:86) | = 𝑁 , and edge set (cid:69) ⊂ (cid:86) × (cid:86) . A weight 𝑤 𝑒 is associated with each edge 𝑒 ∈ (cid:69) of the lattice. We as-sume that the weights 𝑤 𝑒 are independently and identicallydistributed random variables, having an absolutely continu-ous probability distribution density 𝜚 ( 𝑤 ) . We also assumethat the graph admits more than one dimer covering. We canassign a cost 𝐸 [ (cid:68) ] to each covering (cid:68) as 𝐸 [ (cid:68) ] (cid:66) ∑︁ 𝑒 ∈ (cid:68) 𝑤 𝑒 , (1)and a corresponding Gibbs weight e − 𝛽𝐸 [ (cid:68) ] , depending onthe ctitious inverse temperature 𝛽 . The associated partitionfunction is given by 𝑍 ( 𝛽 ) (cid:66) ∑︁ (cid:68) e − 𝛽𝐸 [ (cid:68) ] , (2)In the 𝛽 → 𝑍 ( ) is simply the number of coverings on (cid:71) . The com-putation of 𝑍 ( ) when (cid:71) is a planar graph has been object ofintense study since the seminal results of Kasteleyn, Fisherand Temperley, recently extended to oriented surfaces [6].We say that (cid:68) ∗ is an optimal covering if 𝐸 [ (cid:68) ∗ ] = min (cid:68) 𝐸 [ (cid:68) ] = − lim 𝛽 →+∞ 𝛽 − ln 𝑍 ( 𝛽 ) . (3)The optimal conguration is almost surely unique and cor-responds to the “ground state” of the model. The introducedrandomness can be thought as due to noise, or impurities. Ina dimer deposition picture, for example, randomness in edgeweights might refer to a space-dependent binding energy ofthe diatomic molecules on the substrate.If (cid:71) is a complete graph or a random graph, the RDM re-covers the “random-link models” studied by statistical physi-cists since the 1980s [7–11].As we will show below, the RDM is also related to randomEuclidean matching problems (REMPs) [12, 13]. In a REMP, a r X i v : . [ c ond - m a t . d i s - nn ] D ec a set of 2 𝑁 points { 𝑥 𝑖 } 𝑁𝑖 = is given, uniformly and indepen-dently generated on a given Euclidean domain in 𝐷 dimen-sion (e.g., the unit hypercube). The goal is to nd the optimalpermutation of 2 𝑁 elements 𝜎 that minimizes the cost 𝐸 [ (cid:68) 𝜎 ] (cid:66) 𝑁 ∑︁ 𝑖 = 𝑑 𝑝 ( 𝑥 𝜎 ( 𝑖 − ) , 𝑥 𝜎 ( 𝑖 ) ) , 𝑝 ∈ R + , (4)where 𝑑 ( 𝑥, 𝑦 ) is the Euclidean distance between 𝑥 and 𝑦 . Eachpermutation denes a pairing, i.e., a set of matched points (cid:68) 𝜎 (cid:66) {( 𝑥 𝜎 ( 𝑖 ) , 𝑥 𝜎 ( 𝑖 − ) )} 𝑁𝑖 = . We denote (cid:68) ∗ ≡ (cid:68) 𝜎 ∗ with 𝜎 ∗ = arg min 𝜎 𝐸 [ (cid:68) 𝜎 ] . In a bipartite variation of the problem, therandom Euclidean assignment problem (REAP), two sets of 𝑁 random points of distinct “colors”, { 𝑥 𝑖 } 𝑁𝑖 = and { 𝑦 𝑖 } 𝑁𝑖 = , aregiven. The pairing has to be such that only points of dierentcolor are matched. In other words, we search for an optimalpermutation of 𝑁 elements 𝜎 that minimizes the cost 𝐸 [ (cid:68) 𝜎 ] (cid:66) 𝑁 ∑︁ 𝑖 = 𝑑 𝑝 ( 𝑥 𝑖 , 𝑦 𝜎 ( 𝑖 ) ) , 𝑝 ∈ R + . (5)Here each permutation denes a set of 𝑁 pairs (cid:68) 𝜎 (cid:66) {( 𝑥 𝑖 , 𝑦 𝜎 ( 𝑖 ) )} 𝑁𝑖 = of points from dierent sets that are matched.Computing the average optimal cost in the REMP or in theREAP is a challenging task, due to the presence of Euclideancorrelations between the pair costs [13–17]. The aforemen-tioned random-link models provide the innite-dimensionallimit of REMPs and REAPs [13].The typical properties of the solutions of these problemsare not trivial. It has been shown, for example, that therandom-link model on the complete graph is “critical”, i.e. itsHessian spectrum is gapless [16]. In Ref. [18], evidences oflong range correlations in the solution of the 1 𝐷 REAP weregiven. On the other hand, the RDM in the 𝛽 → uniformly sampled dimer coverings, theirunion generates a set of curves and paths. Kenyon predictedthe convergence of such curves to a Schramm-Löwner evo-lution SLE [19], and proved the conformality of the loops[20] on bipartite lattices. He also showed that the limit mea-sure of the possible dimer coverings has conformal invari-ance properties on the square lattice [21]. We recall here thatan SLE 𝜅 curve 𝛾 ( 𝑡 ) in the upper complex plane H is given by 𝛾 ( 𝑡 ) = 𝑔 − 𝑡 ( 𝜉 𝑡 ) , where 𝑔 𝑡 ( 𝑧 ) satises the Löwner equation,d 𝑔 𝑡 ( 𝑧 ) d 𝑡 = 𝑔 𝑡 ( 𝑧 ) − 𝜉 𝑡 , 𝑔 ( 𝑧 ) 𝑧 = lim 𝑧 →+∞ 𝑔 𝑡 ( 𝑧 ) 𝑧 = . (6)Here 𝜉 𝑡 is the driving function of the process, and in the caseof a SLE 𝜅 it is given by a Brownian process with (cid:104) 𝜉 𝑡 (cid:105) = (cid:104) 𝜉 𝑡 (cid:105) = 𝜅𝑡 [22, 23].Finally, it is known that there is a special correspon-dence between matchings and spin glasses in two dimen-sions [5, 24]. The problem of nding the ground state of thetwo-dimensional Edwards–Anderson (EA) model [25] can bemapped into a planar matching problem [26, 27]. Notably,the 2 𝐷 EA has a glass transition at zero temperature exhibit-ing scale invariance [28–30], and it has been suggested thatconformal invariance might hold as well. In particular, EA domain walls have been found to be consistently describedby a SLE 𝜅 with 𝜅 ≈ . 𝐷 lattices in the 𝛽 → +∞ limit, i.e. on the ground state, and search for cor-respondences with the REMP, the REAP and 2 𝐷 spin glasses.We will show that similar critical properties appear in someRDMs, in the REMP and in the EA model, suggesting the ex-istence of a unique universality class for these models. More-over, such properties depend on the nature of the underlyinggraph only. a. Models and methods. We analyze three types of lat-tices: the honeycomb (H) lattice, the triangular (T) lattice, andthe square (Q) lattice. Each lattice is obtained considering 𝐿 rows of 𝐿 sites, displaced in such a way that the lattice edgelength is xed to 1. The total number of sites is therefore2 𝑁 = 𝐿 . We impose periodic boundary conditions in bothdirections. For 𝐿 =
4, for example, the three lattices areWe associate to each edge 𝑒 a random weight 𝑤 𝑒 , extractedfrom the exponential distribution 𝜚 ( 𝑤 ) = e − 𝑤 . Once allweights have been assigned, Edmond’s blossom algorithm[32, 33] provides us the optimal weighted dimer covering (cid:68) ∗ .Subsequently, we select a random edge ˆ 𝑒 ∈ (cid:68) ∗ and we cut it.Forbidding the edge ˆ 𝑒 plays the role of a local perturbationthat induces an excitation. Re-running the algorithm on thegraph in which ˆ 𝑒 is not present, we obtain a new optimal so-lution (cid:68) ∗ ˆ 𝑒 of higher cost with respect to (cid:68) ∗ . The dierence Δ 𝐸 ˆ 𝑒 (cid:66) 𝐸 [ (cid:68) ∗ ˆ 𝑒 ] − 𝐸 [ (cid:68) ∗ ] ≥ Δ 𝐸 ˆ 𝑒 = 𝑂 ( ) for 𝑁 (cid:29) (cid:68) ∗ ˆ 𝑒 and (cid:68) ∗ , i.e., (cid:83) ˆ 𝑒 = { 𝑒 ∈ (cid:69) : 𝑒 ∈ (cid:68) ∗ (cid:52) (cid:68) ∗ ˆ 𝑒 } . (7)This selects a set of edges in the original graph. Pictorially,e.g., on the H modelObviously, ˆ 𝑒 ∈ (cid:83) ˆ 𝑒 and it is easy to see that (cid:83) ˆ 𝑒 is a self-avoiding single cycle. If the system is critical, we expect thatthe localized perturbation induces, with nite probability, anextensive rearrangement that aects arbitrarily far points.The cycle (cid:83) ˆ 𝑒 is a fractal object. We will denote by 𝑆 ˆ 𝑒 (cid:66) | (cid:83) ˆ 𝑒 | the number of edges of the cycle, thus its size.Alongside with the RDM, we consider the REMP and theREAP with a total of 2 𝑁 points on a square of size 𝐿 = √ 𝑁 with periodic boundary condition. We adopt the cost inEq. (5) with 𝑝 = V a r [ 𝜗 ] HQREAP V a r [ 𝜗 ] TREMP ln 𝐿 Figure 1.
Left.
Cumulative distribution of the size 𝑠 of the excitationfor the dierent models considered in the papers. The representedsizes are 𝐿 = , ,
500 (from left to right in each case) for eachdimer model. For the Euclidean cases, instead, the represented sizesare 2 𝑁 = , , Right.
Variance of thewinding angle as a function of 𝐿 for all models. Lines are t obtainedusing the function 𝑓 ( 𝐿 ) = 𝑎 + 𝜅 ln 𝐿 . problems exactly as in the RDM, exciting the optimal solu-tion (cid:68) ∗ into a new matching (cid:68) ∗ ˆ 𝑒 by forbidding a pair ˆ 𝑒 ∈ (cid:68) ∗ .In the REMP, the dierence between 𝜎 ∗ and 𝜎 ∗ ˆ 𝑒 can be quan-tied as before considering the cycle (cid:83) ˆ 𝑒 as in Eq. (7). In theREAP the cycle (cid:83) ˆ 𝑒 is obtained starting from the edge set inEq. (7) and then joining consecutive points of the same color:this is because in the REAP the typical distance of a pair inthe optimal solution scales as √ ln 𝑁 [35], whereas the dis-tance between points of the same color is (cid:79) ( ) as in all otherconsidered cases. b. Criticality and fractal dimension. We numericallyevaluated the probability distribution function P [ 𝑆 ˆ 𝑒 > 𝑠 ] ofhaving a cycle of length greater than 𝑠 for all models con-sidered above. Scaling theory for critical systems states thatsuch a probability distribution can be written as P [ 𝑆 ˆ 𝑒 > 𝑠 ] = 𝑠 − 𝜁 𝜌 (cid:16) 𝑠 𝜆 𝐿 − (cid:17) , (8)for some scaling function 𝜌 , such that 0 < lim 𝑧 → 𝜌 ( 𝑧 ) < +∞ .The scaling exponents 𝜁 > 𝜆 > 𝐿 increases, a power-law tail develops in all considered mod-els. This implies that local perturbations induce extensive re-arrangements with nite probability in the thermodynamiclimit and the models are indeed on a “critical point”. A nu-merical estimation of 𝜁 by tting the tail in Fig. 1 also showsthat the power-law exponent is the same for the H and the Qmodel and the REAP, whereas the T model has an exponentvery close to the REMP one, see Table I.Let us now evaluate the fractal dimension 𝐷 f of the cycle.Assuming now that (cid:104) 𝑆 ˆ 𝑒 (cid:105) ∼ 𝐿 𝛼 , where (cid:104)•(cid:105) is the average over H Q REAP T REMPlim 𝑁 (cid:10) 𝑁 − 𝐸 [ (cid:68) ∗ ] (cid:11) . ( ) . ( ) ∞ . ( ) . ( )(cid:104) Δ 𝐸 ˆ 𝑒 (cid:105) . ( ) . ( ) . ( ) . ( ) . ( ) 𝛼 . ( ) . ( ) . ( ) . ( ) . ( ) 𝛾 . ( ) . ( ) . ( ) . ( ) . ( ) 𝜁 . ( ) . ( ) . ( ) . ( ) . ( ) 𝜁 (from t) 0 . ( ) . ( ) . ( ) . ( ) . ( ) 𝐷 f . ( ) . ( ) . ( ) . ( ) . ( ) 𝜅 . ( ) . ( ) . ( ) . ( ) . ( ) + 𝜅 . ( ) . ( ) . ( ) . ( ) . ( ) Table I. Asymptotic average optimal cost and of the scaling expo-nents for random dimer covering problems on the torus. Note thatthe REAP has an anomalous scaling in the cost, namely (cid:104) 𝐸 [ (cid:68) ∗ ](cid:105) ∼ 𝑁𝜋 ln 𝑁 for 𝑁 (cid:29) all instances, from Eq. (8) we have 𝛼 = − 𝜁𝜆 . (9)The gyration radius of the cycle is dened as 𝑅 𝑒 (cid:66) 𝑆 𝑒 ∑︁ 𝑖,𝑗 𝑑 ( 𝑟 𝑖 , 𝑟 𝑗 ) (10)where 𝑟 𝑖 is the position of the 𝑖 th node in the cycle and thesum runs over all pairs of vertices of the cycle. Conditioningon cycles of lenght 𝑠 , it satises a scaling law of the form [36] (cid:10) 𝑅 𝑒 (cid:11) 𝑆 ˆ 𝑒 = 𝑠 = 𝑠 𝐷 − 𝑔 (cid:16) 𝑠 𝜆 𝐿 − (cid:17) . (11)where 𝑔 is a scaling function such that 0 < lim 𝑧 → 𝑔 ( 𝑧 ) < +∞ .Assuming that (cid:10) 𝑅 𝑒 (cid:11) ∼ 𝐿 𝛾 , we obtain 𝛾 = 𝐷 − − 𝜁𝜆 . (12)With reference to Eq. (9) a relation between the three expo-nents 𝛼 , 𝛾 and 𝜁 and the fractal dimension 𝐷 f is easily found, 𝐷 f = − 𝛾 + 𝛼, (13a) 𝜁 = − 𝛾 − 𝛾 + 𝛼 . (13b)The fractal dimension and the power-law exponent 𝜁 can beextracted by a careful measurement of 𝛼 and 𝛾 . The values of 𝛼 and 𝛾 are estimated using the method of ratios [37] [38].All our results are collected in Table I. We observe thatthe values of the exponents naturally splits in two groups:the rst one including the Q model, the H model and theREAP, with cycles having 𝐷 f = . ( ) , and the secondone including the T model and the REMP, with cycles hav-ing 𝐷 f = . ( ) . The fact that the REAP and the H and Qmodels share the same exponents, and similarly the T modelshares its exponents with the REMP, suggests that the onlyrelevant feature that determines the scaling is the nature ofthe underlying graph, i.e., the fact of being bipartite or not.In all cases 𝜁 >
0, i.e., the local perturbation induces a long-range rearrangement with nite probability in the large 𝑁 limit. Moreover, the estimated value for 𝜁 is in agreementwith the result of the direct t, conrming the consistency ofour scaling ansatz. c. Conformality. The presence of criticality suggests theinspection of a stronger invariance, namely conformal invari-ance. Given a cycle of length 𝑠 that can be contracted to apoint (i.e., that does not wind around the torus), we can re-move from it 𝑘 < 𝑠 adjacent steps, by mapping their supportto the boundary of a canonical domain (e.g., the unit disc)[39]. The remaining portion of the cycle can be conformallymapped to a curve in the upper half-plane H with one end-point in the origin and the other at innity. The resultingcurve can be analized then with the so-called zipper algo-rithm [40, 41] to extract the driving function 𝜉 𝑡 of the curveappearing in Eq. (6) [42]. In the hypothesis that the curve re-sults from a SLE 𝜅 , the obtained driving function has to be aBrownian process with (cid:104) 𝜉 𝑡 (cid:105) = 𝜅𝑡 . We perform the analysison the Q and on the T model. The results are given in Fig. 2.The driving function is found indeed to be a Gaussian pro-cess with (cid:104) 𝜉 𝑡 (cid:105) = (cid:104) 𝜉 𝑡 (cid:105) ∼ 𝜅𝑡 with 𝜅 ≈ . ( ) . Althoughwe have not been able to resolve the dierence in 𝜅 betweenthe dierent lattices, the results strongly support that the ob-tained curves are SLE 𝜅 with 𝜅 (cid:39) . (cid:83) ˆ 𝑒 . Then, we choose anorientation at random for the cycle (clockwise or anticlock-wise), and we compute a function 𝜗 ( 𝑒 ) of the cycle edges, insuch a way that, if 𝑒 + 𝜗 ( 𝑒 + ) = 𝜗 ( 𝑒 ) + angle ( 𝑒, 𝑒 + ) . Here angle ( 𝑒, 𝑒 + ) is theturning angle from 𝑒 to 𝑒 + (cid:104) 𝜗 ( 𝑒 )(cid:105) =
0. The variance of the angle is found to growwith 𝐿 , according to the lawVar [ 𝜗 ] = 𝑎 + 𝜅 𝐿, (14)see Fig. 1. The values of 𝜅 are given in Table I. Once again, theT model and the REMP are found to have the same value of 𝜅 .This value is dierent from the value of 𝜅 obtained for the Qand the H model and for the REAP. In the hypothesis that theobtained curves are SLE 𝜅 , Var [ 𝜗 ] has to behave as in Eq. (14),and the quantity 𝜅 appearing in Eq. (14) is exactly the param-eter of the Schramm–Löwner evolution [43–45]. The quan-tity 𝑎 is instead a nonuniversal constant. Moreover, given anSLE 𝜅 , then 𝐷 f = min (cid:0) + 𝜅 , (cid:1) [22]. This relation approxi-mately holds for all the considered models, see Table I. d. Conclusions. The obtained values 𝜁 > − − − − − − − 𝜉 𝑡 /√ 𝑡 P d f QT N ( ,𝜅 ) Figure 2.
Left.
A contractible loop with 𝑠 =
522 in the T modelwith 𝐿 =
400 (top) is split into two parts (red, of 𝑘 =
90 steps,and black, of 𝑠 − 𝑘 steps). By a conformal transformation, the blackpart is mapped into a curve in the upper half plane stemming fromthe origin (bottom). Right.
Probability distribution function of thedriving function for 𝐿 =
400 on the Q lattice and on the T lattice. Thenumerical results have been obtained averaging over a wide range oftime 𝑡 and are compared with a normal distribution (cid:78) ( , 𝜅 ) havingzero mean and variance 𝜅 = . matching problem on a bipartite complete graph). Moreover,the fractal dimension of the cycles obtained in the REMP andin the T model is compatible with the one of 2 𝐷 spin-glass do-main walls, which is estimated to be 𝐷 f = . ( ) [46–48].These results show that non-bipartite matching problems andspin glasses in 2 𝐷 belong to the same universality class. Onthe other hand, the H model, the Q model and the REAP havecycles with 𝐷 f = . ( ) , compatible with the fractal dimen-sion of the loop-erased self-avoiding walk in 2 𝐷 , which has 𝐷 f = and is a SLE [49, 50]. Finally, as test of possible con-formal invariance, we mapped the contractible excitations inthe half upper plane and we computed the driving function 𝜉 𝑡 generating the process. We found, for all analyzed cases, that 𝜉 𝑡 is a Brownian process with variance (cid:104) 𝜉 𝑡 (cid:105) ∼ 𝜅𝑡 , as expectedfor an SLE 𝜅 . As additional test, the variance of the windingangle is found to scale asymptotically as Var [ 𝜗 ] ∼ 𝜅 ln 𝐿 with 𝜅 = ( 𝐷 f − ) , as predicted for SLE 𝜅 .Our results undercover a non trivial connection betweenthe RDM on monopartite lattices, the REMP and the EA spin-glass model in 2 𝐷 . These models are all critical at zero tem-perature and, more importantly, share the same universal-ity class, exhibiting signs of conformal invariance. RDMs onbipartite lattices and the REAP, on the other hand, show asimilar critical zero-temperature behavior, but with dierentscaling exponents. The underlying graph topology is there-fore an essential feature for the determination of the univer-sality class of the critical point. e. Acknowledgments. The authors are grateful to CarloLucibello for providing us his preliminary data about theREMP. The authors would like to thank Bertrand Duplantier,Scott Kirkpatrick, Nicolas Macris and Andrea Sportiello foruseful discussions. The research of G.P. and G.S. has beensupported by the Simons Foundation (grant No. 454949, G.Parisi). R.M. has been supported by Swiss National Founda-tion grant No. 200021E 17554. [1] L. Lovász and D. Plummer,
Matching Theory , AMS Chelsea Pub-lishing Series (2009).[2] P. W. Kasteleyn, in
Graph Theory and Theoretical Physics (Aca-demic Press, 1967) pp. 43–110.[3] P. Kasteleyn, Physica , 1209 (1961).[4] H. Temperley and M. E. Fisher, Philos. Mag. , 1061 (1961).[5] P. W. Kasteleyn, J. Math. Phys. , 287 (1963).[6] D. Cimasoni and N. Reshetikhin, Commun. Math. Phys. ,187 (2007); Commun. Math. Phys. , 445 (2008).[7] H. Orland, J. Phys. Lett. (Paris) , 763 (1985).[8] M. Mézard and G. Parisi, J. Phys. Lett. (Paris) , 771 (1985).[9] G. Parisi and M. Ratiéville, Eur. Phys. J. B , 457 (2002).[10] S. Caracciolo, M. P. D’Achille, E. M. Malatesta, and G. Sicuro,Phys. Rev. E , 052129 (2017).[11] G. Parisi, G. Perrupato, and G. Sicuro, J. Stat. Mech.: TheoryExp. , 033301 (2020).[12] M. Mézard and G. Parisi, EPL , 913 (1986).[13] M. Mézard and G. Parisi, J. Phys. (Paris) , 2019 (1988).[14] S. Caracciolo, C. Lucibello, G. Parisi, and G. Sicuro, Phys. Rev.E , 012118 (2014).[15] S. Caracciolo and G. Sicuro, Phys. Rev. Lett. , 230601 (2015).[16] C. Lucibello, G. Parisi, and G. Sicuro, Phys. Rev. E , 012302(2017).[17] S. Caracciolo, M. D’Achille, and G. Sicuro, Phys. Rev. E ,042102 (2017).[18] E. Boniolo, S. Caracciolo, and A. Sportiello, J. Stat. Mech. The-ory Exp. , P11023 (2014).[19] R. W. Kenyon and D. B. Wilson, Trans. Am. Math. Soc. ,1325 (2011).[20] R. Kenyon, Commun. Math. Phys. , 477 (2014).[21] R. Kenyon, Ann. Probab. , 759 (2000).[22] S. Rohde and O. Schramm, Ann. Math. , 883 (2005).[23] J. Cardy, Ann. Physics , 81 (2005).[24] M. E. Fisher, J. Math. Phys. , 1776 (1966).[25] S. F. Edwards and P. W. Anderson, J. Phys. F: Metal Phys. , 965(1975).[26] L. Bieche, J. P. Uhry, R. Maynard, and R. Rammal, J. Phys. A , 2553 (1980). [27] F. Barahona, J. Phys. A , 3241 (1982).[28] A. J. Bray and M. A. Moore, J. Phys. C: Solid State Phys. ,L463 (1984).[29] T. Shirakura and F. Matsubara, Phys. Rev. Lett. , 2887 (1997).[30] A. K. Hartmann and A. P. Young, Phys. Rev. B , 180404(R)(2001).[31] C. Amoruso, A. K. Hartmann, M. B. Hastings, and M. A. Moore,Phys. Rev. Lett. , 267202 (2006).[32] J. Edmonds, Can. J. Math. , 449 (1965).[33] B. Dezsö, A. Jüttner, and P. Kovács, Electron. Notes Theor.Comput. Sci , 23 (2011).[34] In Ref. [8] it is shown that the 𝑝 = 𝐷 case has an exactly solv-able 𝐷 → +∞ limit. We focused therefore on 𝑝 = 𝐷 = , 259(1984).[36] Here a hyperscaling assumption has been made: we assumethat, to get a proper 𝐿 → +∞ limit, the argument of 𝑔 in Eq. (11)has the same scaling properties of the one of 𝜌 in Eq. (8).[37] S. Caracciolo, R. G. Edwards, S. J. Ferreira, A. Pelissetto, andA. D. Sokal, Phys. Rev. Lett. , 2969 (1995).[38] See the Supplementary information.[39] M. Gherardi, Phys. Rev. E , 032128 (2013).[40] D. E. Marshall and S. Rohde, SIAM Journal on Numerical Anal-ysis , 2577 (2007).[41] T. Kennedy, J. Stat. Phys. , 803 (2008).[42] See the Supplementary information.[43] B. Wieland and D. B. Wilson, Phys. Rev. E , 056101 (2003).[44] B. Duplantier and H. Saleur, Phys. Rev. Lett. , 2343 (1988).[45] B. Duplantier and I. A. Binder, Phys. Rev. Lett. , 264101(2002).[46] O. Melchert and A. K. Hartmann, Phys. Rev. B , 174411(2007).[47] W. Wang, M. A. Moore, and H. G. Katzgraber, Phys. Rev. Lett. , 100602 (2017); Phys. Rev. E , 032104 (2018).[48] F. Corberi, L. F. Cugliandolo, F. Insalata, and M. Picco, J. Stat.Mech.: Theory Exp. , 043203 (2019).[49] S. N. Majumdar, Phys. Rev. Lett. , 2329 (1992).[50] G. F. Lawler, O. Schramm, and W. Werner, Ann. Probab. ,939 (2004). Appendix A: Extrapolation of the exponents.
We extracted the values of 𝛼 and 𝛾 using the method of ratios, i.e., considering, for each 𝐿 , a system of size 𝐿 and a systemof size 2 𝐿 . Assuming now that (cid:104) 𝑆 ˆ 𝑒 (cid:105) ∼ 𝐿 𝛼 , the value of the exponent 𝛼 has been estimated aslog (cid:104) 𝑆 ˆ 𝑒 (cid:105) 𝐿 (cid:104) 𝑆 ˆ 𝑒 (cid:105) 𝐿 = 𝛼 + 𝛼 ( ) 𝐿 𝜔 + 𝑜 (cid:18) 𝐿 𝜔 (cid:19) , (A1)where (cid:104)•(cid:105) 𝐿 denote an average at size 𝐿 . The t of our data, given in Fig. 3, has been performed using a tting function 𝑓 ( 𝐿 ) = 𝛼 + 𝛼 ( ) 𝐿 − 𝜔 , with 𝛼 , 𝛼 ( ) and 𝜔 free parameters. Similarly for the gyration radius, assuming that (cid:104) 𝑅 𝑒 (cid:105) ∼ 𝐿 𝛾 , we havelog (cid:10) 𝑅 𝑒 (cid:11) 𝐿 (cid:10) 𝑅 𝑒 (cid:11) 𝐿 = 𝛾 + 𝛾 ( ) 𝐿 𝜔 + 𝛾 ( ) 𝐿 𝜔 + 𝑜 (cid:18) 𝐿 𝜔 (cid:19) . (A2)In this case, the t has been obtained using a tting function 𝑓 ( 𝐿 ) = 𝛾 + 𝛾 ( ) 𝐿 − 𝜔 + 𝛾 ( ) 𝐿 − 𝜔 , with 𝛾 , 𝛾 ( ) and 𝛾 ( ) free parameters,whilst 𝜔 was xed to the same value estimated in the analysis for 𝛼 . The data points have been obtained averaging over 10 –10 dierent instances for each value of 𝐿 . In the case of lattice models, we used 8 ≤ 𝐿 ≤ . . . 𝛼 REMPT . . . 𝛾 REMPT .
05 0 . . . . 𝛼 HQREAP .
05 0 . . . . . 𝛾 HQREAP 𝐿 − 𝜔 𝜔 H 0 . ( ) Q 0 . ( ) REAP 1 . ( ) T 1 . ( ) REMP 0 . ( ) Figure 3. Extrapolation of 𝛼 (left) and 𝛾 (right) for the dierent models. The value of both exponents in the bipartite models is clearly dierentfrom the one obtained for the monopartite ones. REAP we considered 64 ≤ 𝐿 ≤ 𝜔 are given in the table in Fig. 3. Appendix B: Analysis of the driving function a. Mapping to the standard chordal geometry.
To study the driving function of the excitations of the system, we restrictourselves to the lattice models and we consider cycles Γ = ( 𝛾 , 𝛾 , . . . , 𝛾 𝑠 ≡ 𝛾 ) of length 𝑠 on the torus. By means of the notationabove, we denote the piecewise linear curve Γ passing through 𝛾 , 𝛾 , . . . , 𝛾 𝑠 ≡ 𝛾 , 𝛾 𝑖 being a site on the lattice for 𝑖 = , . . . , 𝑠 .We consider contractible cycles only. We say that a cycle is contractible if it does not wind around the torus. Fixing 𝑘 < 𝑠 , inour approach, we map a portion of the cyles into the standard chordal geometry as follows.1. By means of a translation, rotation, and scale transformation that we denote 𝜎 , we map Γ into a curve in C such that 𝜎 ( 𝛾 ) = 𝜎 ( 𝛾 ) = − C of the real segment [− , ] to the complement in C of the unit disc, by 𝜙 ( 𝑧 ) = 𝑧 +√ 𝑧 + √ 𝑧 − 𝜄 ( 𝑧 ) = / 𝑧 .4. We map the interior of the disc to H , by a Möbius transformation 𝜇 ( 𝑧 ) = 𝑖 + 𝑧 − 𝑧 . (B1)5. Given now the curve ˜ Γ = Φ ◦ Γ , where Φ (cid:66) 𝜇 ◦ 𝜄 ◦ 𝜙 ◦ 𝜎 , we perform 𝑘 steps of the standard zipper algorithm in H . Thisuniformizes the rst 𝑘 − Γ , mapping each of the 𝛾 𝑖 for 𝑖 = , . . . , 𝑘 − Γ 𝑘 = ( 𝛾 𝑘 , 𝛾 𝑘 + , . . . , 𝛾 𝑠 ) ⊂ Γ (a chordal curve from the point 𝛾 𝑘 to the point 𝛾 in the original domain), is mapped in a curvein the upper half plane.6. Finally, the zipper algorithm, applied to the remaining 𝑠 − 𝑘 − Γ , gives the set of pairs ( 𝑡 𝑖 , 𝜉 𝑡 𝑖 ) 𝑠 − 𝑘 − 𝑖 = , with 𝜉 𝑡 ≡ 𝜉 𝑡 𝑖 is the driving function at the “time” 𝑡 𝑖 extracted by the algorithm.Pictorially, f i tt ed κ t max k=60k=90k=120 1.9 2 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 0.002 0.004 0.006 0.008 0.01 0.012 0.014 1.9 2 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 0.002 0.004 0.006 0.008 0.01 0.012 0.014Triangular lattice f i tt ed κ t max k=60k=90k=120 1.9 2 2.1 2.2 2.3 2.4 2.5 2.6 2.7 2.8 0.002 0.004 0.006 0.008 0.01 0.012 0.014 Figure 4. Fitted values of 𝜅 from simulations on the square lattice (top) and triangular lattice (bottom) as functions of the upper cuto 𝑡 max ;the three symbols refer to the three values of 𝑘 considered. Shaded bands indicate our nal estimates (dashed line) and their variability. − − −
20 0 20 − − 𝛾 𝛾 Γ − −
50 0 50 − 𝜎 ( 𝛾 ) 𝜎 ( 𝛾 ) 𝜎 ◦ Γ − −
100 0 100 − 𝜙 ◦ 𝜎 ◦ Γ − − . . − − . . 𝜄 ◦ 𝜙 ◦ 𝜎 ◦ Γ − − . . . . 𝜇 ◦ 𝜄 ◦ 𝜙 ◦ 𝜎 ◦ Γ ≡ Φ ◦ Γ . . . . 𝑘 -step zipper − − −
20 0 20 − − 𝛾 𝛾 Γ − −
50 0 50 − 𝜎 ( 𝛾 ) 𝜎 ( 𝛾 ) 𝜎 ◦ Γ − −
100 0 100 − 𝜙 ◦ 𝜎 ◦ Γ − − . . − − . . 𝜄 ◦ 𝜙 ◦ 𝜎 ◦ Γ − − . . . . 𝜇 ◦ 𝜄 ◦ 𝜙 ◦ 𝜎 ◦ Γ ≡ Φ ◦ Γ . . . . 𝑘 -step zipper − − −
20 0 20 − − 𝛾 𝛾 Γ − −
50 0 50 − 𝜎 ( 𝛾 ) 𝜎 ( 𝛾 ) 𝜎 ◦ Γ − −
100 0 100 − 𝜙 ◦ 𝜎 ◦ Γ − − . . − − . . 𝜄 ◦ 𝜙 ◦ 𝜎 ◦ Γ − − . . . . 𝜇 ◦ 𝜄 ◦ 𝜙 ◦ 𝜎 ◦ Γ ≡ Φ ◦ Γ . . . . 𝑘 -step zipper − − −
20 0 20 − − 𝛾 𝛾 Γ − −
50 0 50 − 𝜎 ( 𝛾 ) 𝜎 ( 𝛾 ) 𝜎 ◦ Γ − −
100 0 100 − 𝜙 ◦ 𝜎 ◦ Γ − − . . − − . . 𝜄 ◦ 𝜙 ◦ 𝜎 ◦ Γ − − . . . . 𝜇 ◦ 𝜄 ◦ 𝜙 ◦ 𝜎 ◦ Γ ≡ Φ ◦ Γ . . . . 𝑘 -step zipper b. Numerical analysis of the driving function. We x 𝑘 = , ,
120 and consider 𝐿 = , 𝑠 > 𝑠 min = 𝑘 +
10 so as to ensure that 𝑘 < 𝑠 and that there at least a few points in the remaining curve; however, wechecked that the results do not depend appreciably on the particular choice of 𝑠 min . For each pair 𝐿, 𝑘 , we measure the 𝑠 − 𝑘 pairs ( 𝑡 𝑖 , 𝜉 𝑡 𝑖 ) 𝑠 − 𝑘 − 𝑖 = as described above, for ∼ independent realizations.We rst performed our analysis on the Q model. The leftmost panels in Fig. 5a show the ensemble-averaged mean squaredisplacement of 𝜉 𝑡 as a function of time. For small times, up to 𝑡 max ≈ . (cid:10) 𝜉 𝑡 (cid:11) is approximately linear in 𝑡 . While thetime 𝑡 max does not seem to change appreciably with lattice size, the average number of steps needed to reach it, (cid:104) 𝑠 ( 𝑡 max )(cid:105) ,increases with 𝐿 and with 𝑘 . Note that the quantity (cid:104) 𝑠 ( 𝑡 max )(cid:105) is computed without taking into consideration the rst 𝑘 stepsof the original curve, which get mapped to the boundary of the domain. For 𝐿 =
200 we nd (cid:104) 𝑠 ( 𝑡 max )(cid:105) ≈ / /
188 for 𝑘 = / /
120 respectively; for 𝐿 =
400 we nd (cid:104) 𝑠 ( 𝑡 max )(cid:105) ≈ / / (cid:104) 𝑠 ( 𝑡 max )(cid:105) and the averagenumber of steps of the curves are 0 . , . , .
51 for 𝐿 =
200 (again for 𝑘 = , , . , . , .
61 for 𝐿 = 𝜅 we performed linear ts of (cid:10) 𝜉 𝑡 (cid:11) versus 𝑡 , using only the data at the larger size, 𝐿 = [ , 𝑡 max ] with varying 𝑡 max . Fit results are shown in Fig. 4. Let us rst focus on the square lattice. The tted valuesfor the square lattice stabilize around 𝑡 max ≈ .
01. The three values of 𝑘 that we considered yield compatible estimates at theplateau. We estimate 𝜅 = . ( ) .Besides the scaling of the mean squared displacement, SLE predicts that the normalized process 𝜉 𝑡 /√ 𝑡 is distributed normallywith mean 0 and variance 𝜅 . To have sucient statistics, we considered the data for 𝑘 = , ,
120 combined, and all times from 𝑡 min = .
001 and 𝑡 max = .
012 (for smaller values of 𝑡 min lattice artefacts become apparent, aecting the tails of the distribution).Fig. 5a shows that the normalized process is approximately Gaussian as expected. The small discrepancies, especially apparentin linear scale for 𝐿 = 𝜅 = . ( ) , see also Fig. 5b. However, notice that, in the T model,the three values of 𝑘 give slightly inconsistent estimates, see Fig. 4. The values for 𝑘 =
60 do not reach a well dened plateau;excluding them from the overall estimate gives 𝜅 = . ( ) . ξ t2 > ~ κ tL=400 < ξ t > t 0.10.20.3 -4 -3 -2 -1 0 1 2 3 4N(0, κ ) P D F ξ t / t -3 -2 -1 P D F ( l og sc a l e ) ξ t / t k=60k=90k=120 0 0.01 0.02 0.03 0.04 0.05 0.06 0 0.01 0.02< ξ t2 > ~ κ tL=200 < ξ t > t 0.10.20.3 -4 -3 -2 -1 0 1 2 3 4N(0, κ ) P D F ξ t / t -3 -2 -1 P D F ( l og sc a l e ) ξ t / t k=60k=90k=120 (a) Numerical analysis of the driving function in the Q model. The top panels refer to 𝐿 = 𝐿 = ξ t2 > ~ κ tL=400 < ξ t > t 0.10.20.3 -4 -3 -2 -1 0 1 2 3 4N(0, κ ) P D F ξ t / t -3 -2 -1 P D F ( l og sc a l e ) ξ t / t k=60k=90k=120 (b) Numerical analysis of the driving function in the T model. Figure 5. Driving function analysis on the Q model (a) and on the T model (b).
Leftmost panels.
Mean square displacement as a function oftime, binned in intervals 𝛿𝑡 = . 𝑘 (the number of steps of the walks that get mappedto the boundary of H ): 𝑘 =
60 (squares), 𝑘 =
90 (circles), 𝑘 =
120 (diamonds). Dashed lines are the SLE prediction for 𝜅 = . Central andrightmost panels.
Empirical probability distribution function of the rescaled process 𝜉 𝑡 /√ 𝜅𝑡 in linear and logarithmic scale (center and right,respectively). Dashed lines are the SLE prediction for 𝜅 = .
07, i.e., Gaussian distributions with mean 0 and variance 𝜅𝜅