Swift-Hohenberg

The Swift-Hohenberg equation u_t = r u - (1 + laplacian)^2 u - u^3, which forms rolls of length 2 pi.

• Name: swift_hohenberg
• Field: physics/pattern_formation
• Formalism: ordinary differential equations
• References: Swift and Hohenberg 1977; Cross and Hohenberg 1993; Datseris and Parlitz 2022

Factory

@register_model(
    field='physics/pattern_formation',
    tags=('field', 'rolls', 'convection', 'spectral', 'pde'),
    references=('Swift and Hohenberg 1977', 'Cross and Hohenberg 1993', 'Datseris and Parlitz 2022'),
    defaults={'until': 200.0, 'dt': 5.0, 'solver': ETDRK4Solver(dt=0.5)},
)
def swift_hohenberg(
    shape: Sequence[int] = (64, 64),
    *,
    length: float = 16 * math.pi,
    r: float = 0.3,
    random_state: int | np.random.Generator | None = 0,
) -> System:
    """Return the Swift-Hohenberg equation `u_t = r u - (1 + laplacian)^2 u - u^3`, which forms rolls of length 2 pi.

    Above onset, `r > 0`, small noise grows into rolls whose wavenumber is near 1, the fastest growing one.

    Args:
        shape:
            The number of sites along every axis of the periodic lattice.

        length:
            The side of the domain.

        r:
            The distance above onset.

        random_state:
            The seed or generator of the noise.

    Returns:
        system:
            The model, with the field `cells.u`.

    Examples:
        >>> from skcomplex.models.physics.pattern_formation import swift_hohenberg
        >>> swift_hohenberg((8, 8), r=0.1).interactions[0].r
        0.1
    """
    shape = tuple(shape)
    values = {'u': _noise(shape, 0.01, random_state)}
    space = LatticeSpace(shape, spacing=length / shape[0])
    return System(
        elements=[
            FieldElement('cells', element=PatternFieldElement, space=space, values=values, discretization='spectral'),
        ],
        interactions=[SwiftHohenbergEnergyInteraction('energy', target='cells', r=r)],
    )

Components

class PatternFieldElement(Element):
    """Site of a real field `u` whose patterns a free energy selects.

    Args:
        label:
            The label of the site or the field.

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

        u:
            The initial value of the field.

    Attributes:
        default_visual_:
            The visual of the site when it sets none.

    Examples:
        >>> from skcomplex.models.physics.pattern_formation import PatternFieldElement
        >>> PatternFieldElement('site', u=0.1).u
        0.1
    """

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

class SwiftHohenbergEnergyInteraction(GradientInteraction):
    """Gradient flow `u_t = r u - (1 + laplacian)^2 u + g u^2 - u^3` of a field in the Swift-Hohenberg free energy.

    The free energy is `F = integral -r u^2 / 2 + ((1 + laplacian) u)^2 / 2 - g u^3 / 3 + u^4 / 4`. Above onset,
    `r > 0`, the modes of wavenumber near 1 grow into rolls, or into hexagons for a large enough quadratic term `g`.

    Args:
        label:
            The label of the interaction, unique among the parts of its system.

        target:
            The path of the field.

        state:
            The name of the state of the field.

        r:
            The distance above onset.

        g:
            The strength of the quadratic term, which breaks the symmetry `u -> -u`.

    Examples:
        >>> from types import SimpleNamespace
        >>> from sage.all__sagemath_symbolics import var
        >>> from skcomplex.models.physics.common.fields import SwiftHohenbergEnergyInteraction
        >>> energy = SwiftHohenbergEnergyInteraction('energy', target='cells').check().prepare()
        >>> namespace = SimpleNamespace(target=SimpleNamespace(u=var('u')), r=var('r'), g=0)
        >>> [law.value for law in energy.run(namespace).laws_]
        [-u^3 + r*u - u - 2*laplacian(u) - laplacian(laplacian(u))]
    """

    def __init__(
        self,
        label: str,
        *,
        target: Annotated[str, RoleSpec()],
        state: str = 'u',
        r: Annotated[float, ParameterSpec()] = 0.3,
        g: Annotated[float, ParameterSpec()] = 0.0,
    ) -> None:
        super().__init__(label, target=target, states=state, variational=True)
        self.state = state
        self.r = r
        self.g = g

    def _get_potential(self, namespace: SimpleNamespace) -> Expression:
        """Return the density `-r u^2 / 2 + (u + laplacian(u))^2 / 2 - g u^3 / 3 + u^4 / 4`.

        Args:
            namespace:
                The namespace of the symbols the density can use.

        Returns:
            density:
                The density of the free energy.
        """
        u = getattr(namespace.target, self.state)
        return -namespace.r * u**2 / 2 + (u + laplacian(u)) ** 2 / 2 - namespace.g * u**3 / 3 + u**4 / 4

Equations

cells.u[i]: diff(u_cells(t), t) == g_energy*u_cells(t)^2 - u_cells(t)^3 + r_energy*u_cells(t) - laplacian(laplacian(u_cells(t))) - 2*laplacian(u_cells(t)) - u_cells(t)