Site percolation

Site percolation probed by a fire lit along the top row, which spans the lattice above p_c = 0.5927.

• Name: site_percolation
• Field: physics/statistical_mechanics
• Formalism: lattice dynamics
• References: Broadbent and Hammersley 1957; Stauffer and Aharony 1994; Charbonneau 2017

Factory

@register_model(
    field='physics/statistical_mechanics',
    tags=('lattice', 'percolation', 'phase transition', 'fire'),
    references=('Broadbent and Hammersley 1957', 'Stauffer and Aharony 1994', 'Charbonneau 2017'),
    defaults={'until': 200, 'dt': 1},
)
def site_percolation(
    shape: Sequence[int] = (64, 64),
    *,
    p: float = 0.6,
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return site percolation probed by a fire lit along the top row, which spans the lattice above `p_c = 0.5927`.

    Sites are occupied with probability `p`, and the fire burns the clusters connected to the top row, one layer of
    sites per step.

    Args:
        shape:
            The number of rows and columns.

        p:
            The occupation probability.

        random_state:
            The seed or generator of the occupations.

    Returns:
        system:
            The model, with the covers `lattice.cover`.

    Examples:
        >>> from skcomplex.models.physics.statistical_mechanics import site_percolation
        >>> site_percolation((10, 10)).elements[0].check().prepare().size_
        100
    """
    shape = tuple(shape)
    cover = (np.random.default_rng(random_state).random(shape) < p).astype(float)
    cover[0] *= BURNING
    space = LatticeSpace(shape, boundary='fixed')
    return System(
        elements=[LatticeElement('lattice', element=PercolationSiteElement, space=space, values={'cover': cover})],
    )

Components

class PercolationSiteElement(Element):
    """Site of a percolation lattice: empty at 0, occupied at 1, burning at 2, or burnt at 3.

    Fire spreads from burning sites to their occupied neighbours and burns out in one step, so a fire lit along one
    edge burns exactly the clusters that touch that edge.

    Args:
        label:
            The label of the site or the lattice group.

        visual:
            How the element is drawn, or `None` for `default_visual_`.

        cover:
            The initial cover, 0, 1, 2, or 3.

    Examples:
        >>> from skcomplex.models.physics.statistical_mechanics import PercolationSiteElement
        >>> PercolationSiteElement('site').cover
        0.0
    """

    def __init__(
        self,
        label: str,
        visual: VisualSpec | None = None,
        *,
        cover: Annotated[float, StateSpec(domain=NN)] = 0.0,
    ) -> None:
        super().__init__(label, visual=visual)
        self.cover = cover

    def _get_laws(self, s: SimpleNamespace) -> list[BaseLaw]:
        """Return the spread and burning out of the fire.

        Args:
            s:
                The namespace of the symbols the laws can use.

        Returns:
            laws:
                The laws of the component.
        """
        return [LocalRuleLaw(s.cover, VonNeumannNeighborhood(center=True), PercolationSiteElement._burn)]

    @staticmethod
    def _burn(cover: Array, *around: Array) -> Array:
        """Return the next cover of a site, which burns out, catches fire from a burning neighbour, or stays."""
        catches = (cover == 1) & np.any(np.stack(around) == BURNING, axis=0)
        return np.where(cover == BURNING, BURNT, np.where(catches, BURNING, cover))