On the universality of crossing probabilities in two-dimensional percolation^{†}

R. P. Langlands, C. Pichet, Ph. Pouliot, Y. Saint-Aubin

† Appeared in Jour. Stat. Phys., vol 67, Nos. 3/4.

1. INTRODUCTION

The geometrical explanation of universality in terms of fixed points of renormalization-group transformations has met with enormous success, but its lack of precision continues to present a challenge to the mathematician, even for relatively simple models, such as percolation in two dimensions.

To begin to reflect on the problem, one can assume that the crossing probabilities (cf. [5]) at the critical probability are universal and can therefore serve as coordinates of the fixed point. The current climate imposes a more critical stance. Indeed comments of several physicists and mathematical physicists have made it clear that this hypothesis was not universally shared, or at the very least that its possible significance was not widely appreciated. We could find no mention of it, or its simple consequences, in standard texts ([4], [8]).

It is not our intention here to comment further on the initial reflections, which will be developed further, both numerically and theoretically, elsewhere. It seemed worthwhile none the less, in view of the disparate views encountered, to examine the hypothesis itself numerically in order to establish a concrete basis for confidence in the usefulness of the crossing probabilities. Although very crude in comparison with many of the numerical results on percolation, the evidence that it is the purpose of this paper to present establishes conclusively that the crossing probabilities are universal, and therefore suitable coordinates for the fixed point, and that several basic models, to be described later, fall into the same universality class.

The mathematical consequence is that attention is focussed not on the critical indices, which are from a mathematical viewpoint both literally and figuratively derived objects, since they are given hypothetically by eigenvalues of the jacobian matrix of the renormalization group at its fixed point, but on an object with a more direct mathematical significance, the fixed point itself. The advantage resulting from the change of focus is of course not only mathematical. Since the new object is of lower order it is in most respects of easier numerical access, and the authors, by no means specialists in simulation, have therefore imposed upon themselves standards other than those appropriate for the calculation of critical exponents. The casual, implicit reference to η_{||} at the end of the paper notwithstanding, we shall not be calculating critical indices. We will be comparing functions, and this creates different problems of accuracy. We have tried, in what seems to us an appropriate substitute for the usual error analysis, to explain clearly in Section 2 and 3 the sources of inevitable errors, and to estimate their magnitude.

It is best to formulate the questions not as they first presented themselves, but in the more cogent manner suggested by our experience. Consider, for the sake of precision, percolation by sites or by bonds on a lattice at critical probability, and place on this lattice a large rectangle with sides parallel to the two axes. Take its width to be an and its height to be bn. The positive numbers a and b are fixed for the moment but n will approach infinity. The exact manner in which the rectangles grow is unimportant, but to be definite we place the lower left corner at the origin.

For a given state of occupation of the sites or bonds, the notion of a horizontal crossing (or in the language of [5], an occupied horizontal crossing) of the rectangle necessarily includes an arbitrary element, because the crossing is from a band on the left to a band on the right, but the exact prescription of the band in terms of width or other features is often unimportant. Thus the probability π_{h}^{n}(a, b) is somewhat ill-defined. None the less the limit

lim_{n→∞} π_{h}^{n}(a, b) = π_{h}(a, b) = π_{h}(a/b),

provided it exists, as we assume, is well defined and depends only on the quotient r = \frac{a}{b}.

The function π_{h}(r) is defined for 0 < r < ∞, is monotone decreasing, and approaches 1 as r approaches 0 and 0 as r approaches ∞. A similar function π_{v}(r) is defined by vertical crossings, and is again monotone, but increasing and approaches 0 as r approaches 0 and 1 as r approaches ∞. Granted the continuity of both functions there is a unique value r_{0} of r such that

π_{v}(r_{0}) = π_{h}(r_{0}).

If the lattice is symmetric with respect to permutation of the two axes, as for a square lattice with its usual orientation, then r_{0} is 1. Otherwise r_{0} is an invariant of the lattice, or more generally of the model, whose value appears to be given in the cases considered in this paper by simple formulas that can be explained heuristically, because they are immediate consequences of a symmetry that is almost certainly present, but we have not been able to prove them.

The functions π_{h}(r) and π_{v}(r) are clearly not universal, because by changing the aspect ratio of the lattice we can force π_{h}(1) to take any value between 0 and 1. Our numerical results establish, however, that the functions

η_{h}(r) = π_{h}(rr_{0}), η_{v}(r) = π_{v}(rr_{0})

are universal. We stress that the models discussed in this paper are symmetric with respect to reflections in the two coordinate axes. If this condition is not satisfied, universality continues to manifest itself but differently.

We presume, although our experiments were restricted to very few models, that the pertinent class of universality includes all those for which the assumptions of Kesten’s book [2] are valid. In particular, our conclusions apply to the probabilities π_{h}^{∗}(r) and π_{v}^{∗}(r) associated to the dual model. Since

π_{h}(r) + π_{v}^{∗}(r) = 1, π_{v}(r) + π_{h}^{∗}(r) = 1,

we conclude that r_{0}^{∗} = r_{0}.

The equation π_{v}^{∗}(r) = π_{v}(r) entailed by universality therefore implies that

π_{h}(r) + π_{v}(r) = 1,

and, as a consequence,

η_{h}(r) + η_{v}(r) = 1.

These equations are amply confirmed by our experiments.

One implication is that

π_{h}(r_{0}) = π_{v}(r_{0}) = \frac{1}{2}.

This equation is readily proved for percolation by bonds on a square lattice, but has not been proven for percolation by sites on a square lattice. It is, moreover, to our surprise, not an equation whose validity is immediately recognized, even by specialists, a strong indication that the consequences of the universality of the crossing probabilities have not always been firmly grasped.

Consider, more generally, intervals α_{1}, …, α_{m}, β_{1}, …, β_{m}, γ_{1}, …, γ_{n}, and δ_{1}, …, δ_{n} on the sides of the basic rectangle of width a and height b. We introduce the event E that on the dilated rectangle of width an and height bn there are crossings from the dilation of α_{i} to that of β_{i} for 1 ≤ i ≤ m but no crossing from the dilation of γ_{j} to that of δ_{j} for 1 ≤ j ≤ n.

It is natural to suppose once again that the limits of the probabilities

lim_{n→∞} π_{E}^{n}(a, b) = π_{E}(\frac{a}{b})

exist, and that

η_{E}(r) = π_{E}(rr_{0})

is a universal function, depending only on E. We present some evidence in support of this hypothesis, but it should be examined more extensively.

We have confined ourselves to very few events and to very few models, and have, as yet, made little attempt to examine dilations of curves other than rectangles; nor have we considered percolation in dimensions other than two. Conversations with Michael Aizenman have greatly clarified our views as to the nature of the universality manifested by the crossing probabilities, and our understanding of their invariance under various transformations of the curves defining the event E. In particular, they have suggested a number of conjectures to which we shall return in a later paper, in which the modifications required for models with less symmetry than those treated here will also be discussed.

A good deal of the work on the present paper was carried out at the Centre de recherches mathématiques of the Université de Montréal and the authors are particularly grateful to its systems manager François Lambert for the help and advice he gave them.

2. EXPERIMENTAL SETTING

The numerical evidence for the universality of crossing probabilities will be obtained on lattices of finite size. It is therefore imperative to discuss our conventions and the limitations due to finiteness. This section is devoted to these topics.

2.1 The six finite models and the four crossing functions

We have studied percolation by sites and by bonds on the three regular lattices of the plane: the square, the triangular and the hexagonal lattices. For each of these six models, 81 different values of the ratio r = \frac{a}{b} are considered. The integers a and b were so chosen that the product ab remained as close as possible to 4 × 10^{4} while the numbers ln \frac{r}{r_{0}} distributed themselves uniformly over the interval (−2, 2). The width a is, in these models, the number of sites in a line, and the height b the number of sites in a column. To avoid any confusion as to which direction is horizontal and which vertical, we have included in Figure 1 diagrams of most pertinent lattices.

As explained in the next paragraph, the values of the crossing probabilities in a finite lattice are quite sensitive to the conventions. Our conventions for percolation by sites are immediate, once we agree what points of the lattice belong to a rectangle of size (a, b), for a crossing must then join a point on the extreme left to one on the extreme right. For a square lattice the conventions are clear; for triangular and hexagonal lattices, we refer the reader to Figures 1a (triangular) and 1b. (hexagonal).

(a)

(b)

(c)

Fig. 1. (a) A triangular lattice with (width, height)=(5,4), (b) a hexagonal lattice with (width, height)=(10,5), (c) a triangular lattice with (width, height)=(4,3) for the percolation by bonds.

The conventions for the percolation by bonds are a little more intricate. We used the same dimensions as for percolation by sites on the corresponding lattice. We chose however to add all bonds attached to the sites, thereby creating spurious sites on the edges of the lattice. (For the triangular lattice, see Figure 1.c where the true sites are depicted by larger dots than the ones used for spurious sites.) Crossings are taken to start from spurious sites and to end at them.

For each of the six models and each of the 81 values of the ratio r, four crossing probabilities are to be measured. We denote the horizontal and vertical probabilities by π_{h} and π_{v}. The probability of a horizontal and a vertical crossing occurring simultaneously is denoted by π_{hv}. Finally π_{d} is the probability of a “diagonal” crossing: a diagonal crossing is a crossing starting from the upper half of the left side and reaching the right half of the bottom edge. If either the width or the height is odd, the diagonal crossing may start from a central site or end at one. The probabilities π_{hv} and π_{d} were added as examples of the events E described in the introduction.

For the models on a square lattice (with percolation by bonds or by sites), the four crossing probability functions are related by:

π_{h}(r) = π_{v}(1/r),  π_{hv}(r) = π_{hv}(1/r),  π_{d}(r) = π_{d}(1/r).

It is useful to introduce a second variable

s = ln \frac{r}{r_{0}} .

For the square lattice r_{0} = 1, and as a function of s, ln \frac{π_{h}}{π_{v}} is odd, while the functions π_{hv} and π_{d} are even. As in the introduction, we define r_{0} for percolation by sites and by bonds on the triangular and hexagonal lattices as that value of r for which π_{h}(r_{0}) = π_{v}(r_{0}). Then s = 0 at r = r_{0}. If the crossing probabilities are universal in our sense then the three functions ln \frac{π_{h}}{π_{v}}, π_{hv} and π_{d} will be symmetric with respect to the s = 0 axis. A secondary goal of our numerical work is to determine the invariants r_{0}. (See section 3.2.)

2.2 Sample size and limitations due to finiteness

In addition to the determination of a value for the critical probability, which we discuss in the following section, we have identified two difficulties in comparing the crossing probabilities for the six models: their sensitivity to the choice of conventions and the statistical errors. We discuss first the sensitivity to conventions as the sampling size we choose is partly determined by it.

To understand the sensitivity to conventions, let us consider the bond percolation on a square lattice and let us label the conventions introduced earlier as the set of rules I. Consider a second set of rules, labelled II, for which only the bonds between immediate neighbors among the m × n sites are drawn. In this new convention no spurious sites need to be introduced. It is clear that

π_{h}^{I}(m, n) < π_{h}^{II}(m + 2, n + 2)

since π_{h}^{I} is the crossing probability on a (m + 2) × (n + 2) lattice with rules II where the horizontal bonds on the top and bottom lines and the vertical ones on the left and right columns are blocked. Because of these blocked bonds, the top and bottom lines of the (m + 2) × (n + 2) lattice cannot be used for horizontal crossings and:

π_{h}^{I}(m, n) < π_{h}^{II}(m + 2, n) < π_{h}^{II}(m + 2, n + 2).

This shows that, to first order in \frac{1}{n}

|π_{h}^{II}(\frac{m}{n}) − π_{h}^{I}(\frac{m}{n})| > \frac{2}{n} |(π_{h}^{II})′(\frac{m}{n})| .

Hence a simple change of conventions alter the result by this quantity. Table I gives an idea of the magnitude of this sensitivity of π_{h} for the two extremes of the measured range of r/r_{0} and for r/r_{0} = 1.

The indeterminacy due to the choice of conventions is, as is clear from Table I, substantial and inevitable. For the general events described in the introduction, and for other models than those considered here, it is even more serious, because the curves defining E can be strongly curved, or the symmetry of the model with respect to the coordinate axes severely violated, so that considerable thought has to be given to the necessary corrections. In the present paper it is sufficient to keep statistical errors within this indeterminacy. Because of computational limitations, this was not possible over the whole range of r/r_{0}. Instead we chose to measure each point with a sample size not smaller than 2.5 × 10^{5}. Statistical errors also appear on Table I. (Statistical errors are taken in this paper to represent a 95% confidence interval.) We observe that it is not π_{h} that appears in the graphs or that is analyzed in the next section but ln(\frac{π_{h}}{1−π_{h}}). Since the derivative of this function is (π_{h}(1 − π_{h}))^{−1}, any error in π_{h} is magnified by a factor of approximately 10^{3} at the ends of our range of investigation, so that possible statistical errors at the extremes are much larger than suggested by the last column of Table I. The statistical errors at the extremes are, however, effectively an order of magnitude smaller than those in the middle. As observed in Section 3, there is another improvement of one order of magnitude introduced by scaling, so that the figure of 10^{3} is ultimately reduced to 10. Even so, care has to be exercised with the results for the extreme points.

Table I. Sensitivity to Conventions

r  π_{h}  \frac{2}{n}|π′_{h}(\frac{m}{n})|  Statistical errors
7.3  7 × 10^{−4}  2 × 10^{−5}  10^{−4}
1  0.5  5 × 10^{−3}  2 × 10^{−3}
0.14  1.0  10^{−4}  10^{−4}

We used a random generator of linear congruential type, x_{i+1} = (ax_{i} + c) mod m, with: a = 142412240584757; c = 11; m = 2^{48}. It is of maximal period m.

3. NUMERICAL RESULTS

3.1 The determination of the critical probabilities

Critical probabilities have been studied extensively and in detail in the literature. For several reasons, we decided none the less to calculate again those we use. First of all, what appears to be a standard reference ([1]) considers only dimensions greater than two, and does not provide references to recent work in dimension two. Moreover, for obvious reasons, it gives the probabilities to only four places. So does, in some cases, [8, p. 17], and we preferred five places. Moreover the recent results that we could find ([6,9]) give, even when the intervals of error are taken into consideration, discrepant values. Finally for the investigation of universality of the crossing probabilities for models other than the standard ones, we will need simple, direct methods for calculating critical probabilities. It seemed useful to experiment on the standard models with the obvious ones, and to do so independently, applying clear, easily described principles.

If N = L^{2} is the lattice size then standard ideas (basically the existence and definition of the critical index ν [8, §4.1]) suggest that for an accuracy δ in the crossing probability we need an accuracy in the critical probability of AδL^{− \frac{1}{ν}}. The value of ν is \frac{4}{3}, and, at least, for a lattice of equal width and height, A can be taken to be of order 1. Here N is about 40000 so that if we take δ = .001 as Table I suggests, we need the critical probability to within two parts in 100, 000. This is what we have tried to achieve. For a horizontal crossing on rectangles with large or small aspect ratio r, the value of A = A_{r}^{h} could, however, be much different.

Finite-size scaling suggests that we introduce the function π_{h}(a, b) of the introduction as a function of L and r = \frac{a}{b}, and write it as in formula (55a) of [8] as

π_{h}(a, b) ∼ Φ_{r}((p − p_{c})L^{\frac{1}{ν}}),

where the function Φ_{r} depends on r. It is difficult to study the function Φ_{r} directly since we do not know its value at 0. We may however also write

π_{v}(a, b) ∼ Ψ_{r}((p − p_{c})L^{\frac{1}{ν}}),

and

π_{h}(a, b) + π_{v}(a, b) ∼ Θ_{r}((p − p_{c})L^{\frac{1}{ν}}).

The advantage of Θ_{r} is that we anticipate, as a result of the universality of crossing probabilities, that its value at 0 is 1. Thus, if we assume as in [8, §4.1, (55b)], that Θ_{r}(x) is a differentiable function of x, we have

Θ_{r}((p − p_{c})L^{\frac{1}{ν}}) = 1 + A_{r}(p − p_{c})L^{\frac{1}{ν}} + B,

where B is of the order of ((p − p_{c})L^{\frac{1}{ν}})^{2}, and thus negligible. The constant A_{r} is the sum of A_{r}^{h} and A_{r}^{v}.

Thus to estimate p_{c} and A_{r}, the latter for a given r, we find, by simulations and for a given L, the function Θ_{r}((p − p_{c})L^{\frac{1}{ν}}), treating it as a linear function of p, and then calculate its intercept with the line Θ = 1 and its slope. Whatever value we choose ultimately to take as an approximation for p_{c} , we can expect that the error it causes in the calculation of π_{h}(r) and π_{v}(r) is comparable to that in the difference between the values Θ_{r} and 1.

On the other hand, we were hesitant to anticipate in our calculations a result, the equality π_{h}(r) + π_{v}(r) = 1, that we were trying to establish. So we used a second method to calculate p_{c}. We started once again with the equation (55a) of [8], which asserts that near p_{c}, Φ_{r}((p − p_{c})L^{\frac{1}{ν}}) is a linear function of p. If one takes this seriously, it suffices to calculate the intersection of these two lines for two values of L in order to calculate both p_{c} and Φ_{r}(0) = π_{h}(r). Since the equation (55a) and its variants are by no means to be taken literally, we preferred to begin by making the procedures they entail explicit in a case for which p_{c} can be calculated exactly. This allows us also to verify that A_{r} is not too large for extreme values of r, on the contrary. We recall those values of the critical probabilities that are known exactly. For percolation by sites on a triangular lattice and for percolation by bonds on a square lattice p_{c} = \frac{1}{2}, and for percolation by bonds on a triangular lattice and on its dual, the hexagonal lattice, the critical probabilities are respectively 2 sin \frac{π}{18} and 1 − 2 sin \frac{π}{18}. Moreover the critical probabilities of the two remaining models have been established by several computational experiments (for site percolation on a square lattice see [6,9]).

To make the methods we use for the calculation of p_{c} clear, we consider percolation by sites on a triangular lattice. The choice of this model is easy to justify. Not only is p_{c} known to be \frac{1}{2}, but it is also known that, at p = p_{c}, the relation π_{h} + π_{v} = 1 holds for any value of s, even for finite lattices. We measured π_{h} and π_{v} for the following lattice sizes:

small lattices large lattices
186 × 215  558 × 644
83 × 480  249 × 1440

and for p from 0.49998 to 0.50002 in steps of size 0.00001. (Even though this is not relevant to the present discussion, the sizes were so chosen that the ratios 186/215 and 558/644 would make s as close as possible to 0 and the two other ratios 83/480 and 249/1440 so that s would be near to ln 5. See the next section.) The number of sites in the small lattices is roughly 4 × 10^{4} and the large lattices contain 9 times that number. For this experiment, the same set of random numbers (between 0 and 2^{48} − 1) were used for a given grid at the five different values of p. Hence we were sure from the beginning that the measurements of π_{h} and π_{v} would be increasing functions of p. The sample sizes have been 1005 K for the small lattices and 500 K for the large ones.

We begin with the results for the two lattices with s ∼ 0 thus r ∼ 1. The measured values of π_{h} + π_{v}, thus the function Θ_{r} are plotted on Figure 2.a as functions of p together with linear fits, the line with the largest slope belonging to the 558 × 644 lattice. (a)

(b)

Fig. 2. Numerical determination of p_{c} for the percolation by sites on a triangular lattice: (a) for the lattices with s ∼ 0, (b) for the lattices with s ∈ ln 5.

The first method yields p_{c} = 0.500003 (for the 186 × 215 lattice) and p_{c} = 0.500001 (for the 558 × 644 lattice). The second yields p_{c} = 0.499999. Calculating A_{1} by dividing the slopes of the lines in Figure 2.a by L^{\frac{1}{ν}}, we obtain in both cases approximately 1.5. Thus our initial estimates of the accuracy to be expected in the values of π_{h} and π_{v} were too generous, but by a factor of only 1.5.

What are the statistical errors for these numbers? A straightforward analysis using the linear regression hypothesis gives for the first method an error of 1.2 × 10^{−7} for the small lattice and of 2 × 10^{−7} for the large one. For the second method one gets an error of 3 × 10^{−7}. These are not to be taken seriously however as we failed to satisfy the independence of the measurements of the five points along the line. Because we used the same set of random numbers at the five values of p, we must face a possible shift in the intercept of the linear fit. This shift cannot be assumed to be less than the the accuracy of one of the five points. Since the sample size for the small lattice was 10^{6} and for the large one 5 × 10^{5}, the two methods give the following estimates for p_{c}:

first method  \{ \begin{matrix}p_{c} = 0.500003 ± 0.000018 & for the 186 × 315 lattice \\ p_{c} = 0.500001 ± 0.000011 & for the 558 × 644 lattice \end{matrix}
second method p_{c} = 0.499999 ± 0.000019.

Hence we can claim to have obtained the value of p_{c} to five places, the last one having an indeterminacy of ±2. The values coincide with the theoretical p_{c} = \frac{1}{2} in the limit of the error.

The values of π_{h}, π_{v} and π_{hv} at p = 0.50000 have been measured to be:

\begin{matrix}lattice & π_{h} & π_{v} & π_{hv} \\ 186 × 215 & 0.5010 & 0.4987 & 0.3220 \\ 558 × 644 & 0.4995 & 0.5002 & 0.3215 \end{matrix}

which coincide to four parts in 1000, as we expected. Hence the small lattice will be sufficient for the purpose at hand.

We turn now to the elongated lattices 83 × 430 and 249 × 1440. The estimates of p_{c} are:

first method  \{ \begin{matrix}p_{c} = 0.49998 ± 0.00006 & for the 83 × 430 lattice \\ p_{c} = 0.50002 ± 0.00004 & for the 249 × 1440 lattice \end{matrix}
second method p_{c} = 0.50005 ± 0.00008.

Moreover, calculating the slopes of the lines in Figure 2.b we obtain for A_{r}, r now being either about .2 or about 5, the values .084 and .083. This means that for the same accuracy in p_{c} we gain an additional figure at the ends of the interval. Given the increase in error when we pass to ln(\frac{π_{h}}{1−π_{h}}) that was mentioned in Section 2, this improvement is certainly welcome.

For the percolation by sites on square and hexagonal lattices, we chose p_{c} with the help of both methods described above. These methods agree fairly well to the first four signicant digits. For the square lattice, we used a 200 × 200 and a 600 × 600 grid. The sample sizes were 1.5 × 10^{6} and 8 × 10^{5} respectively. The results are:

first method  \{ \begin{matrix}p_{c} = 0.592712 ± 0.000014 & for the 200 × 200 lattice \\ p_{c} = 0.592740 ± 0.000009 & for the 600 × 600 lattice \end{matrix}
second method p_{c} = 0.592762 ± 0.000019.

The runs to calculate the functions π_{h}, π_{v}, π_{hv} and π_{d} were started before the final sample size for the present experiment was reached; p_{c} was set to 0.59272. The above results indicate that 0.59273 (or even 0.59274) might have been a better choice. The errors introduced by the early choice of p_{c}, if any, are smaller or equal to the statistical errors.

For the hexagonal lattice, we used a 265 × 153 and a 989 × 571 grid. The sample sizes were respectively 5 × 10^{6} and 10^{6} and the results are

first method  \{ \begin{matrix}p_{c} = 0.697014 ± 0.000007 & for the 265 × 153 lattice \\ p_{c} = 0.697034 ± 0.000006 & for the 989 × 571 lattice \end{matrix}
second method p_{c} = 0.697049 ± 0.000011.

We measured the functions π_{h}, π_{v}, π_{hv} and π_{d} at p_{c} = 0.69703.

3.2 The constants r_{0} and the four functions π_{h}, π_{v}, π_{hv} and π_{d}

Figure 3 shows the results for the percolation by bonds on a hexagonal lattice (dots) together with polynomial fits for the percolation by sites on a square lattice (curve); the functions plotted are ln \frac{π_{h}}{1−π_{h}}, ln \frac{π_{v}}{1−π_{v}}, ln π_{hv} and ln π_{d}. (See below for the discussion of the fits.) Since all models are visually indistinguishable, we present diagrams only for this comparison.

(a)

(b)

Fig. 3. The four crossing probability functions for the percolation by bonds on a hexagonal lattice (dots) and by sites on a square lattice (curve): (a) ln[π_{h}/(1 − π_{h})], (b) ln[π_{v}/(1 − π_{v})], (c) ln π_{hv},

(d) ln π_{d}.

(c)

(d)

Fig. 3. (Continued)

To plot the above curves for the hexagonal lattice, we had to fix the constant r_{0}. The constants r_{0} were also sought for the three other models not on the square lattice. Guided by the numerical evidence, we chose:

r_{0} = \sqrt{3} for the hexagonal lattice
r_{0} = \frac{\sqrt{3}}{2} for the triangular lattice

for the models of percolation by sites and by bonds. These values of r_{0} are those suggested by the hypothesis that r_{0} = 1 when the triangular and hexagonal lattices are represented in their usual symmetric forms, and the fundamental domains not deformed to rectangles as in our programs. Our numerical simulations strongly confirm these values.

As an example, we fitted a curve through the points of the function ln \frac{π_{h}}{1−π_{h}} of the model of percolation by sites on a hexagonal lattice. Since this function is thought to be odd, we used polynomials with terms (s − a), (s − a)^{3} and (s − a)^{5}, varying the parameter a to get the best fit. The best a was 0.0010 which corresponds to a difference between \sqrt{3} and the measured value of r_{0} of 0.1%. (For this fit we excluded the three points at both extremities of the range of s because of their low accuracy.) Similar results were obtained for the other functions and the other models. Because of the simplicity of their interpretation and the close agreement with those obtained by computation, we henceforth use the exact values and not the numerical estimates.

To compare the six models, we chose to measure one of them with a good accuracy and to fit the four curves ln \frac{π_{h}}{1−π_{h}}, ln \frac{π_{v}}{1−π_{v}}, ln π_{hv} and ln π_{d} with polynomials of the proper parity; this allowed us to compare points of different models with neighboring but distinct values of s. The easiest model to study was the model of percolation by sites on a square lattice since, by symmetry, only 41 of the 81 values of s in the range [−2, 2] need to be measured. For this model we pushed the sample size to over 10^{6}. For this size, the errors on the estimates of π_{h} and π_{v} vary from 10^{−3} for s ∼ 0 to 6 × 10^{−5} for s ∼ ±2. (Note that Δπ_{h}(s ∼ 2) ≃ 6 × 10^{−5} is a rather large relative error since π_{h}(s ∼ 2) ≃ 8 × 10^{−4}.) The results for this model are gathered in Appendix A. (Only the first 41 points are given, the others being obtainable by the permutation width ↔ height.) We have added, for convenience, two columns with the ratio r and its inverse.

To fit ln \frac{π_{h}}{1−π_{h}}, we tried odd polynomials with 3 and 4 terms. We tried also to fit the measurements, excluding the 3 points at both extremities of the range (s ∼ ±2). As the residual sum of squares [7] for 4 terms is almost equal to the the residual sum for 3 terms when the whole set of data is considered but larger when the extreme points are deleted, we conclude that our numerical study cannot give a proper estimate of the coefficient of the fourth term. Similar methods were used for ln π_{hv} and ln π_{d}. The results are:

fit of ln \frac{π_{h}}{1 − π_{h}} = p_{h}(s) = −2.062s − 0.305s^{3} − 0.022s^{5}
fit of ln \frac{π_{v}}{1 − π_{v}} = p_{v}(s) = p_{h}(−s)
fit of ln π_{hv} = p_{hv}(s) = −1.139 − 1.300s^{2} − 0.035s^{4} − 0.005s^{6}
fit of ln π_{d} = p_{d}(s) = −1.122 − 0.618s^{2} − 0.018s^{4} − 0.004s^{6}.

To compare each of the remaining models with the above one, we calculated the root-mean-square deviations from these fits; for example:

Δ_{h} = \sqrt{\frac{1}{n} ∑ [(ln \frac{π_{h}}{1 − π_{h}} − p_{h}) (s_{i})]^{2}},

the sum ranging over the n points measured (41 or 81). We also computed rms deviations Δ′, given by a similar expression but with the 3 points at both extremities of the range of s deleted.

(As the reader will see, these points are the main source of errors, because of their low accuracy.) The results are contained in Table II. The first line gives the rms deviations of the measurements used to obtain the fits and the fits themselves. The others are the deviations of the other models from the above fits for the percolation by sites on a square lattice. What are the acceptable rms deviations Δ for the sample size (∼ 2.5 × 10^{−5}) that we used for the other five models? The 7th line (labelled “statistics”) gives the rms deviations for the measured quantities assuming that the errors are of purely statistical origin.

TABLE II. Deviations from the Fits for the Percolation by Sites on a Square Lattice

Δ_{h}  Δ′_{h}  Δ_{v}  Δ′_{v}  Δ_{hv}  Δ′_{hv}  Δ_{d}  Δ′_{d}
square-sites  0.018 0.012 0.015 0.013 0.013 0.012 0.0047 0.0045
hexagonal-sites 0.023 0.012 0.025 0.022 0.032 0.015 0.0080 0.0074
triangular-sites 0.044 0.024 0.035 0.023 0.041 0.029 0.0186 0.0175
square-bonds 0.030 0.028 0.030 0.022 0.040 0.035 0.042 0.040
hexagonal-bonds 0.043 0.028 0.056 0.035 0.074 0.038 0.027  0.024
triangular-bonds 0.031 0.024 0.036 0.029 0.044 0.030 0.030  0.028
statistics  0.023 0.015 0.023 0.015 0.023 0.015 0.008 0.007
conventions  0.063 0.052 0.063 0.052  —  —  —  —

Note that the models of percolation by bonds are slightly more off than the ones by sites. This appears to be a consequence of the limitations arising from the convention used, and underlines the need for the careful choices of §2.1. Recall that the two conventions discussed in section 2.2 led to a systematic deviation of π_{h} given by

\frac{2}{n} |π′_{h} (\frac{m}{n})| .

This, by itself, produces a Δ_{h} = 0.063 (or Δ′_{h} = 0.052). Hence, we conclude that in the limit of our analysis and in the range of s studied, the four crossing probability functions coincide for the six models.

Systematic errors are most easily detected by examining 1−π_{h} −π_{v}. In Figure 4.c, for a model of bond percolation, a systematic positive error is clear. It is also clear in Figure 4.a but there it is negative, so that π_{h} + π_{v} tends to be greater than 1, even though we used a value for p_{c} that was slightly too small. These errors are presumably the result of the finite size of our lattices. For a triangular lattice, which is self-dual even at a finite size and for which we used a known value of p_{c}, Figure 4.b shows no systematic error. The left side of Figure 4.a is obtained from the right by reflection, and the one point on the right far above the horizontal axis appears to be a failing of our random-number generator In general, the error in π_{h} + π_{v} is seen to remain within the five parts in 1000, that has been our implicit goal, and the systematic errors, due to the conventions, dominate the statistical errors, especially in 4.c. The differences are smaller for large values of |s| because, for these points, one of the crossing probabilities is essentially 1 and absolute statistical errors are then minute on the scale chosen for the graph.

(a)

(b)

Fig. 4. The difference (1 − π_{h} − π_{v}) as a function of s: (a) percolation by sites on a square lattice,

(b) percolation by sites on a triangular lattice, (c) percolation by bonds on a triangular lattice.

As a last remark, we compare the results with a prediction of Cardy [2] about the asymptotic behavior of η_{h}(r) = π_{h}(r/r_{0}). Using finite-size scaling, he suggests that:

π_{h}(r) ∼ Ce^{−πr/3},

as r → ∞. In [2] he takes the constant to be 1, but that was an oversight and on the basis of more recent work ([3]) it appears that it should be 1.42635. Figure 5 displays the points

(r, −πr/3 + ln(C) − ln π_{h}(r/r_{0})

obtained for site-percolation on a square lattice. It shows clearly the limitations in the accuracy of our results for the verification of the prediction in this form. They do however permit the verification of a stronger prediction, a formula for the function η_{h} (cf. [3]).

Fig. 5. Asymptotic behavior of ln π_{h}.

APPENDIX

Table III contains the four crossing probability functions π_{h}, π_{v}, π_{hv} and π_{d} measured for the percolation by sites on a square lattice. The numbers were obtained with a sample size of over 10^{6}.

Table III. The Crossing Probability Functions π

width  height  r  r^{−1}  π_{h}  π_{v}  π_{hv}  π_{d}

200  200  1.000  1.0000  0.50072  0.50036  0.32250  0.32480
205  195  1.051  0.9512  0.47499  0.52645  0.32160  0.32637
210  190  1.105  0.9048  0.44926  0.55181  0.31878  0.32299
216  186  1.161  0.8611  0.42149  0.57657  0.31143  0.31892
221  181  1.221  0.8190  0.39796  0.60303  0.30597  0.31835
227  176  1.290  0.7753  0.37071  0.63016  0.29562  0.31297
232  172  1.349  0.7414  0.34925  0.65291  0.28640  0.30787
238  168  1.417  0.7059  0.32400  0.67676  0.27328  0.30142
244  164  1.488  0.6721  0.30177  0.69927  0.26060  0.29459
250  160  1.562  0.6400  0.27850  0.72155  0.24601  0.28656
257  156  1.647  0.6070  0.25556  0.74512  0.23031  0.27934
263  152  1.730  0.5779  0.23474  0.76561  0.21537  0.27009
270  148  1.824  0.5481  0.21223  0.78836  0.19800  0.25864
277  145  1.910  0.5235  0.19495  0.80668  0.18381  0.25170
284  141  2.014  0.4965  0.17496  0.82583  0.16668  0.23970
291  137  2.124  0.4708  0.15609  0.84461  0.15019  0.22912
298  134  2.224  0.4497  0.14094  0.86028  0.13679  0.21713
306  131  2.336  0.4281  0.12482  0.87531  0.12180  0.20637
314  128  2.453  0.4076  0.11120  0.89023  0.10919  0.19447
322  124  2.597  0.3851  0.09624  0.90563  0.09491  0.18224
330  121  2.727  0.3667  0.08367  0.91724  0.08276  0.17080
338  118  2.864  0.3491  0.07208  0.92765  0.07146  0.15863
347  115  3.017  0.3314  0.06177  0.93820  0.06141  0.14842
355  113  3.142  0.3183  0.05484  0.94623  0.05456  0.13935
364  110  3.309  0.3022  0.04591  0.95498  0.04576  0.12732
374  107  3.495  0.2861  0.03811  0.96278  0.03804  0.11628
383  104  3.683  0.2715  0.03091  0.96962  0.03087  0.10561
393  102  3.853  0.2595  0.02612  0.97431  0.02608  0.09662
403  99  4.071  0.2457  0.02118  0.97927  0.02117  0.08662
413  97  4.258  0.2349  0.01731  0.98299  0.01730  0.07906
423  94  4.500  0.2222  0.01357  0.98674  0.01356  0.07010
434  92  4.717  0.2120  0.01075  0.98945  0.01074  0.06190
445  90  4.944  0.2022  0.00864  0.99180  0.00864  0.05559
456  88  5.182  0.1930  0.00654  0.99356  0.00654  0.04862
468  85  5.506  0.1816  0.00468  0.99550  0.00468  0.04173
480  83  5.783  0.1729  0.00351  0.99650  0.00351  0.03598
492  81  6.074  0.1646  0.00273  0.99745  0.00273  0.03082
504  79  6.380  0.1567  0.00190  0.99813  0.00190  0.02687
517  77  6.714  0.1489  0.00142  0.99868  0.00142  0.02245
530  75  7.067  0.1415  0.00098  0.99908  0.00098  0.01845
544  74  7.351  0.1360  0.00077  0.99930  0.00077  0.01587

References

[1] J. Adler, Y. Meir, A. Aharony, and A. B. Harris, Series study of percolation moments in general dimension, Phys. Rev. B, 41, 9183-9206 (1990)

[2] J.L. Cardy, Finite-size scaling in strips: Antiperiodic boundary conditions, J. Phys. A, 17, L961-L964 (1984)

[3] J.L. Cardy, Critical percolation in finite geometries, Preprint (1991)

[4] G. Grimmett, Percolation, Springer (1989)

[5] H. Kesten, Percolation Theory for Mathematicians, Birkhaüser (1982)

[6] M. Rosso, J.F. Gouyet, B. Sapoval, Determination of percolation probabilities from the use of a concentration gradient Phys. Rev. B, 32, 6053-6054 (1985).

[7] G.A.F. Seber, Linear Regression Analysis, Wiley (1977)

[8] D. Stauffer, Introduction to Percolation Theory, Taylor Francis (1985)

[9] R.M. Ziff,Test of scaling exponents for percolation-cluster experiments, Phys. Rev. Lett., 56, 545-548 (1986)