Percolation properties of the classic Sierpinski carpet and sponge

TL;DR

This paper demonstrates that the iterative construction of Sierpinski carpets and sponges exhibits critical phenomena similar to percolation, deriving critical exponents and using finite-size scaling to determine transport properties.

Contribution

It introduces a novel analogy between fractal constructions and percolation, deriving critical exponents and applying finite-size scaling for transport properties.

Findings

Critical exponents are derived for Sierpinski structures.

Finite-size scaling accurately predicts scalar transport properties.

The percolation analogy applies at any iteration stage.

Abstract

Iterative construction of a Sierpinski carpet or sponge is shown to be a critical phenomenon analogous to uncorrelated percolation. Critical exponents are derived or calculated (by random walks over the carpet or sponge at infinite iteration) that are related by equations identical to those obtained from percolation theory. Finite-size scaling then gives accurate values for the scalar transport properties (e.g., effective conductivity) of the carpet or sponge at any stage of iteration.

Tables

$d$	$d_w$	$D_0$	$t/\nu$
$2$	$2.09675$	$0.77376$	$0.203957$
$3$	$2.02026$	$0.935312$	$0.0546151$

steps/walk	${\delta }_{\text{2D}}$	${\delta }_{\text{3D}}$
${10}^2$	$2.11382\left(175\right)$	$2.02364\left(107\right)$
${10}^3$	$2.10737\left(103\right)$	$2.02219\left(89\right)$
${10}^4$	$2.10526\left(101\right)$	$2.02150\left(60\right)$
${10}^5$	$2.10357\left(82\right)$	$2.02127\left(41\right)$

Equations

$$\mathcal{H}=\frac{\ln \left(b^2-m\right)}{\ln b}=\frac{\ln 8}{\ln 3}\approx 1.89279$$

$$\sigma =⟨\sigma (\mathbf{r})⟩D_w$$

$$D_w=\frac{⟨R(t{)}^2⟩}{2dt}$$

$$D_w=D_0{\left(\frac{\xi }{{\xi }_0}\right)}^{2-d_w}$$

$${\xi }^{(i)}=a_i^{-\nu }$$

$$\nu =\frac{\ln 3}{\ln \left(9/8\right)}={\left(2-\mathcal{H}\right)}^{-1}.$$

$${\sigma }^{(i)}={\sigma }_1{a}_i{D}_w^{(i)}={\sigma }_1{a}_i{D}_0{\left(\frac{{\xi }^{(i)}}{{\xi }_0}\right)}^{2-d_w}$$

$${\sigma }^{(i)}={\sigma }_1{a}_i{D}_0{\left({\xi }^{(i)}\right)}^{2-d_w}={\sigma }_1{D}_0{a}_i^t$$

$$t=1+\nu \left(d_w-2\right).$$

$$\sigma (\xi )\sim {\sigma }_1{D}_0{\xi }^{-t/\nu }$$

$$\sigma (L)={\sigma }_1{D}_0{L}^{-t/\nu }$$

$$\mathcal{H}=\frac{\ln \left(b^3-m\right)}{\ln b}=\frac{\ln 26}{\ln 3}\approx 2.96565$$

$$\nu =\frac{\ln 3}{\ln \left(27/26\right)}={\left(3-\mathcal{H}\right)}^{-1}.$$

$$D_w(\xi )\sim {\xi }^{2-d_w}$$

$$a(\xi )\sim {\xi }^{-1/\nu }.$$

$$d_w<2+\frac{1}{\nu }=2+(d-\mathcal{H}).$$

$$D_w(L)=D_0{L}^{2-d_w}.$$

$$\frac{⟨R(t{)}^2⟩}{2dt}=D_0{⟨R(t{)}^2⟩}^{1-d_w/2}$$

$$⟨R(t{)}^2⟩={\left(2dt{D}_0\right)}^{2/d_w}.$$

$$\ln ⟨R(t{)}^2⟩=\frac{2}{d_w}\ln t+\frac{2}{d_w}\ln \left(2d{D}_0\right).$$

$$\frac{\xi }{{\xi }_0}=n^{1/d_w}.$$

$$\frac{R}{{\xi }_0}=n^{1/\delta }$$

$$\delta (n)=\frac{2\ln n}{\ln ⟨R(n{)}^2⟩}$$

$$\varphi =\frac{{\left(\frac{\xi }{{\xi }_0}\right)}^D}{{\left(\frac{\xi }{{\xi }_0}\right)}^3}={\left(\frac{\xi }{{\xi }_0}\right)}^{D-3}$$

$$\sigma (L)={\sigma }_e{D}_0{\left(\frac{L}{{\xi }_0}\right)}^{-t/\nu }={\sigma }_e\varphi D_0{\left(\frac{L}{{\xi }_0}\right)}^{2-d_w}.$$

$$P_{ij}=\frac{{\sigma }_j}{{\sigma }_i+{\sigma }_j}{\left[\sum _{k=1}^{2d}\left(\frac{{\sigma }_k}{{\sigma }_i+{\sigma }_k}\right)\right]}^{-1}.$$

$$T_i={\left[2\sum _{k=1}^{2d}\left(\frac{{\sigma }_k}{{\sigma }_i+{\sigma }_k}\right)\right]}^{-1}.$$

Full text

I Introduction

The classic Sierpinski carpet (Mandelbrot, ) is a recursive, self-similar fractal embedded in two-dimensional (2D) Euclidean space. The generator is shown in Fig. 1 (the center square is removed, leaving the surrounding eight squares); the carpet after the third iteration is shown in Fig. 2. This center-hole $(3,1)$ Sierpinski carpet has Hausdorff (fractal) dimension

[TABLE]

given the scaling factor $b=3$ and number $m=1$ of eliminated squares in the generator.

With each iteration one-ninth of the mass is removed while the carpet remains self-similar. The correlation length $\xi$ for this system is the length of the carpet, and so increases to infinity   when measured by comparison to the size* $\mathrm{h}_{i}$* of the smallest hole   as the number i of iterations goes to infinity. Thus the iterative construction of this recursive, self-similar fractal resembles the approach to a percolation threshold.

In this paper, relations between critical exponents are derived that are identical to those characterizing uncorrelated percolation. This suggests the term “recursion percolation” for the critical phenomena exhibited by recursive fractals.

The following section briefly describes the Walker Diffusion Method (WDM) by which the analytical and numerical results for the 2D Sierpinski carpet and 3D Sierpinski sponge are obtained. The percolation analogy is then developed in Sec. III. Values for parameters specific to the transport properties of the carpet and sponge cannot be derived: those are calculated in Sec. IV by random walks over a carpet or sponge of infinite size (infinite iteration). Relevant papers in the literature are discussed in Sec. V. Concluding remarks are made in Sec. VI.

II Walker Diffusion Method

This application of the WDM (CVS99, ) utilizes the relation

[TABLE]

between the effective conductivity $\sigma$ of a composite material and the (dimensionless) diffusion coefficient $D_{w}$ obtained from walkers diffusing through a digital representation of the composite. The factor $\left\langle\sigma(\mathbf{r})\right\rangle$ is the volume average of the constituent conductivities.

In effect, the phase domains that make up the composite are host to walker populations, where the walker density of a population is proportional to the conductivity value of its host domain. The principle of detailed balance ensures that the population densities are maintained, by providing the following rule for walker diffusion over the digitized composite: a walker at site (or pixel/voxel) i attempts a move to a randomly chosen adjacent site j during the time interval $\tau=\left(4d\right)^{-1}$, where d is the Euclidean dimension of the system; this move is successful with probability $p_{ij}=\sigma_{j}/\left(\sigma_{i}+\sigma_{j}\right)$, where $\sigma_{i}$ and $\sigma_{j}$ are the conductivities of sites i and j, respectively. The path of the walker thus reflects the composition and morphology of the domains that are encountered.

The diffusion coefficient $D_{w}$ is calculated using the equation

[TABLE]

where the set $\left\{R\right\}$ of walker displacements, each occurring over the time interval $t$, must have a Gaussian probability distribution that is necessarily centered well beyond $\xi$. The correlation length $\xi$ is identified as the length scale above which a composite material attains the “effective”, or macroscopic, value of a scalar transport property (electrical conductivity, for example).

For displacements $R<\xi$, the walker diffusion is anomalous rather than Gaussian due to the heterogeneity of the composite at length scales less than $\xi$. However there is an additional characteristic length $\xi_{0}<\xi$, below which the composite is effectively homogeneous. Then a walker displacement of $\xi$ requiring a travel time $t_{\xi}=\xi^{2}/\left(2\,d\,D_{w}\right)$ is produced by a walk comprising $\left(\xi/\xi_{0}\right)^{d_{w}}$ segments of length $\xi_{0}$, each requiring a travel time of $t_{0}=\xi_{0}^{2}/\left(2\,d\,D_{0}\right)$, where $D_{0}$ is the walker diffusion coefficient calculated from displacements $R\leq\xi_{0}$. Thus $t_{\xi}=\left(\xi/\xi_{0}\right)^{d_{w}}t_{0}$, which gives the relation

[TABLE]

between the walker diffusion coefficient $D_{w},$ the fractal dimension $d_{w}$ of the walker path, and the correlation length $\xi$.

III Recursion percolation

For convenience the Sierpinski carpet or sponge at iteration i is denoted by $\mathrm{S}_{i}$. Further, all lengths, areas, and volumes are in units of the smallest feature, which is the smallest hole $\mathrm{h}_{i}$. That is necessarily the size of a single pixel or voxel, so $\mathrm{h}_{i}$ has length/area/volume equal to 1.

As these fractals are self-similar at all length scales, the characteristic length $\xi_{0}=h_{i}=1$.

With each iteration i, the areal fraction $a_{i}$ of the $\mathrm{S}_{i}$ carpet that is conducting (i.e., not lost to cut-outs) decreases according to $a_{i}=\left(8/9\right)^{i}$, while the correlation length $\xi^{(i)}$ increases according to $\xi^{(i)}=3^{i}\,\xi_{0}=3^{i}$. These two equalities can be written $\ln a_{i}=i\ln\left(8/9\right)$ and $\ln\xi^{(i)}=i\ln 3$, respectively, so producing the power-law relation

[TABLE]

with the exponent

[TABLE]

This result makes the case that recursive construction of the Sierpinski carpet is a percolation-like phenomenon (S=000026A, ). Of course, the “percolation threshold” is approached (from above) as iteration $i\rightarrow\infty$, causing $a_{i}\rightarrow 0$.

According to Eqs. (2) and (4), the effective conductivity of $\mathrm{S}_{i}$ is

[TABLE]

where $a_{i}$ is the areal fraction of the carpet that is conducting, and $\sigma_{1}$ is the conductivity of that material. The walker diffusion coefficient $D_{0}<1$ because walkers near the non-conducting cut-outs tend to linger there (an effect of the walker diffusion rule stated in Sec. II). Both $D_{0}$ and the walker path dimension $d_{w}$, which pertain to the transport properties of the carpet, must be calculated, not derived.

Use of Eq. (5) in Eq. (7) gives

[TABLE]

with the exponent

[TABLE]

Then use of Eq. (5) in Eq. (8) produces the asymptotic relation

[TABLE]

giving the finite-size scaling relation

[TABLE]

for all $L=3^{i}$. [Note that $\sigma(L)$ is the effective conductivity of an infinite 2D array of carpets of size $L$. The length $L$ plays the role of the correlation length $\xi$.]

The formalism above is straightforwardly applied to the Sierpinski sponge as well. All relations for the sponge are obtained by replacing $a_{i}$ with $v_{i}$ in the equations above. Here $v_{i}=\left(26/27\right)^{i}$ is the volume fraction of the $\mathrm{S}_{i}$ sponge that is conducting. The Hausdorff dimension of the sponge is

[TABLE]

giving the exponent

[TABLE]

As the exponent $d_{w}$ is central to the development above, it is interesting to consider whether that development imposes any limits on its value. In fact an analytic bound is obtained by comparing the asymptotic behavior (meaning: as iteration $i\rightarrow\infty$) of $D_{w}^{(i)}$ with that of the conducting areal fraction $a_{i}$ of the carpet or conducting volume fraction $v_{i}$ of the sponge. From Eq. (4),

[TABLE]

and from Eq. (5),

[TABLE]

The value $D_{w}^{(i)}$, reflecting the walker behavior, is responsive to the value $a_{i}$ or $v_{i}$ (rather than vice versa), suggesting that $1/\nu>d_{w}-2$. Thus

[TABLE]

Then $d_{w}<2.10721$ in the case of the 2D carpet, and $d_{w}<2.03435$ in the case of the 3D sponge.

Note that Eq. (16), giving an analytic upper bound for $d_{w}$, should apply to similar recursive fractals as well.

The relations in this section are recognizable from standard percolation theory (CVSperc, ). In particular, the latter gives $\nu/\beta=\left(d-D\right)^{-1}$ and $t=\beta+\nu\left(d_{w}^{*}-2\right)$, where $D$ is the fractal dimension of the incipient infinite cluster of conducting sites, and the exponent $\beta$ is less than 1, reflecting the fact that conductor sites not belonging to the percolating cluster cannot contribute to the conductivity of the system. In contrast, the Sierpinski carpet and sponge do not have such “stranded” conductor sites (thus the “percolation thresholds” $a_{\infty}=0$ and $v_{\infty}=0$).

A more important distinction is that in standard percolation the exponent relations and values are obtained in the limit $\xi\rightarrow\infty$; that is, at the percolation threshold. Because $d_{w}$ increases from $2$ to $d_{w}^{*}$ as the correlation length $\xi\rightarrow\infty$, the derivation that leads to Eq. (16) is not applicable: that derivation relies on the constancy of the $d_{w}$ value over all iterations of the carpet and sponge.

IV Numerical methods and results

The value of the walker path dimension $d_{w}$, and the value of the diffusion coefficient $D_{0}$ associated with the length scale $\xi_{0}$, must be calculated. For a fractal system of finite size $L$, Eq. (4) may be written

[TABLE]

This relation can be expressed in terms of the computable variable $\left\langle R(t)^{2}\right\rangle$:

[TABLE]

which simplifies to

[TABLE]

Thus walks over the fractal system will produce points $(\ln t,\ln\left\langle R(t)^{2}\right\rangle)$ that satisfy the equation

[TABLE]

A linear fit to the points produces a plot from which the values $d_{w}$ and $D_{0}$ can be ascertained.

This graphical approach is taken for the Sierpinski carpet and sponge. These may be created using the subroutine given in Appendix A. Note that in both cases the smallest hole, being indivisible, is the size of a single site (pixel or voxel).

Walks over the carpet or sponge are accomplished by use of the variable residence time algorithm (CVS99, ), described in Appendix B. The algorithm takes advantage of the statistical nature of the diffusion process to eliminate (while accounting for) unsuccessful attempts by the walker to move to a neighboring site.

To allow very long walks, all walks are actually taken over a carpet or sponge at infinite iteration. Note that the subroutine locates conducting and non-conducting sites with respect to an origin, which in the case of the sponge is the $\left(i,j,k\right)$ site with coordinates $\left(0,0,0\right)$. Thus the sponge occupies all space with site index values $i,j,k$ greater than or equal to [math]. A move by a walker at $\left(i,j,k\right)$ is of course determined by the conductivities of the adjacent sites: those values (1 or 0) are obtained by calls to the subroutine.

Figure 3 is the plot of points obtained from walks over the infinite Sierpinski carpet. The slope $2/d_{w}$ of the fitted line gives the value $d_{w}=2.09675$. The y-intercept $(2/d_{w})\ln(2\,d\,D_{0})$ of the line gives the value $D_{0}=0.77376$.

Each point in Fig. 3 is obtained from $40$ sequences of $10^{6}$ walks of time $t$. [A sequence of $10^{6}$ walks is actually a single walk of time $10^{6}\times t$. During that long walk every displacement $R(t)$ is recorded, for a total of $10^{6}$ displacements.] The plotted value $\left\langle R(t)^{2}\right\rangle$ is the average of all walks of time $t$ (that is, the average of all sequences). In every case the number of sequences is sufficient that additional sequences would change the value $\left\langle R(t)^{2}\right\rangle$ by only an insignificant amount (far less than the point size in the figure).

A sequence of walks is initiated by placing a walker at a randomly chosen conducting site $\left(i\gg 0,j\gg 0\right)$ of the infinite carpet.

Similar calculations are made for the 3D Sierpinski sponge.

Figure 4 is the plot of points obtained from walks over the infinite Sierpinski sponge. The slope $2/d_{w}$ of the fitted line gives the value $d_{w}=2.02026$. The y-intercept $(2/d_{w})\ln(2\,d\,D_{0})$ of the line gives the value $D_{0}=0.935312$.

Each point in Fig. 4 is obtained from $40$ sequences of $10^{6}$ walks of time $t$.

For convenience the numerical results obtained above are presented together in Table I.

Interestingly, an upper bound on the value $d_{w}$ can be obtained by considering the number of steps in a walk. Recall from Sec. II that a walker displacement $\xi$ is produced by $(\xi/\xi_{0})^{d_{w}}$ steps, each of length $\xi_{0}$. Thus the relation between the number of steps $n=(\xi/\xi_{0})^{d_{w}}$ and displacement $\xi$ is

[TABLE]

This relation applies to displacements $R<\xi$ as well. However, calculations for short walks may be affected by the finite size of the sites that compose the carpet and sponge. Therefore consider the relation

[TABLE]

where both the displacement $R$ and the exponent $\delta$ are determined by the number of steps $n$. In any case, the calculated value $\delta\rightarrow d_{w}$ as $R\rightarrow\xi$.

These Sierpinski fractals have the characteristic length $\xi_{0}=1$. Then it is computationally convenient to square both sides of Eq. (22), since $R^{2}$ has an integer value. With these changes, Eq. (22) produces the relation

[TABLE]

where $\left\langle R(n)^{2}\right\rangle$ is the average of all the $R^{2}$ values obtained from a very large number of walks of $n$ steps.

Thus the exponent $\delta(n)\rightarrow d_{w}$ as the number $n$ of steps in a walk increases. This is apparent in Table II, where values $\delta(n)>d_{w}$ are recorded. Each value is obtained from $40$ sequences of walks, each sequence comprising $10^{5}$ walks, over the infinite carpet or sponge. The average value $\left\langle R(n)^{2}\right\rangle$ used in Eq. (23) is taken from all walks of $n$ steps.

Note that a value given as $1.234(5)$ means $1.234$ with standard deviation $0.005$, and so indicates the range $1.229$ to $1.239$, centered on $1.234$. The standard deviation for $\delta(n)$ is calculated from the $40$ values obtained by the $40$ sequences of walks.

V Prior research

A number of papers report calculations (and in one case an experiment) that directly or indirectly produce a value for $d_{w}$. Note that a value given as $1.234(5)$ indicates the range $1.229$ to $1.239$, centered on $1.234$.

Gefen et al. (Gefen, ) overlay a resistor network on the classic Sierpinski carpet, and obtain a finite-size scaling relation for the resistance by a renormalization method. That produces the resistance exponent $\widetilde{\zeta}_{R}=0.194$, giving $d_{w}=\mathcal{H}+\widetilde{\zeta}_{R}=2.087$. (Note, by the way, that $t/\nu=d-2+\widetilde{\zeta}_{R}$ (S=000026A, ).)

Barlow et al. (Barlow, ) consider resistor network approximations of the classic Sierpinski carpet, and apply electrical circuit theory to calculate the spectral dimension $d_{s}=1.80525$, giving $d_{w}=2\,\mathcal{H}/d_{s}=2.09698$.

Kim et al. (Kim, ) perform random walks on the $\mathrm{S}_{7}$ carpet, and obtain $d_{w}=2.106(16)$ from their relation $\left\langle R(N)^{2}\right\rangle\propto N^{2/d_{w}}$, where $N$ (with values up to $10^{4}$) is the number of steps in a walk.

Zhuang et al. (Zhuang, ) measure the resistivity for carpets $\mathrm{S}_{1}$, $\mathrm{S}_{2}$, and $\mathrm{S}_{3}$ cut out of copper and aluminum sheets. Those data points lie on a line $\ln\left(1/R\right)\propto\ln L$, where $R(L)$ is the resistance and $L$ is the size scaling of the samples, indicating finite-size scaling with conductivity exponent $t/\nu=0.22(1)$.

Aar$\tilde{\mathrm{a}}$o Reis (AaraoReis, ) performs random walks on carpets $\mathrm{S}_{2}$, $\mathrm{S}_{3}$, $\mathrm{S}_{4}$, and $\mathrm{S}_{5}$, and overlays the four plots of $\left\langle R(N)^{2}\right\rangle^{1/2}/L$ versus $L\,N^{-\nu_{w}}$ ($L=3^{i}$, and $N$ up to $8\times 10^{4}$ steps). The best data collapse (coincidence of the four curves) occurs for $\nu_{w}=0.476(5)$, giving $d_{w}=1/\nu_{w}=2.101(22)$. This approach is taken as well for sponges $\mathrm{S}_{2}$, $\mathrm{S}_{3}$, and $\mathrm{S}_{4}$, producing $\nu_{w}=0.492(6)$, giving $d_{w}=1/\nu_{w}=2.033(25)$.

Suwannasen et al. (Suwannasen, ) perform random walks with $2^{15}$ walkers initially distributed at random on the $\mathrm{S}_{15}$ carpet. Their log-log plot of $\left\langle R(t)^{2}\right\rangle$ versus $C\,t^{2/d_{w}}$ gives $d_{w}=2.10(1)$.

VI Concluding remarks

A notable aspect of this paper is the appearance of the factor $D_{0}$. Most importantly, it appears in the finite-size scaling relation Eq. (11) giving the effective conductivity of a finite carpet or sponge.

Random walk calculations that relate the number of steps in a walk to the consequent walker displacement $R$ are shown to give an upper bound for the value of the exponent $d_{w}$. The longer the walk (that is, as $R\rightarrow\xi$), the tighter the upper bound.

In addition, an analytic expression giving an upper bound for $d_{w}$ is derived by comparing the asymptotic ($i\rightarrow\infty$) behavior of $D_{w}^{(i)}$, and $a_{i}$ or $v_{i}$, governed by the exponents $d_{w}$ and $\nu$, respectively. Thus Eq. (16) should apply to similar recursive fractals as well.

The Sierpinski sponge has been used as a heuristic model for porous rock, which is typically found to have a fractal character over some range of length scales. In fact the relations derived here carry over to isotropic porous rock, with the exception that the conducting volume fraction $v_{i}$ of the $\mathrm{S}_{i}$ sponge is replaced by the connected porosity $\phi$ of the rock, given by the relation

[TABLE]

where the right-hand side of the equality is the volume fraction of the rock that is occupied by the connected pore space of (fractal) mass dimension $D<3$. The corresponding exponent is $\nu=(3-D)^{-1}$.

A porous rock sample of size $L$ saturated with an electrolyte solution having conductivity $\sigma_{e}$ would exhibit an effective conductivity given by the finite-size scaling relation Eq. (11):

[TABLE]

Acknowledgements.

I thank Professor Indrajit Charit (Department of Nuclear Engineering & Industrial Management) for arranging my access to the resources of the University of Idaho Library (Moscow, Idaho).

Appendix A: Sierpinski fractal construction subroutine

This subroutine determines whether element $i,j$ (carpet) or $i,j,k$ (sponge) of the array representing the Sierpinski fractal is conducting or insulating, and returns the value $1$ or [math], respectively. Note that a corner of the array is an element with $i=j=k=0$. This implementation is written in C.

int SierpinskiFractal(int i, int j, int k)

{

while (i>0 || j>0 || k>0) {

if (i%3 == 1 && j%3 == 1 && k%3 == 1) return 0;

i /= 3;

j /= 3;

k /= 3;

}

return 1;

}

Appendix B: Variable residence time algorithm

According to this algorithm (CVS99, ), the actual behavior of the walker is well approximated by a sequence of moves in which the direction of the move from a site $i$ is determined randomly by the set of probabilities $\left\{P_{ij}\right\}$, where $P_{ij}$ is the probability that the move is to adjacent site $j$ (which has conductivity $\sigma_{j}$) and is given by the equation

[TABLE]

The sum is over all sites adjacent to site $i$. The time interval over which the move occurs is

[TABLE]

Note that this version of the variable residence time algorithm is intended for orthogonal systems (meaning a site in a 3D system has six neighbors, for example).