Percolation in lattices¶
Bond Percolation
This is an implementation of bond percolation. See Chapter 3 of [POG]. See also my blog post Percolation and self-avoiding walks related to this code.
AUTHORS:
Sebastien Labbe (2012-12-17): initial version, for pog lecture group
REFERENCES:
- POG
Geoffrey Grimmett, Probability on Graphs, http://www.statslab.cam.ac.uk/~grg/books/pgs.html
EXAMPLES:
Bond percolation sample¶
We construct a bond percolation sample in dimension d=2 with probability of open edges p=0.3:
sage: from slabbe import BondPercolationSample
sage: S = BondPercolationSample(p=0.3, d=2)
sage: S
Bond percolation sample d=2 p=0.300
An edge is defined uniquely as a starting point in Z^d and an axis direction given by an integer i such that 1 <= i <= d. One may ask if a given edge is in the sample S:
sage: ((34,56), 2) in S # random
False
The result is cached so the same answer is returned again:
sage: ((34,56), 2) in S # random
False
sage: ((34,56), 2) in S # random
False
The cluster containing the point zero is returned as an iterator:
sage: S.cluster()
<generator object at ...>
It may be finite of infinite. If you believe it is finite, you may compute its cardinality. If the cluster is infinite, it will not halt:
sage: S.cluster_cardinality() # not tested, might not halt
13
For larger values of p, the cluster might be larger if not infinite. In this case you may want to stop the computation at a certain point determined in advance. The following method does this. And it returns the cardinality if it is smaller than the stop value:
sage: S = BondPercolationSample(p=0.45, d=2)
sage: S.cluster_cardinality_stop_at(stop=10) # random
'>=10'
sage: S.cluster_cardinality_stop_at(stop=100) # random
'>=100'
sage: S.cluster_cardinality_stop_at(stop=1000) # random
625
Bond percolation samples¶
Construction of 20 bond percolation samples. For each of them, compute the cardinality of the open cluster containing zero:
sage: from slabbe import BondPercolationSamples
sage: S20 = BondPercolationSamples(p=0.4, d=2, n=20)
sage: S20.cluster_cardinality(stop=100) # random
[4, 2, 1, 4, 1, 10, 62, 71, 1, 25, 19, 2, 2, 42, '>=100', 1, 18, 2, '>=100', 20]
By considering “larger than 100” to be an infinite cluster, this gives a value of 2/20 = 0.10 for the percolation probability:
sage: S20.percolation_probability(stop=100) # random
0.100
By increasing the stop value, the computations can be redone again on the same samples. In this case, by using a stop value of 1000, we get a value of 0 for the percolation probability of p=0.4:
sage: S20.cluster_cardinality(stop=1000) # random
[4, 2, 1, 4, 1, 10, 62, 71, 1, 25, 19, 2, 2, 42, 176, 1, 18, 2, 186, 20]
sage: S20.percolation_probability(stop=1000) # random
0.000
Percolation probability¶
One can define the percolation probability function for a given dimension d. It will generate n samples and consider the cluster to be infinite if its cardinality is larger than the given stop value:
sage: from slabbe import PercolationProbability
sage: T = PercolationProbability(d=2, n=10, stop=100)
sage: T
Percolation Probability $\theta(p)$
d = dimension = 2
n = # samples = 10
stop counting at = 100
Compute the value for a certain probability p:
sage: T(0.4534) # random
0.300
Of course, this value will change for another equal percolation probability since it depends on the samples:
sage: T = PercolationProbability(d=2, n=10, stop=100)
sage: T(0.4534) # random
0.600
Anyway, it is useful to draw the plot of the percolation probability:
sage: T = PercolationProbability(d=2, n=10, stop=100)
sage: T.return_plot((0,1), adaptive_recursion=4, plot_points=4) # optional long
Graphics object consisting of 2 graphics primitives
Here we use Sage adaptative recursion algorithm for drawing plots which finds the particular important intervals to ask for more values of the function. See help section of plot function for details. Because T might be long to compute we start with only 4 points
TODO¶
Make it 100% doctested (presently 21/24 = 87%)
Base it on DiscreteSubset code
Fix tikz2pdf use
Do we want to use?:
sage: from sage.sets.set_from_iterator import EnumeratedSetFromIterator
sage: from itertools import count
sage: S = EnumeratedSetFromIterator(count)
sage: S
{0, 1, 2, 3, 4, ...}
and ?:
sage: M = FiniteSetMaps(["a", "b"], [3, 4, 5]); M
Maps from {'a', 'b'} to {3, 4, 5}
Methods and classes¶
-
class
slabbe.bond_percolation.
BondPercolationSample
(p, d=2)¶ Bases:
sage.structure.sage_object.SageObject
Let $L^d = (Z^d,E^d)$ be the hypercubic lattice.
A sample contained in the set {0,1}^{E^d}.
An edge e in E is open (=1) in the sample with probability p.
Cached __contains__ method below does the job of memory.
-
children
(pt)¶ Return an iterator over open neighbors of the point pt.
INPUT:
pt
- tuple, point in Z^dm
- integer, limit
EXAMPLES:
The result is consistent:
sage: from slabbe import BondPercolationSample sage: S = BondPercolationSample(0.5) sage: list(S.children((0,0))) # random [(1, 0), (-1, 0), (0, -1)]
Might be different for another sample:
sage: S = BondPercolationSample(0.5) sage: list(S.children((0,0))) # random [(-1, 0)]
In dimension 3:
sage: S = BondPercolationSample(0.5,3) sage: list(S.children((0,0,0))) # random [(1, 0, 0), (0, -1, 0), (0, 0, -1)] sage: S = BondPercolationSample(0.5,3) sage: list(S.children((0,0,0))) # random [(1, 0, 0), (-1, 0, 0)] sage: list(S.children((0,0,0))) # random [(1, 0, 0), (-1, 0, 0)]
sage: S = BondPercolationSample(1,2) sage: list(S.children((0,0))) # random [(1, 0), (-1, 0), (0, 1), (0, -1)] sage: S = BondPercolationSample(1,3) sage: list(S.children((0,0,0))) # random [(1, 0, 0), (-1, 0, 0), (0, 1, 0), (0, -1, 0), (0, 0, 1), (0, 0, -1)]
-
cluster
(pt=None)¶ Return an iterator over the open cluster containing the point pt.
INPUT:
pt
- tuple, point in Z^d. If None, pt=zero is considered.
EXAMPLES:
sage: from slabbe import BondPercolationSample sage: S = BondPercolationSample(0.5) sage: it = S.cluster() sage: next(it) (0, 0)
-
cluster_cardinality
(pt=None)¶ INPUT:
pt
- tuple, point in Z^d. If None, pt=zero is considered.
EXAMPLES:
sage: from slabbe import BondPercolationSample sage: BondPercolationSample(0.01).cluster_cardinality() # random 1 sage: BondPercolationSample(0.4).cluster_cardinality() # random 28
-
cluster_cardinality_stop_at
(stop, pt=None)¶ Return the cardinality of the cluster or the strin “>=STOP” if the size is larger than stop value.
INPUT:
stop
- integerpt
- tuple, point in Z^d. If None, pt=zero is considered.
EXAMPLES:
sage: from slabbe import BondPercolationSample sage: S = BondPercolationSample(0.3,2) sage: S.cluster_cardinality() # random 13 sage: S.cluster_cardinality_stop_at(1000) # random 13 sage: S.cluster_cardinality_stop_at(100) # random 13 sage: S.cluster_cardinality_stop_at(10) # random '>=10'
-
cluster_in_box
(m, pt=None)¶ Return the cluster (as a list) in the primal box [-m,m]^d containing the point pt.
INPUT:
m
- integerpt
- tuple, point in Z^d. If None, pt=zero is considered.
EXAMPLES:
sage: from slabbe import BondPercolationSample sage: S = BondPercolationSample(0.3,2) sage: S.cluster_in_box(2) # random [(-2, -2), (-2, -1), (-1, -2), (-1, -1), (-1, 0), (0, 0)]
-
edges_in_box
(m)¶ Return an iterator over all edges in the primal box [-m,m]^d.
INPUT:
m
- integer
EXAMPLES:
sage: from slabbe import BondPercolationSample sage: S = BondPercolationSample(1,2) sage: for a in S.edges_in_box(1): print(a) ((-1, -1), (0, -1)) ((-1, -1), (-1, 0)) ((-1, 0), (0, 0)) ((-1, 0), (-1, 1)) ((0, -1), (1, -1)) ((0, -1), (0, 0)) ((0, 0), (1, 0)) ((0, 0), (0, 1))
-
neighbor
(pt, d)¶ Return the neighbors of the point pt in direction d.
INPUT:
pt
- tuple, point in Z^ddirection
- integer, possible values are 1, 2, …, d and -1, -2, …, -d.
EXAMPLES:
sage: from slabbe import BondPercolationSample sage: S = BondPercolationSample(0.5,2) sage: S.neighbor((2,3),1) (3, 3) sage: S.neighbor((2,3),2) (2, 4) sage: S.neighbor((2,3),-1) (1, 3) sage: S.neighbor((2,3),-2) (2, 2)
-
plot
(m, pointsize=100, thickness=3, axes=False)¶ Return 2d graphics object contained in the primal box [-m,m]^d.
INPUT:
pointsize
, integer (default:100
),thickness
, integer (default:3
),axes
, bool (default:False
),
EXAMPLES:
sage: from slabbe import BondPercolationSample sage: S = BondPercolationSample(0.5,2) sage: S.plot(2) # optional long
It works in 3d!!:
sage: S = BondPercolationSample(0.5,3) sage: S.plot(3, pointsize=10, thickness=1) # optional long Graphics3d Object
-
save_pdf
(m)¶ EXAMPLES:
sage: from slabbe import BondPercolationSample sage: BondPercolationSample(0.3, d=2).save_pdf(20) # optional long Creation du fichier tikz_sample_d2_p300_m20.tikz Using template '/Users/slabbe/.tikz2pdf.tex'. tikz2pdf: calling pdflatex... tikz2pdf: Output written to 'tikz_sample_d2_p300_m20.pdf'.
-
tikz
(m)¶ Return tikz code.
EXAMPLES:
sage: from slabbe import BondPercolationSample sage: S = BondPercolationSample(0.5,2) sage: S.tikz(2) \begin{tikzpicture} [inner sep=0pt,thick, reddot/.style={fill=red,draw=red,circle,minimum size=5pt}] \clip (-2.4, -2.4) rectangle (2.4, 2.4); \draw (..., ...) -- (..., ...); ... \node[circle,fill=none,draw=red,minimum size=0.8cm,ultra thick,inner sep=0pt] at (0,0) {}; \node[above right] at (0,0) {$(0, 0)$}; \end{tikzpicture}
-
zero
()¶ EXAMPLES:
sage: from slabbe import BondPercolationSample sage: S = BondPercolationSample(0.5,3) sage: S.zero() (0, 0, 0)
sage: S = BondPercolationSample(0.5,5) sage: S.zero() (0, 0, 0, 0, 0)
-
-
class
slabbe.bond_percolation.
BondPercolationSamples
(p, d, n)¶ Bases:
sage.structure.sage_object.SageObject
Return a list of n BondPercolationSample of given parameter p and dimension d.
EXAMPLES:
sage: from slabbe import BondPercolationSamples sage: BondPercolationSamples(0.2,2,3) <slabbe.bond_percolation.BondPercolationSamples object at ...>
-
cluster_cardinality
(stop)¶ Return the list of cardinality of the cluster for each sample or +Infinity if the size is larger than stop value.
INPUT:
stop
- integer
EXAMPLES:
sage: from slabbe import BondPercolationSamples sage: S20 = BondPercolationSamples(p=0.2, d=2, n=20) sage: S20.cluster_cardinality(100) # random [1, 4, 1, 2, 2, 6, 5, 1, 5, 9, 2, 2, 2, 1, 2, 4, 4, 3, 2, 1]
sage: d = 2 sage: n = 5 sage: for p in srange(0,1,0.1): print(p,BondPercolationSamples(p,d,n).cluster_cardinality(100)) # optional long 0.000000000000000 [1, 1, 1, 1, 1] 0.100000000000000 [5, 1, 2, 1, 2] 0.200000000000000 [3, 1, 4, 1, 1] 0.300000000000000 [11, 7, 2, 1, 3] 0.400000000000000 [3, 3, 35, 10, '>=100'] 0.500000000000000 ['>=100', '>=100', '>=100', '>=100', 26] 0.600000000000000 ['>=100', '>=100', '>=100', '>=100', '>=100'] 0.700000000000000 ['>=100', '>=100', '>=100', '>=100', '>=100'] 0.800000000000000 ['>=100', '>=100', '>=100', '>=100', '>=100'] 0.900000000000000 ['>=100', '>=100', '>=100', '>=100', '>=100']
-
ntimes_over_size
(stop)¶ EXAMPLES:
sage: from slabbe import BondPercolationSamples sage: S = BondPercolationSamples(0.2,2,20) sage: S.ntimes_over_size(100) # random 0 sage: S = BondPercolationSamples(0.4,2,20) sage: S.ntimes_over_size(100) # random 1 sage: S = BondPercolationSamples(0.5,2,20) sage: S.ntimes_over_size(100) # random 17
-
percolation_probability
(stop)¶
-
-
class
slabbe.bond_percolation.
PercolationProbability
(d, n, stop, verbose=False)¶ Bases:
sage.structure.sage_object.SageObject
EXAMPLES:
sage: from slabbe import PercolationProbability sage: f = PercolationProbability(d=2, n=10, stop=100) sage: f Percolation Probability $\theta(p)$ d = dimension = 2 n = # samples = 10 stop counting at = 100
-
return_plot
(interval=0, 1, adaptive_recursion=4, plot_points=4, adaptive_tolerance=0.1)¶ Return a plot of percolation probability using basic sage plot settings.
INPUT:
interval
, default=(0,1)adaptive_recursion
, default=0plot_points
, default=10adaptive_tolerance
default=0.10
EXAMPLES:
sage: from slabbe import PercolationProbability sage: T = PercolationProbability(d=2, n=10, stop=100) sage: T.return_plot() # optional long Graphics object consisting of 1 graphics primitive
-
-
slabbe.bond_percolation.
compute_percolation_probability
(range_p, d, n, stop)¶ EXAMPLES:
sage: from slabbe.bond_percolation import compute_percolation_probability sage: compute_percolation_probability(srange(0,0.8,0.1), d=2, n=5, stop=100) # random d = 2, n = number of samples = 5 stop counting at = 100 p=0.0000, Theta=0.000, if |C|< 100 then max|C|=1 p=0.1000, Theta=0.000, if |C|< 100 then max|C|=1 p=0.2000, Theta=0.000, if |C|< 100 then max|C|=5 p=0.3000, Theta=0.000, if |C|< 100 then max|C|=6 p=0.4000, Theta=0.000, if |C|< 100 then max|C|=31 p=0.5000, Theta=1.00, if |C|< 100 then max|C|=-Infinity p=0.6000, Theta=1.00, if |C|< 100 then max|C|=-Infinity p=0.7000, Theta=1.00, if |C|< 100 then max|C|=-Infinity
sage: range_p = srange(0,0.8,0.1) sage: compute_percolation_probability(range_p, d=2, n=5, stop=100) # not tested d = 2, n = number of samples = 5 stop counting at = 100 p=0.0000, Theta=0.000, if |C|< 100 then max|C|=1 p=0.1000, Theta=0.000, if |C|< 100 then max|C|=1 p=0.2000, Theta=0.000, if |C|< 100 then max|C|=5 p=0.3000, Theta=0.000, if |C|< 100 then max|C|=6 p=0.4000, Theta=0.000, if |C|< 100 then max|C|=31 p=0.5000, Theta=1.00, if |C|< 100 then max|C|=-Infinity p=0.6000, Theta=1.00, if |C|< 100 then max|C|=-Infinity p=0.7000, Theta=1.00, if |C|< 100 then max|C|=-Infinity
sage: range_p = srange(0.45,0.55,0.01) sage: compute_percolation_probability(range_p, d=2, n=10, stop=1000) # not tested d = 2, n = number of samples = 10 stop counting at = 1000 p=0.4500, Theta=0.000, if |C|< 1000 then max|C|=378 p=0.4600, Theta=0.000, if |C|< 1000 then max|C|=475 p=0.4700, Theta=0.000, if |C|< 1000 then max|C|=514 p=0.4800, Theta=0.100, if |C|< 1000 then max|C|=655 p=0.4900, Theta=0.700, if |C|< 1000 then max|C|=274 p=0.5000, Theta=0.700, if |C|< 1000 then max|C|=975 p=0.5100, Theta=0.700, if |C|< 1000 then max|C|=16 p=0.5200, Theta=0.700, if |C|< 1000 then max|C|=125 p=0.5300, Theta=0.900, if |C|< 1000 then max|C|=4 p=0.5400, Theta=0.700, if |C|< 1000 then max|C|=6
sage: range_p = srange(0.475,0.485,0.001) sage: compute_percolation_probability(range_p, d=2, n=10, stop=1000) # not tested d = 2, n = number of samples = 10 stop counting at = 1000 p=0.4750, Theta=0.200, if |C|< 1000 then max|C|=718 p=0.4760, Theta=0.200, if |C|< 1000 then max|C|=844 p=0.4770, Theta=0.200, if |C|< 1000 then max|C|=566 p=0.4780, Theta=0.500, if |C|< 1000 then max|C|=257 p=0.4790, Theta=0.200, if |C|< 1000 then max|C|=566 p=0.4800, Theta=0.300, if |C|< 1000 then max|C|=544 p=0.4810, Theta=0.300, if |C|< 1000 then max|C|=778 p=0.4820, Theta=0.500, if |C|< 1000 then max|C|=983 p=0.4830, Theta=0.300, if |C|< 1000 then max|C|=473 p=0.4840, Theta=0.500, if |C|< 1000 then max|C|=411
sage: range_p = srange(0.47,0.48,0.001) sage: compute_percolation_probability(range_p, d=2, n=20, stop=2000) # not tested d = 2, n = number of samples = 20 stop counting at = 2000 p=0.4700, Theta=0.0500, if |C|< 2000 then max|C|=1666 p=0.4710, Theta=0.100, if |C|< 2000 then max|C|=1665 p=0.4720, Theta=0.000, if |C|< 2000 then max|C|=1798 p=0.4730, Theta=0.0500, if |C|< 2000 then max|C|=1717 p=0.4740, Theta=0.150, if |C|< 2000 then max|C|=1924 p=0.4750, Theta=0.150, if |C|< 2000 then max|C|=1893 p=0.4760, Theta=0.150, if |C|< 2000 then max|C|=1458 p=0.4770, Theta=0.150, if |C|< 2000 then max|C|=1573 p=0.4780, Theta=0.200, if |C|< 2000 then max|C|=1762 p=0.4790, Theta=0.250, if |C|< 2000 then max|C|=951
-
slabbe.bond_percolation.
percolation_graphics_array
(range_p, d, m, ncols=3)¶ EXAMPLES:
sage: from slabbe.bond_percolation import percolation_graphics_array sage: percolation_graphics_array(srange(0.1,1,0.1), d=2, m=5) # optional long sage: P = percolation_graphics_array(srange(0.45,0.55,0.01), d=2, m=5) # optional long sage: P.save('array_p45_p55_m5.png') # not tested sage: P = percolation_graphics_array(srange(0.45,0.55,0.01), d=2, m=10) # optional long sage: P.save('array_p45_p55_m10.png') # not tested