Percolation thresholds on high dimensional D_n and dense packing lattices
PPercolation thresholds on high dimensional D n and dense packing lattices Yi Hu and Patrick Charbonneau
1, 2 Department of Chemistry, Duke University, Durham, North Carolina 27708, USA Department of Physics, Duke University, Durham, North Carolina 27708, USA (Dated: February 22, 2021)The site and bond percolation problems are conventionally studied on (hyper)cubic lattices, whichafford straightforward numerical treatments. The recent implementation of efficient simulation algo-rithms for high-dimensional systems now also facilitates the study of D n root lattices in n dimensionas well as E -related dense packing lattices. Here, we consider the percolation problem on D n for n = 3 to 13 and on E relatives for n = 6 to 9. Precise estimates for both site and bond per-colation thresholds obtained from invasion percolation simulations are compared with dimensionalseries expansion on D n lattices based on lattice animal enumeration. As expected, the bond perco-lation threshold rapidly approaches the Bethe lattice limit as n increases for these high-connectivitylattices. Corrections, however, exhibit clear yet unexplained trends. Interestingly, the finite-sizescaling exponent for invasion percolation is found to be lattice and percolation-type specific. I. INTRODUCTION
Percolation being one of the simplest critical phenom-ena, its models play particularly important roles in sta-tistical physics [1]. Minimal models – lattice-based ones,in particular – have thus long been used to test notionsof universality as well as the relationship between mean-field and renormalization group predictions. On lattices,two covering fractions p can be defined: (i) the probabil-ity that a vertex is occupied, and (ii) the probability thatan edge between nearest-neighbor vertices is occupied.As p increases, a percolating cluster forms at a threshold p sitec or p bondc , depending on the covering choice [1]. Be-cause precise thresholds values are prerequisite for strin-gently assessing criticality [2–5], (and because thresholdsare lattice specific and lack an analytical expression [6],)substantial efforts have been directed at estimating themthrough numerical simulations [2, 7–10] and graph-basedpolynomial methods [11–14]. The strong dependence ofcriticality on spatial dimension n , especially above andbelow its upper critical dimension, n u = 6, motivates ex-panding these efforts over an extended range of n [1, 15].In this context, the invasion percolation algorithm re-cently introduced by Mertens and Moore [2, 16] is partic-ularly interesting. In short, the algorithm directly growsa percolating cluster, and thus provides both the univer-sal asymptotic critical behavior and the lattice-specificfinite-size scaling correction. Most crucially, by avoid-ing the explicit construction of a lattice grid, the schemepreserves a polynomial space complexity as n increases.Threshold values up to ten significant digits of precisionhave thus been obtained on hypercubic lattices ( Z n ) upto n = 13 [17].Hypercubic lattices, although geometrically straight-forward, are in some ways not natural systems to consideras dimension increases. Recall that lattices can be seenas discretizations of Euclidean space R n in which eachlattice point is a vertex of a cell in that tessellation. As n increases, the cubic cells that tile Z n become increas-ingly dominated by spiky corner sites. The cells of root lattices, D n , are relatively less spiky in n ≥
3. A wayto quantify this effect is to compare the sphere packingfraction for different lattices. In this sense, D n packingsare exponentially denser than their Z n counterpart bya factor of 2 n/ − [18]. Similarly, the eight-dimensional E lattice corresponds to a sphere packing fraction twicethat of D (16 times that of Z ); E -related lattices, E , E and Λ are also the densest known sphere packingsin their corresponding dimension (see Appendix A). Thisadvantage has motivated the recent consideration of D n and E -related periodic boundary conditions for high-dimensional numerical simulations [19, 20]. For a samecomputational cost, these periodic boxes have a larger in-scribed radius than hypercubes and thus present less pro-nounced finite-size corrections. Considering these latticesmay thus help suppress obfuscating pre-asymptotic cor-rections to percolation criticality [2, 4], especially aroundthe upper critical dimension n u . Further interest in D n lattices also stems from its inclusion of the canonicalthree-dimensional face-centered cubic lattice, D .In this work, we investigate the two canonical latticepercolation thresholds in D n for n = 3 to 13 as well as for E -related in n = 6 ∼
9. In Section II, we first derive theseries expansion for both p sitec and p bondc on D n latticesbased on lattice animal enumeration. We then describethe invasion percolation algorithm in Section III, and an-alyze the numerical percolation results in Section IV. Webriefly conclude in Section V. II. SERIES EXPANSION
In this section we derive high-dimensional series expan-sions for both site and bond percolation thresholds on D n lattices by counting lattice animals embedded on the lat-tice [17, 21]. For site percolation a site animal of size v isa cluster of v lattice vertices connected after connectingall neighboring vertex pairs. Similarly, for bond percola-tion a bond animal of size e consists of a connected setof e lattice edges. In both cases, the perimeter t is thenumber of incident vertices (or edges) for the lattice ani- a r X i v : . [ c ond - m a t . s t a t - m ec h ] F e b mal but not part of it. Two lattice animals are distinct ifthey do not overlap through translation. We here denotethe number of site and bond animals of perimeter t ona n -dimensional lattice as g tv ( n ) and g te ( n ), respectively,where both t and g are functions (polynomials in D n lat-tice) in terms of n with a functional mapping t (cid:55)→ g . A. Site percolation
We first consider the site percolation threshold usingsite animals. Following Mertens et al. [17] we define thepolynomial A v ( q ) = (cid:88) t g tv q t , (1)in terms of q = 1 − p . In particular, A v (1) ≡ A v givesthe total number of lattice animals of size v in an n -dimensional lattice. At covering fraction p , the expectedsite cluster size on the lattice is then S = (cid:88) v v p v − A v (1 − p ) ≡ ∞ (cid:88) (cid:96) =0 b (cid:96) ( n ) p (cid:96) , (2)where we have expanded S as a power series in p . Because A v ( q ) has a factor of p v − , obtaining b (cid:96) only requires A , ..., A (cid:96) +1 , i.e., counting g t to g t(cid:96) +1 . Once these termsare known, p sitec can be approximated by re-summing theterms using a Pad´e approximant b (cid:96) − /b (cid:96) .The objective is thus to count g tv and express it asa polynomial in n . On hypercubic lattices the compu-tational cost of this enumeration is greatly simplifiedby the introduction of proper dimension for lattice an-imals [17, 21], but this approach is not obviously gener-alizable for D n lattices. We instead implement a moregeneric, brute-force algorithm [22], which traverses everypossible lattice animal via a breadth first search (BFS)of the lattice vertices.Starting at the origin, we add every nearest-neighborsite (as described in Appendix A) to the perimeter set.In that set, we then choose one site and add it to the siteanimal set according to the following criteria:1. if the coordinates lexicographically greater than theorigin;2. if the site is newly added to the perimeter set atthe previous iteration, or lexicographically greaterthan all sites in the site animal set.These two conditions guarantee that a site animal – af-ter properly accounting for translational invariance – iscounted exactly once. Once a new site is selected, theperimeter set is updated with nearest neighbors of thissite, with new sites being selected until the pre-assignedsize v is reached. The perimeter t of each generated lat-tice animal is also calculated. Therefore, by running thealgorithm once with assigned v and n , a series of integervalues of ( t ( n ) , g tv ( n )) can be obtained. The next step entails obtaining the analytical polyno-mial form, t ( n ) and g tv ( n ). We first consider t . Knowingthat the polynomial form of t ( n ) has the same leading or-der as the vertex connectivity, i.e., t = 2 n ( n −
1) + O ( n )which is quadratic, we can relate the site animals in dif-ferent dimensions by a linear fit of t , t = 2 n ( n − v + c t n + c t . (3)The coefficients c t and c t can then be extracted with t results from two different dimensions.The polynomial g tv ( n ) is also obtained by solving alinear system. The (upper bound of the) order ofthis polynomial must, however, be determined in ad-vance. Because the total number of lattice animals is ∼ [2 n ( n − v − , the order of g tv ( n ) is also at most n v − . And because the orientational degeneracy un-der D n symmetry requires that g tv ( n ) always has roots n ( n − n v − . There-fore, we require the numerical g tv ( n ) results for at most2 v − g tv ( n ) / [ n ( n − v − (cid:88) i =0 c gi n i , (4)to obtain c gi . Results for t and g polynomials are avail-able in Ref. 23. While the validity of this fitting form hasyet to be mathematically demonstrated, the correctnessof g tv polynomials can be empirically tested by checkingthat the residual vanishes when fitting the results of a(larger-than-necessary) number of dimensions. We haveevaluated site animals up to dimension n = 15, which issufficient to solve Eq. (4) in v ≤
6. However, because thetotal number of site animals, A v ∼ n v − , grows expo-nentially with v , obtaining results for v > t and g tv with v ≤
6, we obtain the first sixterms in the expansion for S (Eq. (2)) b = 1 ,b = 2 n ( n − ,b = 2 n ( n − n − n + 7) ,b = 2 n ( n − n − n + 57 n − n + 12) ,b = 2 n ( n − n − n + 272 n − n + 804 n − n + 857) ,b = 2 n ( n − n − n + 1004 n − n +6018 n − n − n + 2840756 n − . (5)The approximant b /b agrees with the threshold of theBethe lattice, i.e., a branching tree of degree z = 2 n ( n − p c , Bethe = 1 z − ≡ σ . (6)In the following we denote 1 /σ the Bethe lattice limit ofthe percolation threshold on a lattice.In general, b (cid:96) has a leading order of n (cid:96) , and b (cid:96) − /b (cid:96) provides an approximation for p c with an error that van-ishes asymptotically as O ( n − ( (cid:96) +2) ). For comparison, b (cid:96) ∼ n (cid:96) for a hypercubic lattice and the approximant b (cid:96) /b (cid:96) +1 have the same order of error O ( n − ( (cid:96) +2) ) [17].Expanding b /b , in particular, gives p sitec = 1 σ + 1 n + 238 n + 172 n + 99932 n + O ( n − ) . (7)The accuracy of this series is evaluated in Sec. IV B. B. Bond percolation
For the bond percolation, we similarly define the bondpolynomial A e ( q ) = (cid:88) t g te q t , (8)which gives the expected bond cluster size S = (cid:88) e e p e − A e (1 − p ) ≡ ∞ (cid:88) (cid:96) =0 b (cid:96) ( n ) p (cid:96) (9)as a polynomial in p . The enumeration scheme for bondanimals is essentially the same as for site animals, withthe exception that we now maintain bonds, which are in-dexed as the coordinates of the lexicographically smallervertex on this bond, in addition to the orientation index– from 1 to n ( n −
1) – of the bond. The bond animalenumeration is then used to obtain a series of numericalvalues ( t ( n ) , g te ( n )). The perimeter polynomial t ( n ) forbond animal is also quadratic with n , but the leadingprefactor is not fixed. The polynomial is thus obtainedby fitting t ( n ) in at least three dimensions. A bond ani-mal of size e includes at most e + 1 sites, hence the orderof g te ( n ) is at most n e , including roots n ( n − e − g te ( n ), similar to Eq. (4), g te ( n ) / [ n ( n − e − (cid:88) i =0 c gi n i . (10)Bond animals can thus be evaluated up to dimension n =12 and Eq. (10) can be solved up to e = 5. Results for t and g polynomials are also available in Ref. 23. Here aswell, because the total number of bond animals A e ∼ n e grows exponentially with e , results for e > b = n ( n − ,b = 2 n ( n − n − n − ,b = 2 n ( n − n − n + 9) ,b = 2 n ( n − n − n + 12 n − n + 27 n + 131 n − ,b = 2 n ( n − n − n + 64 n − n + 56 n + 328 n + 1534 n − n + 7499) . (11) Note that we have b (cid:96) ∼ O ( n (cid:96) +2 ) which is two orders(in n ) higher than for site percolation. Note also thatunlike for site percolation, b /b here has yet to convergeto the Bethe lattice limit at leading order. For (cid:96) ≥ p bondc ≈ b (cid:96) − /b (cid:96) has anerror of O ( n − ( (cid:96) +3) ) one order smaller than for p sitec . Inparticular, expanding b /b gives p bondc = 1 σ + 1 n + 8116 n + O ( n − ) . (12)The accuracy of this series is also evaluated in Sec. IV B. III. INVASION PERCOLATION
In this section we briefly describe the invasion percola-tion algorithm by Mertens and Moore [2] (derived fromRef. 24) for an arbitrary lattice structure, and then ana-lyze its complexity for the considered lattices.As stated in the introduction, the algorithm grows asingle cluster without explicitly storing the lattice grid.Two data structures are then used: a set (collection ofunique elements) S to maintain all sites (or bonds) be-longing to the cluster as well as those incident to them;and a priority queue Q for the stepwise growth of thecluster. For site percolation, starting from the origin ev-ery neighboring vertex is inserted (following Appendix A)into S . For each of these new vertices, a random weight w i ∈ [0 ,
1) is assigned and the vertex is inserted into Q with w i as the key. For the next step, the vertex of min-imum weight in Q is popped, incrementing the clustersize N . The previous steps are repeated until the pre-assigned cluster size N = N is attained. The expectedset size at a certain N , denoted B ( N ) = (cid:104)|S ( N ) |(cid:105) , iscomputed by averaging the set size among independentrealizations. For the bond percolation, we start with anarbitrary bond incident to the origin and otherwise thesame procedure is used.The cluster obtained by invasion percolation process si-multaneously approaches the giant component at p c withthe scaling form [2, 16] NB ( N ) ≈ p c (1 − cN − δ ) (13)where δ is the correction exponent and c a fitting con-stant.For each instance, the space complexity is n |S| + |Q| ∼ nN/p c ∼ O ( nσN )where the factor of n accounts for the size of an n -dimensional vector. For D n lattice the space complex-ity is thus O ( n N ). Although the space complexity islarger, by a factor of n , than for Z n lattices ( O ( n N )),the memory requirement is still moderate for contem-porary computers. The time complexity depends onthe complexity of the insertion to Q which is at most -8 -7 -6 -5 -4 -3 FIG. 1. Convergence of
N/B ( N ) to site (diamonds) and bond(asterisks) percolation thresholds on D n lattices in n = 4(blue), 8 (red) and 12 (yellow). Finite-size correction scalesas O ( N − δ ) (solid lines), where δ → n → ∞ . O ( n + log |Q| ) ≈ O (log N ), and thus O ( N log N ) in to-tal. In practice, we can grow clusters up to N = 1 . × in n = 3 and up to 2 × in n = 13, within a memory us-age of less than 10 GB. At least 10 independent clustersare obtained for each lattice, which results in each real-ization is usually taking less than a minute on an AMDRyzen 3900x processor. IV. RESULT AND DISCUSSION
In this section we compare the numerical thresholdvalues obtained from the invasion percolation describedin Sec. III with the series expansion results obtained inSec. II.
A. Numerical thresholds
Table I reports both site and bond percolation thresh-olds for D n as well as for E -related lattices obtained byfitting the numerical N/B ( N ) results with Eq. 13. Thevalues for n = 3-6 are consistent with published values,and, except for n = 3, our results are at least an orderof magnitude more accurate. For 6 < n ≤
13 no priorresult is known. As in Ref. 2 for Z n lattices, p c results for D n lattices are obtained with higher precision – for com-parable computational efforts – as n increases (Fig. 1).Because the correction exponent δ (Eq. (13)) increaseswith n , finite-size corrections then decay faster. For therange of n considered, this advantage compensates forthe decrease in N imposed by the growing memory cost.As a result, a relative uncertainty of 10 − to 10 − is ob-tained for all investigated dimensions.Because δ controls the convergence rate of invasion TABLE I. Site and bond percolation threshold on D n and E related latticesLattice p sitec p bondc D D D D D D D D D D D E E E FIG. 2. Finite-size scaling exponent δ for site and bond per-colation on Z n , D n and E -related lattices, E , E , E andΛ . (Results for Z n site percolation in n = 4 to 13 are fromRef. 2.) Error bars from fitting are smaller than (comparableto) the marker size in n ≤ n > δ generically grows with n but its value isnot universal. percolation, it is interesting to compare its behavior fordifferent lattices. As a first glance, δ increases with n for both Z n and D n lattices and tends to 1 as dimen-sion increases, as expected from the Bethe lattice analy-sis [2, 16]. While for site percolation on Z n , D n and E -related lattices δ appears similar, the exponent evolvesdifferently for bond percolation on different lattices aswell as for either type of percolation on a same lattice.Because the exact value of δ depends on the type of per-colation as well as on lattice geometry, we conclude thatthe exponent is not universal. As a corollary, δ may bea useful quantity for selecting a lattice for studying criti-cality; a greater δ indeed implies a faster decay of certainfinite-size corrections. B. Comparison with series expansion
Our precise numerical thresholds for D n lattices canbe compared with the series prediction obtained for bothsite and bond percolation in Sec. II. The relative error ofthe expansion up to n − (cid:96) term, defined as η ( (cid:96) ) p = (cid:12)(cid:12)(cid:12) p c , simulation − p ( (cid:96) )c , series (cid:12)(cid:12)(cid:12) /p c , simulation , (14)is shown in Fig. 3(a). As expected, these thresholds con-verges gradually to the Bethe lattice value, 1 /σ , in thelarge n limit. For site percolation, this convergence rateis fairly slow – a ∼
10% deviation persists even in n = 13– but introducing higher-order terms in the series dra-matically reduces that error. Including terms of order upto n − leads to a relative error of ∼ .
1% in n = 13. Forbond percolation, because the prefactors for both n − and n − in the expansion form are zero, the deviation isalready down to ∼ .
1% in n = 13. Including two moreterms in Eq. (12) further divides the error by a factor ∼ n . The series expansion in Eqs. (7) and (12) is thusexpected to predict percolation thresholds with very highaccuracy for n > Z n , D n and E -related lat-tices are compared with the Bethe lattice result inEq. (6). (Although a dimensional series expansion isnot available for E -related lattices, the site connec-tivity, z = 72 , , ,
272 for E , E , E and Λ lat-tices [18], respectively, alone suffices for this comparison(see Appendix A).) In all three cases, the Bethe latticeprediction better matches the bond than the site per-colation threshold (Fig. 3(b)). For Z n and D n latticesthis result is expected from the series expansion. Inthe large n limit, the deviation of p bondc from the Bethelattice limit is of O ( n − ) , O ( n − ) , O ( n − ) and O ( n − )for p sitec ( Z n ) , p bondc ( Z n ) , p sitec ( D n ) and p bondc ( D n ), respec-tively. For E -related lattices, for which no such seriesexist, the same trend is observed. More specifically, thedeviation is < ∼
1% for bond percolation and < ∼
40% forsite percolation. This concordance suggests that the ef-fect might be more than a mere coincidence. Yet it lacksa physical explanation. A generic scaling form for thepercolation threshold beyond the Bethe lattice approxi-mation might be informative in this respect, but is stillfound lacking. -4 -3 -2 -1 -3 -2 -1 (a)(b) FIG. 3. (a) Relative error for the site (diamonds) and bond(asterisks) percolation thresholds on D n lattices predicted byseries expansion for various highest-order terms. Note that forsite percolation the high-order lines are truncated in small n because the relative error then changes sign. Lines are guidesto the eye. (b) Percolation thresholds on Z n , D n and E -related lattices (markers with dotted line) compared to theBethe lattice limit 1 /σ (solid line), which matches the bondpercolation threshold well in all three lattice types. V. CONCLUSION
We have reported the series expansion and numericalpercolation thresholds for D n lattices as well as the nu-merical thresholds for E -related lattices from n = 6 to9. The excellent agreement between the two independentapproaches cross-validates their results. Remarkably,bond percolation presents much faster decaying finite-size corrections than site percolation for invasion percola-tion in n >
6. This finding suggests that pre-asymptoticcorrections might be most efficiently suppressed in theformer. The Bethe lattice approximation to the perco-lation threshold also presents a markedly higher preci-sion for bond percolation than for site percolation for D n , due to the vanishing of the first subleading ordercoefficients in the series expansion. This feature thus ap-pears to be generic for lattices other than Z n , for whichit was first reported [17, 26, 27]. Our finding identifyunresolved features of percolation and set the stage forinvestigating percolation criticality on high-dimensionallattices beyond the conventional hypercubic geometry. ACKNOWLEDGMENTS
We thank R. M. Ziff for carefully maintaining the Per-colation threshold Wikipedia page, which has greatly fa-cilitated our literature search. This work was supportedby a grant from the Simons Foundation (
Appendix A: Lattice packing
In this appendix we briefly review the structure of thehigh-dimensional lattices considered in this study, follow-ing the construction in Ref. 30. As reference, the con-ventional n -dimensional hypercubic lattices, Z n , is de-fined as a set of n -dimensional vectors of integer compo-nents. The nearest-neighbor vector ( a , a , ..., a n ) in Z n are ( ± , n − ) (this notation means ( ± , , ... ) and theirpermutations). The number of nearest neighbors (kissingnumber) is thus 2 n . D n lattices can be viewed as a sub-set of Z n in which the coordinates have even sum. Thenearest neighbor vectors are ( ± , n − ), thus resultingin 2 n ( n −
1) nearest neighbors in total. In n = 3, for ex-ample, the 12 nearest-neighbor vectors for the D ≡ fcc lattice read(1 , , , (1 , − , , ( − , , , ( − , − , , (1 , , , (1 , , − , ( − , , , ( − , , − , (0 , , , (0 , , − , (0 , − , , (0 , − , − .D , D and D lattices are the densest packings ofequal spheres in the corresponding dimensions. Thedensest sphere packings for n = 6 to 9 are E , E , E and Λ lattices, respectively. In particular, the E lat-tice consists of two D lattice points with offset (
12 8 ).The nearest-neighbor vectors of E can be viewed as fourgroups, ± (
12 8 ) , , (
12 4 , −
12 4 ) ,
70 vectors , (
12 2 , −
12 6 ) and ( −
12 2 ,
12 6 ) ,
56 vectors , ( ± , ) ,
112 vectors , (A1)and thus each vertex has 240 nearest neighbors in total. E lattice is a cross-section of E in n = 7. One of thechoices to generate nearest neighbor vectors in E is thesubset of Eq. (A1) with zero sum, which results in 126vectors in total. Further constraining a + a = 0 leadsto 72 nearest neighbor vectors in E . The Λ lattice isnot unique, but one of its forms can be constructed sim-ilarly to E . It consists of two D lattice points, offsetby ( , ..., , d = 10 the (presumed) densest packingis a non-lattice [31], and thus n = 9 offers a natural endto our consideration of dense packing lattices.Finally, we note that in some dimensions there existstructures comparable to lattices considered here. Forexample, in n = 3 face-centered cubic ( D ) is stronglyrelated to the hexagonal closed-packed structure. Asa result, their site percolation thresholds are close al-though not identical [32]. In n = 5 , ,
7, four similarstructures are known [31] and their p c values may alsodiffer marginally. In n = 9 a continuum of structures canbe constructed. In the current study these alternativelattices were not considered, hence their constructionsare omitted from this appendix. [1] D. Stauffer and A. Aharony, Introduction To PercolationTheory (Taylor & Francis, 1994).[2] S. Mertens and C. Moore, Phys. Rev. E , 022120(2018).[3] W. Huang, P. Hou, J. Wang, R. M. Ziff, and Y. Deng,Phys. Rev. E , 022107 (2018).[4] G. Biroli, P. Charbonneau, and Y. Hu, Phys. Rev. E ,022118 (2019).[5] Z. Xun and R. M. Ziff, Phys. Rev. Research , 013067(2020).[6] J. C. Wierman, Phys. Rev. E , 027105 (2002).[7] C. D. Lorenz and R. M. Ziff, Phys. Rev. E , 230 (1998). [8] X. Xu, J. Wang, J.-P. Lv, and Y. Deng, Front. Phys. ,113 (2014).[9] M. Kotwica, P. Gronek, and K. Malarz, Int. J. Mod.Phys. C , 1950055 (2019).[10] Z. Xun and R. M. Ziff, Phys. Rev. E , 012102 (2020).[11] C. R. Scullard and R. M. Ziff, Phys. Rev. Lett. ,185701 (2008).[12] C. R. Scullard and R. M. Ziff, J. Stat. Mech. Theory Exp. , P03021 (2010).[13] J. L. Jacobsen, J. Phys. A , 454003 (2015).[14] C. R. Scullard and J. L. Jacobsen, Phys. Rev. Research , 012050 (2020). [15] S. Kirkpatrick, Phys. Rev. Lett. , 69 (1976).[16] S. Mertens and C. Moore, Phys. Rev. E , 042116(2017).[17] S. Mertens and C. Moore, J. Phys. A , 475001 (2018).[18] J. H. Conway and N. J. A. Sloane, “Certain importantlattices and their properties,” in Sphere Packings, Lat-tices and Groups (Springer New York, New York, NY,1988) pp. 94–135.[19] L. Berthier, P. Charbonneau, and J. Kundu, Phys. Rev.Lett. , 108001 (2020).[20] G. Biroli, P. Charbonneau, E. I. Corwin, Y. Hu, H. Ikeda,G. Szamel, and F. Zamponi, arXiv preprint (2020),2003.11179.[21] W. F. Lunnon, Comput. J , 366 (1975).[22] S. Mertens, J. Stat. Phys. , 1095 (1990).[23] “Duke digital repository,”https://doi.org/10.7924/xxxxxxxxx.[24] D. Wilkinson and J. F. Willemsen, J. Phys. A , 3365(1983). [25] S. C. Van der Marck, Int. J. Mod. Phys. C , 529 (1998).[26] D. Gaunt, M. Sykes, and H. Ruskin, J. Phys. A , 1899(1976).[27] D. Gaunt and H. Ruskin, J. Phys. A , 1369 (1978).[28] R. Pordes, D. Petravick, B. Kramer, D. Olson,M. Livny, A. Roy, P. Avery, K. Blackburn, T. Wenaus,F. W¨urthwein, I. Foster, R. Gardner, M. Wilde, A. Blate-cky, J. McGee, and R. Quick, in J. Phys. Conf. Ser. , 78,Vol. 78 (2007) p. 012057.[29] I. Sfiligoi, D. C. Bradley, B. Holzman, P. Mhashilkar,S. Padhi, and F. Wurthwein, in , 2, Vol. 2 (2009) pp. 428–432.[30] J. H. Convay and N. J. A. Sloane, IEEE Trans. Inf. The-ory , 227 (1982).[31] J. H. Conway and N. J. A. Sloane, Discrete Comput.Geom. , 383 (1995).[32] C. D. Lorenz, R. May, and R. M. Ziff, J. Stat. Phys.98