Stochastic fold with a ramped control parameter

A stochastic fold whose control is ramped through the tipping point of its lower branch.

• Name: ramped_stochastic_fold
• Field: physics/climate
• Formalism: stochastic differential equations
• References: Fang, Wang, Kurths, and Fan 2025; Scheffer et al. 2009

Factory

@register_model(
    field='physics/climate',
    tags=('tipping points', 'fold', 'ramp', 'sde', 'early warning signals', 'critical slowing down'),
    references=('Fang, Wang, Kurths, and Fan 2025', 'Scheffer et al. 2009'),
    defaults={'until': 500.0, 'dt': 0.5, 'solver': EulerMaruyamaSolver(dt=0.01)},
)
def ramped_stochastic_fold(
    *,
    rate: float = 0.004,
    sigma: float = 0.02,
    a: float = -1.0,
    b: float = 1.0,
    phi0: float = -1.0,
    phi1: float = 1.0,
) -> System:
    """Return a stochastic fold whose control is ramped through the tipping point of its lower branch.

    The state follows `dx = (a x^3 + b x + phi(t)) dt + sigma dW` with `phi = phi0 + rate t` until it reaches `phi1`,
    Eq. 4 of Fang et al. with `phi = -1 + 0.004 t`. It starts on the lower branch, which vanishes at the static fold
    `phi_c = sqrt(-4 b^3 / (27 a))`, `2 / (3 sqrt(3))` for `a = -1` and `b = 1`. Near the fold, recovery slows down,
    so the variance and the lag-1 autocorrelation of the state rise. A ramp faster than the relaxation carries the
    control past `phi_c` before the state tips, by a lag that grows with the rate.

    Args:
        rate:
            The rate at which the control rises.

        sigma:
            The amplitude of the noise.

        a:
            The cubic coefficient, negative.

        b:
            The linear coefficient, positive.

        phi0:
            The initial control.

        phi1:
            The final control.

    Returns:
        system:
            The model, with the state `system.x`, starting on the lower branch at `phi0`.

    Examples:
        >>> from skcomplex.models.physics.climate import ramped_stochastic_fold
        >>> ramped_stochastic_fold(rate=0.01).elements[0].rate
        0.01
    """
    fold = RampedFoldElement(
        'system',
        x=_lower_branch(a, b, phi0),
        a=a,
        b=b,
        sigma=sigma,
        phi0=phi0,
        phi1=phi1,
        rate=rate,
    )
    return System(elements=[fold])

Components

class RampedFoldElement(Element):
    """Stochastic fold `dx = (a x^3 + b x + phi(t)) dt + sigma dW` whose control rises linearly from `phi0` to `phi1`.

    For `a < 0 < b` the fold has two stable branches for `|phi| < phi_c = sqrt(-4 b^3 / (27 a))`, and the lower
    branch vanishes at `phi = phi_c`. The control rises at the rate `rate` and stays at `phi1` once it reaches it.

    Args:
        label:
            The label of the system.

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

        x:
            The initial state.

        a:
            The cubic coefficient, negative.

        b:
            The linear coefficient, positive.

        sigma:
            The amplitude of the noise.

        phi0:
            The initial control.

        phi1:
            The final control.

        rate:
            The rate at which the control rises.

    Examples:
        >>> from skcomplex.models.physics.climate import RampedFoldElement
        >>> RampedFoldElement('system').rate
        0.004
    """

    def __init__(
        self,
        label: str,
        visual: VisualSpec | None = None,
        *,
        x: Annotated[float, StateSpec()] = -1.3247,
        a: Annotated[float, ParameterSpec()] = -1.0,
        b: Annotated[float, ParameterSpec()] = 1.0,
        sigma: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.02,
        phi0: float = -1.0,
        phi1: float = 1.0,
        rate: float = 0.004,
    ) -> None:
        super().__init__(label, visual=visual)
        self.x = x
        self.a = a
        self.b = b
        self.sigma = sigma
        self.phi0 = phi0
        self.phi1 = phi1
        self.rate = rate

    def _get_laws(self, s: SimpleNamespace) -> list[BaseLaw]:
        """Return the fold with its ramped control and its noise.

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

        Returns:
            laws:
                The laws of the component.
        """
        control = self.phi0 + ramp(s.t, 0.0, (self.phi1 - self.phi0) / self.rate, self.phi1 - self.phi0)
        return [RateLaw(s.x, s.a * s.x**3 + s.b * s.x + control), NoiseLaw(s.x, s.sigma)]

Equations

system.x: diff(x_system(t), t) == a_system*x_system(t)^3 + b_system*x_system(t) + where(t < 0.0, 0, where(t < 500.0, 0.004*t, 2.0)) - 1.0