The density of critical percolation clusters touching the boundaries of strips and squares

TL;DR

This paper combines conformal field theory and high-precision simulations to analyze the density of critical percolation clusters touching boundaries in various 2D geometries, providing new theoretical and numerical insights.

Contribution

It derives theoretical formulas for cluster densities touching boundaries in 2D percolation and validates them with high-precision numerical simulations.

Findings

Density of crossing clusters in an infinite strip is proportional to a specific sine and cosine function.

Contours for cluster densities in squares and rectangles are numerically determined and compared with theory.

Theoretical predictions match well with high-precision numerical results.

Abstract

We consider the density of two-dimensional critical percolation clusters, constrained to touch one or both boundaries, in infinite strips, half-infinite strips, and squares, as well as several related quantities for the infinite strip. Our theoretical results follow from conformal field theory, and are compared with high-precision numerical simulation. For example, we show that the density of clusters touching both boundaries of an infinite strip of unit width (i.e. crossing clusters) is proportional to $(\sin \pi y)^{-5/48}\{[\cos(\pi y/2)]^{1/3} +[\sin (\pi y/2)]^{1/3}-1\}$. We also determine numerically contours for the density of clusters crossing squares and long rectangles with open boundaries on the sides, and compare with theory for the density along an edge.

Figures

Equations

$${\rho }_b(y)=(\sin \pi y{)}^{-5/48}\left[{\left(\cos \frac{\pi y}{2}\right)}^{1/3}+{\left(\sin \frac{\pi y}{2}\right)}^{1/3}-1\right].$$

$${\rho }_e(y)=(\sin \pi y{)}^{-5/48}.$$

$${\rho }_0(y)$$

$${\rho }_1(y)$$

$$\rho(z,x_{a},x_{b})\propto(z-\bar{z})^{-5/48}F\bigg{(}\frac{(x_{b}-x_{a})(\bar{z}-z)}{(\bar{z}-x_{a})(x_{b}-z)}\bigg{)}\ ,$$

$$F_{\pm}(\eta)=\bigg{(}\frac{2-\eta}{2\sqrt{1-\eta}}\pm 1\bigg{)}^{1/6}\ .$$

$${\rho }_{+}(r,\theta ,x_a\to 0,x_b\to \infty )$$

$${\rho }_{-}(r,\theta ,x_a\to 0,x_b\to \infty )$$

$${\rho }_a(r,\theta )\propto (r\sin \theta {)}^{-5/48}.$$

$$w(z)=\frac{1}{\pi }\log (z).$$

$${\rho }_{0\bar{1}}(y)$$

$${\rho }_{1\bar{0}}(y)$$

$${\rho }_b(y)={\rho }_0(y)+{\rho }_1(y)-{\rho }_e(y),$$

$$y=(Y+a)/(512+2a).$$

$$\rho (x)={\left(\frac{|z^{\prime }(x)|}{1-\lambda z(x)}\right)}^{1/3}.$$

$$z(w)=1-\wp \left(\left(\frac{1}{2}+\frac{i}{2}\right)(iw+1)K(2);4,0\right),$$

$${\rho }_b(x,y=0)=(\sin \pi x{)}^{1/3}.$$

$${\rho }_b(x,y\gg 0)=(\sin \pi x{)}^{11/48}.$$

$${\rho }_b(x=0,y)=(\sin \pi y{)}^{-1/3}\left[{\left(\cos \frac{\pi y}{2}\right)}^{2/3}+{\left(\sin \frac{\pi y}{2}\right)}^{2/3}-1\right],$$

$${\rho }_0(y{)}^6+{\rho }_1(y{)}^6={\rho }_e(y{)}^6.$$

$${\int }_0^1{\rho }_b(y)dy=1.$$

$$B$$

$$=$$

$${\rho }_0\sim y^{-5/48},$$

$${\rho }_b\sim y^{11/48}.$$

Full text

percolation, cluster density, crossing

I Introduction

Percolation in two-dimensional systems is an area that has a long history, but remains under very active current study. A number of very different methods has been applied to critical 2-D percolation, including conformal field theory (CFT) [1], other field-theoretic methods [2], modular forms [3], computer simulation [4], Stochastic Loẅner Evolution (SLE) processes [5] and other rigorous methods [6]. (Because the literature is so very extensive we have cited only a few representative works.)

More specifically, there is a great deal of work of recent work studying universal properties of crossing problems in critical percolation in two dimensions (i.e., [7, 1, 8, 9, 10, 5, 11, 3, 12]). Another interesting and also practically important universal feature of percolation at the critical point is the density, defined as the number of times a site belongs to clusters satisfying some specified boundary condition (such as clusters touching certain parts of the boundary) divided by the total number of samples $N$, in the limit that $N$ goes to infinity. This problem has been addressed for clusters touching any part of the boundary of a system in various geometries, including rectangles, strips, and disks, via conformal field theory [13] and by solving the problem for a Gaussian free field and then transforming to other statistical mechanical models, including percolation [14]. In a recent Letter [30], we considered the problem of clusters simultaneously touching one or two intervals on the boundary of a system, and considered cases where those intervals shrink to points (anchor points). These results exhibit interesting factorization, are related to two-dimensional electrostatics, and highlight the universality of percolation densities. Note that the density at a point $z$ of clusters which touch specified parts of the boundary is proportional to the probability of finding a cluster that connects those parts of the boundary with a small region around the point $z$.

In this paper we consider the problem of the density $\rho_{b}$ of critical percolation clusters in various geometries where the clusters simultaneously touch both of the boundaries (i.e. crossing clusters), and several related quantities.

The first case we consider is an infinite strip, with boundaries parallel to the $x$-axis at $y=0$ and $y=1$, so the crossing is in the vertical direction. (All our models are defined so that the crossing is vertical. Fig. 4 below illustrates the geometries that we consider.) For the infinite strip we find (leaving out an arbitrary normalization constant here and elsewhere) that

[TABLE]

This may be compared to a previous result [13, 14] for clusters touching either the upper or lower boundary (or both) which is simply given by

[TABLE]

We also show that the density of clusters touching one boundary irrespective of touching the other is given by

[TABLE]

where $\rho_{0}$ corresponds to those clusters touching the lower boundary at $y=0$, and $\rho_{1}$ corresponds to those clusters touching the upper boundary at $y=1$. (Note that $\rho_{0}$ is the analog of the order parameter profile $\langle\sigma\rangle_{+,f}$ for the Ising case (see Eq. (16) in [19].) We also find expressions for clusters touching one boundary and not the other, which are combinations of the above results.

Perhaps the main new theoretical result in the above is (3) (or equivalently (4)), which follows straightforwardly from the results in [30]. The derivation is given in section II.

A second type of theoretical prediction gives the density variation along a boundary (the general expressions is in (15)). This is used to predict the density along the edge in several geometries (see (17), (18), and (19) below).

The above theoretical results are found to be consistent with numerical simulations to a high degree of accuracy. We include the results of numerical simulations of the density contours of clusters crossing square and rectangular systems vertically, with open boundaries on the sides, and compare with theory along the boundaries.

Our theoretical treatment is related to previous use of conformal field theory to study order parameter profiles in various 2-D critical models with edges and similar research [13, 15, 16, 17, 18, 19, 20, 21, 22, 23, 24, 25, 26, 27, 28]. (Note that the density which we consider is the expectation value of the spin operator in percolation, which is the order parameter in this setting.) Many of these prior results make use of the original Fisher-de Gennes proposal [29] for the behavior of the order parameter at criticality near a straight edge. In this paper, we limit ourselves to critical percolation. We also include the results of high-precision computer simulations. In addition, the formula (15) for the density along the edge of a system is new, to our knowledge.

Note that, because of the fractal nature of critical percolation clusters, the density of clusters is, strictly speaking, zero everywhere in the system. However, if we properly renormalize the density as the lattice mesh size goes to zero, the density can remain finite. At the boundaries, for some quantities, this results in the density diverging, as for $\rho_{e}(y)$ at $y=0$ and $y=L$ (but remaining integrable). For $\rho_{b}(y)$ of equation (1), on the other hand, the renormalized density remains finite everywhere. When comparing to numerical simulations, one has to normalize the data so that the densities coincide with the theoretical prediction using whatever normalization convention is chosen for the theoretical results. The resulting normalization constant, which must be applied to the numerical data, is specific for each system and is non-universal.

In the following sections, we first give the theoretical derivation of our infinite strip formulas above. Then we present the numerical results. This is followed by numerical results on square and (long) rectangular systems. These are compared with theory for the density along the edges of these systems. We end with a few concluding remarks.

II Theory for the infinite strip

We first consider the density of critical percolation clusters which span the sides of an infinite 2-D strip. We can find the density predicted by conformal field theory [31] using the results of [30]. In that article we showed that in the upper half plane the density $\rho$ of clusters connected to an interval $(x_{a},x_{b})$ is

[TABLE]

where the function $F(\eta)$ was determined by conformal field theory and takes on one of two forms,

[TABLE]

If we parameterize $z$ as $z=re^{i\theta}$ and let $x_{a}\to 0$ and $x_{b}\to\infty$, then $\eta=1-e^{2i\theta}$ and using (6) we can rewrite (5) as

[TABLE]

For the positive real axis $\theta\to 0$ and $\rho_{+}\sim\theta^{-5/48}$ while for the negative real axis $\theta\to\pi$ and $\rho_{+}\sim(\pi-\theta)^{11/48}$. The powers here arise from the fixed and free boundary exponents, respectively, in the bulk-boundary operator product expansion of the magnetization operator $\psi$ [30]. (More specifically, as it approaches the boundary, $\psi\sim\mathbf{1}$ or $\psi\sim\phi_{1,3}$, which have conformal dimensions [math] and $1/3$, respectively.) This shows that $\rho_{+}$ is the density of clusters attached to the positive axis. Because $\rho_{-}(r,\theta)=\rho_{+}(r,\pi-\theta)$, it follows that $\rho_{-}$ is the density of clusters attached to the negative real axis.

The final density that we need is that of clusters connected arbitrarily to the axis. This is given by $\langle\psi(z,\bar{z})\rangle_{\rm fixed}\propto(z-\bar{z})^{-5/48}$ [30, 13, 14] which may also be written

[TABLE]

These densities are unnormalized. However for points that are short distances above the positive (negative) axis, but very far from the origin, the relation $\rho_{+(-)}\approx\rho_{a}$ holds. This condition holds since the points are far from the free boundary, and thus dominated by the fixed boundary. It is satisfied by our expressions for $\rho_{+},\rho_{-},$ and $\rho_{a}$, so they are properly normalized relative to one another.

We next map these densities into the infinite strip w\in$$\{x+iy| $x\in(-\infty,\infty)$, $y\in(0,1)\}$ using

[TABLE]

This leads to the expressions for $\rho_{0}(y)$, $\rho_{1}(y)$, and $\rho_{e}(y)$ given by equations (3), (4), and (2), respectively.

Using these functions we can also determine the densities of clusters that touch one side but not the other,

[TABLE]

In a similar manner we can find the density of clusters touching both sides, $\rho_{b}(y)$. Adding $\rho_{0}$ and $\rho_{1}$ includes all configurations that touch either side, but double counts the clusters that touch both sides. Subtracting $\rho_{e}$ leaves only those clusters that touch both sides of the strip. Thus

[TABLE]

which gives equation (1).

III Simulations for the infinite strip

To approximate the infinite strip, we considered rectangular systems with periodic boundary conditions in the horizontal direction, for both site and bond percolation on square lattices. Here we report our results for site percolation on the square lattice, for a system of 511 (vertical) $\times$ 2048 (horizontal) sites, at $p=p_{c}=0.5927462$ [32]. We generated 500,000 samples to compute the average densities, using a Leath-type of algorithm to find all clusters touching the upper and lower boundaries.

The aspect ratio of the rectangle used in our simulations was $2048/511=4.008\ldots$, which might seem a bit small. However, since the correlation length of the system is governed by the width of the rectangle ($511$), the effect of the finite ratio drops off exponentially with the length of the rectangle, so our results should be very close to those of an infinite strip, an expectation which is born out by the results described below. (Furthermore, the probability of finding a horizontally wrapping cluster drops exponentially with the aspect ratio, so if a longer rectangle were used, very few wrapping clusters would be found and the data would have poor statistics.)

The agreement of predicted and simulated values is excellent. In Fig. 1 we show the results for $\rho_{b}(y)$ with no adjustments, other than an overall normalization (in particular, the extrapolation length $a$, described below, is set to zero). Fig. 2 (upper set of points) shows the ratio of simulation to theoretical results, equation (13). Most points fall within 1%, except those right near the boundary. However, we can do better.

When comparing simulations on a (necessarily) finite lattice with the results of a continuum theory, there is the question, due to finite-size and edge effects, of what value of the continuous variable $y$ should correspond to the lattice variable $Y$ — specifically, where the boundary of the system should be placed. On the lattice, the density goes to zero at $Y=0$ and $Y=512$. A naïve assignment of the continuum coordinate would therefore be $y=Y/512$. However, we can get a better fit to the data near the boundaries by assuming that the effective boundary is a distance $a$ (in units of the lattice spacing) beyond the lattice boundaries — that is, at $Y=-a$ and $Y=512+a$. The distance $a$ is an effective “extrapolation length” where the continuum density far from the wall extrapolates to the value zero [33]. This is accomplished by defining the continuum variable $y$ by

[TABLE]

Then, $Y=-a$ corresponds to $y=0$ and $Y=512+a$ corresponds to $y=1$. Note, $Y=0$ corresponds to $y=a/(512+2a)$ and $Y=512$ corresponds to $y=(512+a)/(512+2a)=1-a/(512+2a)$, and so the theoretical extrapolated density $\rho_{b}(y)$ will be greater than zero at these points on the actual boundary. The spacing between all points is stretched by a small amount because of the denominator in equation (14), but this stretching does not have much effect on the behavior of $\rho_{b}$ near the center. The main effect of $a$ on the shape of the theoretical curves of $\rho_{b}$ is near the boundaries.

By choosing an extrapolation length of $a=0.26$, we can get a much better fit of the data, as can be seen in Fig. 2 (lower set of points). A plot of the data analogous to Fig. 1 puts most of the data points right in the center of the theoretical circles, but would barely be visible when plotted on this scale. With $a=0.26$, the error is now reduced to less than 0.1%, except right near the boundaries, as can be seen in Fig. 2. A more thorough study of the extrapolation coefficient $a$ would require the study of different sized lattices, and the demonstration that $a$ is independent of the lattice size. We have not carried this out. Note, however, that a constant $a$ implies that if one keeps the physical size of the lattice fixed (so that the increasing number of lattice points makes the mesh size go to zero), the extrapolation length, measured in physical units, will also go to zero.

Note also that the distance in the $y$ directions was 511 rather than 512 because we used one row at $Y=0$, in conjunction with the periodic boundary conditions, to represent both horizontal boundaries of the system on the lattice. That is, $Y=1$ and $Y=511$ are the lowest and highest rows where we occupy sites in the system, where $Y$ represents the lattice coordinate in the vertical direction.

We note that these simulations were carried out before the theoretical predictions were made.

We have also carried out simulations measuring the density of clusters touching one edge, $\rho_{0}(y)$. The results are shown in Fig. 3. We also plot the results of the theoretical prediction, equation (3), on the same plot, and find agreement within 0.5% without using an extrapolation length $a$.

IV Percolation on a square and semi-infinite strip

In this section, we consider the density of crossing clusters on a square with open boundaries and also on a (long) rectangle. We compare the numerical results with the predictions of conformal field theory for the density along an edge. The various cases considered are illustrated in Fig. 4.

Note that, as mentioned, the crossing is always in the vertical direction. In a slight abuse of notation, we use $\rho_{b}$ for the density of a clusters that touch both anchoring intervals in all cases. The different situations may be distinguished by the arguments of $\rho_{b}$, e.g. $\rho_{b}(y)$ for the infinite strip, where there is no $x$ dependence, and $\rho_{b}(x,y=0)$ along the bottom of the semi-infinite strip as in (17) below.

Our simulations of percolation densities on an open square examined clusters that cross in the vertical direction, with open boundary conditions on the left- and right-hand sides. We considered site percolation on a square lattice of $511\times 511$, with 2,000,000 samples generated. The resulting contours are shown in Fig. 5. As in the periodic case, the density goes to zero at the upper and lower boundaries because, compared to an infinite system, these boundaries intersect many possible crossing paths, leading to large holes in the clusters near the boundaries. As a consequence, relatively few points on the boundary will be part of the crossing clusters, and in the limit that the mesh size goes to zero, their density evidently goes to zero. The density also goes to zero on the sides because of the open conditions there. Interestingly, the contour curves are almost symmetric in the horizontal and vertical directions, indicating that the anchoring and open boundaries have a similar effect on the density.

We have not carried out a field-theory calculation to find the density inside the square or rectangular systems. To do so requires a six-point function whose calculation would be unwieldy.

It is however relatively easy to calculate the variation of the density along the bottom edge of the square, now normalized so that the density remains non-zero as the mesh size goes to zero. To do this, we consider the density of crossings from one of the anchoring intervals to a single point on the other interval, using the boundary spin operator. Now crossing from the top edge to one point $x$ on the bottom edge (which is given by a three-point function, depending only on $x$ and $\lambda$) automatically implies crossing from the top to bottom (which is given by a four-point function, depending only on $\lambda$). Therefore the density at $x$ is proportional to the ratio of the former to the latter. It follows generally that, up to a ($\lambda$-dependent) multiplicative constant, one has [34]

[TABLE]

Here $z(w)$ maps the region of interest ($w$) onto the 1/2-plane ($z$), $x$ is the $w$-coordinate along the anchoring interval of interest, and $\lambda$ is the conformally invariant cross-ratio for the anchoring points. For a rectangle, this depends on the aspect ratio $r=K(\sqrt{1-\lambda})/K(\sqrt{\lambda})$ [35, 1].

The mapping for the square is

[TABLE]

with $\wp(u;g_{2},g_{3})$ the Weierstrass elliptic function and $K(z)$ the elliptic integral function. This mapping takes a unit square $x,y\in(0,1)$ into the upper half plane. For the square $\lambda=1/2$, and we can take $x\in(0,1)$. In Fig. 6 we compare the measurement and theory. Clearly the agreement is excellent.

In the case of a half-infinite strip $x\in(0,1),y\in(0,\infty)$, the density of sites along the $x$-axis belonging to clusters crossing vertically is found from (15) using $z(w)=\sin^{2}(\pi w/2)$, and $\lambda=0$. This gives

[TABLE]

Of course, for an infinite strip, the probability of crossing (in the long direction) is zero, so one must consider the limit of a large system, and calculate the density given that crossing takes place, and take the limit that the length of the rectangle goes to infinity. It turns out that numerically, the density at the edge for the square, equations (15,16), differs only very slightly from that of the half-infinite strip, given by equation (17). From the point of view of the density along an anchoring edge, the square is not much different from the half-infinite strip.

For $y\gg 0$ in the above half-infinite strip (or equivalently for any $y$ for a fully infinite strip in the vertical direction) one can also find the density along the $x$-direction of the vertically crossing clusters in closed form

[TABLE]

This function may be found by transforming the density of clusters connecting two boundary points, derived in [30], into the infinite strip. We then take the limit as the two anchoring points move infinitely far away in opposite directions, while normalizing the density so that it remains finite. (Related order-parameter profiles for the Ising case are studied in [22].) Interestingly, a plot of this density profile (written as a function of $y$ rather than $x$) is very similar in appearance to that of the vertically crossing clusters $\rho_{b}(y)$ given in equation (1) and plotted in Fig. 1. When normalized so that they have the same value at $y=1/2$, the maximum difference between the two is at $y=0$ and $y=1$, where (18) is only 1.5 % below (1). This small difference indicates that open boundaries and anchoring boundaries have similar effects on the density of the crossing percolation clusters, and is consistent with the near symmetry seen in the contours in Fig. 5.

We can also find the density along the left and right (open) edges. For the case of a half-infinite strip rotated by $90^{\circ}$ with respect to the one above, ($y\in(0,1),x\in(0,\infty)$), we find

[TABLE]

which is similar in form to $\rho_{b}(y)$ of equation (1) (which in fact corresponds to $x\gg 0$ for the geometry considered here) but with different exponents. This similarity arises because the derivations of (19) and (13) are virtually identical, except that for (19) we leave a free interval between the two anchoring intervals (where the boundary spin operator sits) which is mapped to the end of the half-infinite strip using sine functions. Comparison with numerical data (not shown) for the density at the short edge of a rectangle of aspect ratio 8 to approximate the infinite strip shows excellent agreement with (19).

V Further comments and conclusions

If we consider densities raised to the sixth power (compare Ref. [30]), we find a Pythagorean-like relation involving $\rho_{e}(y)^{3}$, $\rho_{0}(y)^{3}$, and $\rho_{1}(y)^{3}$ (which we present without interpretation),

[TABLE]

There seems to be no simple relation involving $\rho_{b}(y)$ other than its basic definition given in equation (13). For the corresponding quantities at the edge of the strip as in (19), the power in equation (20) is 3 rather than 6.

Although the overall normalization of a density such as $\rho_{b}(y)$ (for the infinite strip, see (1)) is arbitrary, we can fix its value by requiring that

[TABLE]

That is, we define $\rho_{b}(y)=B\sin(\pi y)^{-5/48}(\cos(\pi y/2)^{1/3}+\sin(\pi y/2)^{1/3}-1)$ where $B$ is a constant. Then equation (21) yields

[TABLE]

Another choice of $B$ is to make $\rho_{b}(1/2)=1$, which yields $B=(2^{5/6}-1)^{-1}=1.27910371\ldots$

In many problems of percolation density profiles, the density goes to infinity at a boundary point, such as occurs for $\rho_{0}(y)$ when $y\to 0$. Interestingly, in the case considered here of clusters touching both boundaries, $\rho_{b}(y)$, the density goes to zero at those boundaries and remains finite everywhere.

To highlight the difference between densities of all clusters touching a boundary vs. the densities of crossing clusters touching one boundary, we consider the limit that the strip width becomes infinite, so that the system becomes a half-plane. Because we have written our results for a strip of fixed (unit) width, this density is given by the behavior of $\rho$ for small $y$. The density of all clusters touching the $x$-axis is thus found from equation (3) (or [13, 14]) to be

[TABLE]

where we have left off the coefficient because the normalization is arbitrary. This density diverges at $y=0$ because it is much more likely to find sites belonging to these clusters near the $x$-axis. On the other hand, the density of crossing clusters touching the $x$-axis is found from equation (1) in the limit $y\to 0$,

[TABLE]

In this case, the density increases as $y$ increases, in contrast to (24), and goes to zero at the anchoring boundary, for the reason mentioned above.

Behavior identical to (25) also follows from the small-$x$ expansion of $\rho_{b}(x,y\gg 0)$ given by (18), which represents the behavior of the density of the vertically spanning clusters at the open boundary. Thus, near the boundaries (but away from the corners), the open and anchoring boundary conditions have identical effects on the density of the crossing clusters.

In conclusion, we have studied the density of vertically percolating clusters in a square system, as well as for half-infinite and infinite strips extending in either the horizontal or vertical direction. The various cases considered are illustrated in Fig. 4.

For the infinite strip, Fig. 4(a), the density for crossing clusters is given by (19) and compared with numerical simulations in Figs. 1 and 2. For the square, Fig. 4(b), we find theoretical results for the anchoring edge densities (see (15) and (16) and Fig. 6), as well as numerical results for the density in the interior (Fig. 5), which exhibit an interesting near symmetry. For the half-infinite strip in the horizontal direction, Fig. 4(e), the density $\rho_{b}(x=0,y)$ at the left open boundary is given by (19), and at $x\gg 0$ (or equivalently, for an infinite strip), the density is given by (1). For a half-infinite strip in the vertical direction, the density along the lower anchoring boundary (Fig. 4(c)) is given by (17) while for $y\gg 0$ (Fig. 4(d)) the density is given by (18). For the half-infinite systems, the densities near a wall are given by the same power-law (25) regardless of whether it anchors the crossing clusters or is open. Note that all of our theoretical predictions were confirmed by computer simulation.

For the future, it would be interesting to study analogous properties for Fortuin-Kasteleyn (FK) clusters of the critical Ising and Potts models.

VI Acknowledgments

This work was supported in part by the National Science Foundation under grants numbers DMS-0553487 (RMZ) and DMR-0536927 (PK).