Fractional Van der Pol

The fractional Van der Pol oscillator D^alpha x = y, D^alpha y = mu (1 - x^2) y - x.

• Name: fractional_van_der_pol
• Field: mathematics/fractional_dynamics
• Formalism: fractional differential equations
• References: Podlubny 1999; Trzaska 2023

Factory

@register_model(
    field='mathematics/fractional_dynamics',
    tags=('fractional', 'memory', 'limit cycle'),
    references=('Podlubny 1999', 'Trzaska 2023'),
    defaults={'until': 40.0, 'dt': 0.1},
)
def fractional_van_der_pol(x: float = 2.0, y: float = 0.0, *, order: float = 0.9, mu: float = 1.0) -> System:
    """Return the fractional Van der Pol oscillator `D^alpha x = y`, `D^alpha y = mu (1 - x^2) y - x`.

    For orders close to 1 the trajectory approaches a closed cycle, and the memory of lower orders shrinks it toward
    the fixed point at the origin.

    Args:
        x:
            The initial displacement.

        y:
            The initial value of the displacement's derivative of order `alpha`.

        order:
            The order `alpha` of the Caputo derivatives, in `(0, 1]`.

        mu:
            The strength of the nonlinear damping.

    Returns:
        system:
            The model, with the states `oscillator.x` and `oscillator.y`.

    Examples:
        >>> from skcomplex.models.mathematics.fractional_dynamics import fractional_van_der_pol
        >>> fractional_van_der_pol(mu=2.0).elements[0].mu
        2.0
    """
    return System(elements=[van_der_pol_type(float(order))('oscillator', x=x, y=y, mu=mu)])

Components

    class FractionalElement(VanDerPolElement):
        def __init__(
            self,
            label: str,
            visual: VisualSpec | None = None,
            *,
            x: Annotated[float, StateSpec(order=order)] = 2.0,
            y: Annotated[float, StateSpec(order=order)] = 0.0,
            mu: float = 1.0,
        ) -> None:
            super().__init__(label, visual=visual, x=x, y=y, mu=mu)

Equations

oscillator.x: caputo(x_oscillator(t), t, 9/10) == y_oscillator(t)
oscillator.y: caputo(y_oscillator(t), t, 9/10) == -(x_oscillator(t)^2 - 1)*mu_oscillator*y_oscillator(t) - x_oscillator(t)