The percolation critical polynomial as a graph invariant

0.13 MB·English

by  Christian R. Scullard

The percolation critical polynomial as a graph invariant
∗
Christian R. Scullard
Lawrence Livermore National Laboratory, Livermore CA 94550, USA
(Dated: November3, 2011)
Everylattice forwhich thebondpercolation critical probability canbefoundexactly possesses a
critical polynomial, with the root in [0,1] providing the threshold. Recent work has demonstrated
thatcriticalpolynomialscanbedefinedonanyperiodiclattice. Ingeneral,thepolynomialdepends
on the lattice and on its decomposition into identical finite subgraphs, but once these are speci-
1 fied, the polynomial is essentially unique. On lattices for which the exact percolation threshold is
1 unknown, the polynomials provide approximations for the critical probability with the estimates
0 appearingtoconvergetotheexactanswerwithincreasingsubgraphsize. Inthispaper,Ishowhow
2 the critical polynomial can be viewed as a graph invariant like the chromatic and Tutte polyno-
mials. Like these, the critical polynomial is computed on a finite graph and may be found using
v
the deletion-contraction algorithm. This allows efficient calculation on a computer, and I present
o
such results for the kagome lattice using subgraphs of up to 36 bonds. For one of these, I find the
N
prediction pc = 0.52440572..., which differs from the numerical value, pc = 0.52440503(5), by only
4 6.9×10−7.

]
n Percolation is the study of the formation of random
n clusters on lattices. Given an infinite lattice, we declare
-
s eachbondtobeopenwithprobabilityp,andclosedwith
i probability 1−p. The resulting random clusters grow
d
. in average size with p until we reach the critical point,
t
a pc, at which an infinite cluster first appears. The de-
m termination of pc is an unsolved problem in all but one
- dimension and on two-dimensional lattices that are self-
d
dual 3-uniform hypergraphs[1]. An example is shownin
n
Figure 1a, where the shaded triangle can represent any
o
c configuration of bonds or sites. Critical thresholds on
[ these lattices are given by the Ziff criterion [2],

1 P(A,B,C)=P(A¯,B¯,C¯) (1) FIG. 1. A 3-uniform hypergraph. The shaded region need
v not be a simple triangle.
1 where P(A,B,C) is the probability that all three ver-
6 tices are connected through open paths, and P(A¯,B¯,C¯)
0 the martini-A lattice, with the probability assignments
is the probability that none are connected. Application
1
of (1) to find a bond threshold results in a polynomial shown in Figure 1a. If we set p1 = 1, this bond is con-
.
1 in p, with order at most equal to the number of bonds tracted and its end vertices merged, the result being the
1 honeycomb lattice with two bonds doubled in parallel.
in the unit triangle. We may consider each bond of an
1
1 n−bondtriangletohaveadifferentprobability,resulting This doubled bond can be replaced by a single effective
bond, so we get,
: in a critical surface,
v
Xi f(p1,...,pn)=0 (2) H(1−[1−p2][1−p3],p4,p5). (5)
r where the function f is at most first order in any of its
a By setting p1 = 0, we delete the bond, and the result is
arguments. Irefertothispropertyas“linearity”. Exam-
thesquarelatticewithbondsdoubledinseries. Thus,we
ples of such critical surfaces include the square lattice,
have,
S(p1,p2)=1−p1−p2, (3)
S(p2p4,p3p5). (6)
and the the honeycomb lattice,
H(p1,p2,p3)=1+p1p2p3−p1p2−p1p3−p2p3. (4) The A lattice is a 3-uniform hypergraph and its critical
surface given by (1). The only way to reduce to the
Note that critical surfaces may be multiplied by an ar- correctdeleted and contracted surfaces and preserve the
bitrary constant, or even another function with no root
required linearity property is to set
in [0,1], without changing the predicted bond thresh-
old. Removing the former freedom by requiring the con- A(p1,p2,p3,p4,p5)=
stanttermtobe+1,andthelatterbydemandinglinear- p1H(1−[1−p2][1−p3],p4,p5)
ity, makes the critical polynomial unique. Consider now +(1−p1)S(p2p4,p3p5), (7)
2
whichhastheformalappearanceofanaverageofthetwo 3. Conversely, if the critical polynomial gives the ex-
specialcases. Thisformisgeneralforlatticesthatcanbe act threshold for a base of a single unit cell, as is
solvedusing(1),andconstitutesthedeletion-contraction the caseforself-dual3−uniformhypergraphs,then
formula for critical surfaces of such systems. As another any critical polynomial using a larger base makes
example, we find for the martini-B lattice (Figure 1b), the same prediction, i.e. the original polynomial
alwaysfactorsout. Thisisalsoaconjecturedprop-
B(p1,p2,p3,p4)=p4S(p1,1−[1−p3][1−p4])
erty but it would seem to be necessary for consis-
+(1−p4)S(p1p3,p4). (8)
tency, and I have found no counter-examples.
Thedeletion-contractionalgorithmcanbeusedtoextend
4. In many cases, if the single-cell polynomial does
the definition of the critical polynomial to lattices for
not givethe exactanswer,then it canbe shown[8]
which the critical threshold is not known exactly. For
that no critical polynomial derived from a finite-
the kagome lattice of Figure 2a, we may write,
sized base will give the correct percolation thresh-
K(p1,p2,p3,p4,p5,p6)= old. Thisistrueforthekagomelattice,asdiscussed
p4B(1−[1−p5][1−p6],p1,p2,p3) below.
+(1−p4)A(p1,p2,p3,p5,p6), (9)
Critical polynomials have been found for all the
where A and B are given by (7) and (8). Critical poly-
Archimedean lattices [8, 10, 11]. The polynomial (10)
nomials are found by setting all probabilities equal, and
was originally found by Wu using his “homogeneity” [3]
in this case K(p,p,p,p,p,p)=0, i.e.,
assumption and recently, he extended the method to in-
1−3p2−6p3+12p4−6p5+p6 =0 (10) clude kagome subnets [9], with excellent results [12]. In
fact, the homogeneity approximation gives a prediction
which gives the well-known [3, 4] approximation pc = for the full q−state Potts criticalfrontier, but at present
0.52442971..., compared to the numerical [5] pc = itappearslimitedtothekagomelatticeanditssubnetsas
0.52440502(5). Thus, a critical polynomial may be de- it relies on a transformation from the triangular lattice
finedonanylatticeforwhichknownlatticesappearupon to these kagome-type graphs. Although I focus on the
deletion andcontractionofa bond, but there is no guar- kagomelatticeinthispaper,mostlybecauseoftheinter-
antee that the resulting predictionwill be exact. Infact, est it has attracted over the years (e.g. [3, 13, 14]), the
the definition extends to any lattice, for if either of the method presented here may be applied on any periodic
latticesappearinginthedeletion-contractionformulahas lattice for either site or bond percolation.
anunknownthreshold,onesimplyappliestheformulare- In previous work, critical polynomials were found “by
cursivelyandcontinues this procedureuntil knowncases hand”, with the (4,6,12) lattice and its 18−bond unit
appear. This method is used to compute other graph cellprobablyrepresentingthelimitofwhatonewouldbe
invariants such as the chromatic and Tutte polynomials inclined to do this way (see the appendix of [8]). How-
[6, 7]. ever, the recursive nature of the algorithm makes it an
In the following, I call the subgraph of a lattice, on idealproblemfora computer,andinthis paper,I report
which different probabilities are assigned, the “base” of the results of using a program to calculate critical poly-
theprocess,whichisthentiledtoformtheinfinitelattice. nomials on the kagome lattice for bases containing 12,
For example, in Figures 2a, 2b, 3, 4 are bases for the 24, and 36 bonds. At each step, the program chooses a
kagome lattice consisting of one, two, four, and six unit bond, and finds the two graphsthat result from its dele-
cells. Some known and conjectured properties of critical tionandcontraction. It knowsasmallnumber ofgraphs
polynomials are as follows: (e.g., triangular, honeycomb, and square lattices), and
it repeats the deletion-contraction algorithm recursively
1. They are unique (up to a trivial multiplication by
untilitrecognizesallthelatticesithasfound. Theoutput
a constant) once the base and tiling are specified
isasetoffunctiondefinitions,likeequations(7),(8)and
[8, 9]. In particular, they are independent of the
(9)whichcanbeevaluatedinacomputeralgebrapackage
bonds chosen to delete and contract at each step.
to get the critical surface and then the critical polyno-
In this sense, the critical polynomial is a property
mial. Of course, there are many issues to overcome in
of a finite graph.
programming this scheme, and a more full account will
2. If the single-cell prediction is not exact, it is gen- be given elsewhere. Here, I show some examples of its
erally very close, usually within 10−5 of the nu- operation on the kagome lattice.
merically determined threshold [10]. In addition,
proper choices of larger bases lead to predictions
closer to the exact answer, with accuracy increas- Base of 2 unit cells
ing with base size. This is not proven, but the
argument for it is the wide variety of examples in The first extension beyond a single unit cell base is
which it is clearly seen to be true. thatshowninFigure2, inwhichwe employa base using
3
FIG. 3. A 24−bond base for the kagome lattice with two
inequivalent embeddings.
FIG. 2. a) a base for thekagome lattice using two unit cells, −31068p22+4248p23−259p24 =0, (13)
shown in two different shades. Each bond on the base has a
different probability; b) a tiling for a); another tiling that is and the solution on [0,1] is pc = 0.52440672..., which
equivalent to b). differs from the numerical result by 1.7×10−6, a great
improvement over the single-cell 6−bond case. There
does not appear to be a way to get a better prediction
two unit cells which are indicated by different colors in using a 24−bond base.
Figure2a. Thepolynomialcanbewritteninthefactored
form
Bases of 6 unit cells
−(1−3p2−6p3+12p4−6p5+p6)
×(−1−p2−2p3+10p4−10p5+3p6). (11)
There are many options for bases of 36 bonds, and I
willdiscussonlythreehere. Ifwetakethe6-unitcellver-
Werecognizethefirstterminbracketsasthepolynomial
forthe 6−bondbase. Thesecondtermhasnorealroots, sion of the base shown in Figure 3a with the embedding
analogousto Figure 3b, we once againfind a polynomial
and thus the prediction here is the same as the one we
found previously, pc = 0.52442971... . There is no other in which (10) appears as a factor. Thus, for this base
and embedding, we getthe 6−bondestimate again. It is
way to tile this base to give a different result.
probably safe to assume that this trend continues with
larger bases and embeddings of this type.
Bases of 4 unit cells Turning now to the base shown in Figure 4a, in which
the embedding is indicated by the matching shapes on
Usingthebaseconsistingof4unitcellsshowninFigure the external vertices, we get the polynomial,
3a, we may tile it in two different ways as indicated in 1−3p4−12p5−20p6−60p7−132p8+56p9
Figures 3b and 3c. Starting with Figure 3b, the critical 10 11 12 13
+684p +1440p +2108p +2052p
polynomial can be written in the factored form −10452p14−68708p15−82980p16+280152p17
−(1−3p2−6p3+12p4−6p5+p6)× +1316026p18−49980p19−12878976p20
(−1−p2−2p3+10p4−10p5+3p6)× +5124684p21+90816816p22−199458252p23
(1−2p2−4p3+7p4+24p5−28p6−64p7+ −12979085p24+816398808p25−1939348056p26
172p8−184p9+110p10−36p11+5p12) (12) [phone]p27−2575935942p28
[phone]p29−984362272p30
Onceagain,the firstterminbracketsis justthe 6−bond [phone]p31−121897767p32+26954680p33
polynomial (10), and the others have no roots in [0,1]. −4096134p34+382956p35−16617p36=0, (14)
Thepredictionisagainthesameasthe6−bondestimate.
Things finally become more interesting when we use which has solution on [0,1] pc = 0.52440607..., differing
the staggered embedding in Figure 3c. In this case, the fromthe numericalvalue by 1.1×10−6. A different base
polynomial is andembeddingisshowninFigure4bandhaspolynomial,
1−6p4−24p5−24p6−24p7+27p8 1−6p4−24p5−14p6+36p7+39p8−100p9
+552p9+1056p10−1224p11−8548p12 −462p10+780p11+4583p12+4812p13
−4872p13+68568p14−50664p15 +−391227363p164p1−7+711600901p11546−p1885−62560p916340p19
−226650p16+643944p17−843684p18 −9675936p20+5297340p21+66607704p22
+684384p19−368886p20+133152p21 −151097304p23−5319734p24+610494828p25
4
have a significant effect on the rate of reduction to the
known simple cases. I have hardly explored this issue
and presently use what amounts to a random bond se-
FIG. 4. 36−bond bases for the kagome lattices. Matching
shapesontheexternalverticesindicatehoweachbaseistiled
to create thelattice. FIG.5. a)thesubgraphandembeddingresultingfromdelet-
ing 22 bonds in Figure 4a, which is a single cell of kagome
lattice and some connecting bonds; b) contracting the con-
−1461237180p26+2022998000p27 necting bonds gives the kagome lattice partitioned into unit
−1949295060p28+1387593528p29 cells.
−745850356p30+303533928p31
−92388675p32+20427736p33−3103578p34
+290052p35−12579p36=0, (15) lection, rejecting only those choices that lead to undue
complications. Moreover, this algorithm is perfect for a
with solution in [0,1], pc = 0.52440572...,slightly better parallel implementation as it would require little inter-
than the previous estimate and within 6.9×10−7 of the
processor communication. It remains to be determined
numerical value.
how much extra performancecan be wroughtfromthese
Clearly,weareobservingthepredictionsconvergingto
considerations,buttheproblemwillhopefullybe seenas
the exact answer with increasing base size. From these
an interesting computational challenge.
few examples it also appears that the estimates are con-
I thank Robert Ziff for helpful advice on this
verging from above, as they all seem to be greater than
manuscript. This work was performed under the aus-
the numerical value. However, I know of no argument
picesoftheU.S.DepartmentofEnergybyLawrenceLiv-
that guarantees this trend will continue. Note also that
ermore National Laboratory under Contract DE-AC52-
no polynomial for the kagome lattice derived from any
07NA27344.
finite-size base will give the exact answer. This can be
seen as follows. Using the base in Figure 4a, we delete
many of the bonds to leave a single unit cell and a few
connecting bonds, as in Figure 5a. Now, we may con-
tracttheconnectingbondstorecoverthekagomesystem ∗ [email protected]
in which the base is a single cell, as in Figure 5b. By [1] B. Bollob´as and O. Riordan,arXiv:1001.4674 (2010).
uniqueness,the predictionforthiscasemustbeequation [2] R. M. Ziff, Phys. Rev.E 73, 016134 (2006).
(10), which contradicts the prediction found by setting [3] F. Y.Wu,J. Phys.C 12 (1979).
all bonds equal in the full critical surface, i.e. equation [4] C.R.ScullardandR.M.Ziff,Phys.Rev.E73,045102(R)
(14). Thus,ifthesingle-cellpolynomialdoesnotprovide (2006).
[5] X. Feng, Y. Deng, and H. W. J. Bl¨ote,
the exact threshold, and for the kagome lattice it does
Phys. Rev.E 78, 031136 (2008).
not, then, although we can get arbitrarily close by using
[6] B. Bollob´as, Modern Graph Theory (Springer,1998).
ever larger bases, no finite critical polynomial derived in
[7] G. Haggard, D. J. Pearce, and G. Royle,
this way will ever solve the problem. ACM Trans. Math. Softw. 37, 24:1 (2010).
I have presented critical polynomials for the kagome [8] C. R.Scullard, J. Stat. Mech. , P09022 (2011).
lattice, up to bases of 36 bonds. This is the limit of [9] F. Y.Wu,Phys. Rev.E 81, 061110 (2010).
the current implementation of the algorithm, as it be- [10] C. R. Scullard and R. M. Ziff, J. Stat. Mech. , P03021
(2010).
comes progressively more difficult to add bonds due to
[11] C. R.Scullard and R. M. Ziff, Phys. Rev.Lett. , 185701
the exponential complexity. It is not uncommon for a
(2008).
large calculation to produce over a million function def- [12] C. Ding, Z. Fu, W. Guo, and F. Y. Wu,
initions. Nevertheless, there is significant room for im- Phys. Rev.E 81, 061111 (2010).
provement in the efficiency, as the rule used to choose [13] R. M. Ziff and H. Gu, Phys.Rev.E 79, 020102 (2009).
the bond for deletion and contraction at each step can [14] C.Tsallis,J.Phys.C:SolidStatePhys.15,L757(1982).