Modified Chua with mixed-mode oscillations

A modified Chua circuit on three time scales, whose oscillations mix small and large amplitudes.

• Name: mixed_mode_chua
• Field: engineering/circuits
• Formalism: ordinary differential equations
• References: Krupa, Popović, and Kopell 2008; Trzaska 2023

Factory

@register_model(
    field='engineering/circuits',
    tags=('circuit', 'mixed-mode oscillations', 'canards', 'slow-fast', 'ode'),
    references=('Krupa, Popović, and Kopell 2008', 'Trzaska 2023'),
    defaults={'until': 200.0, 'dt': 0.5, 'solver': ScipyODESolver(method='Radau')},
)
def mixed_mode_chua(
    x: float = 0.01,
    y: float = 0.0,
    z: float = 0.0,
    *,
    epsilon: float = 0.01,
    a: float = 1.5,
    b: float = -1.0,
    alpha: float = 1.0,
    R: float = 0.0,
    beta: float = 0.005,
    c: float = 0.00053,
) -> System:
    """Return a modified Chua circuit on three time scales, whose oscillations mix small and large amplitudes.

    The circuit follows `ε x' = -y + (a + b x) x^2`, `y' = alpha (x - R y - z)`, and `z' = c - β y`. For `ε = 0.01`,
    `a = 1.5`, `b = -1`, and `β = 0.005`, the circuit has small oscillations only for `c = 0`, large ones only for
    `c = 0.00105`, and mixed-mode oscillations for `c = 0.00053`.

    Args:
        x:
            The initial voltage of the first capacitor.

        y:
            The initial current of the inductor.

        z:
            The initial voltage of the second capacitor.

        epsilon:
            The first capacitance, which sets the fast time scale.

        a:
            The quadratic coefficient of the diode.

        b:
            The cubic coefficient of the diode.

        alpha:
            The inverse inductance.

        R:
            The resistance.

        beta:
            The gain of the controlled source over the second capacitance.

        c:
            The offset of the controlled source.

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

    Examples:
        >>> from skcomplex.models.engineering.circuits import mixed_mode_chua
        >>> mixed_mode_chua(c=0.0).elements[0].c
        0.0
    """
    circuit = MixedModeChuaCircuitElement(
        'circuit',
        x=x,
        y=y,
        z=z,
        epsilon=epsilon,
        a=a,
        b=b,
        alpha=alpha,
        R=R,
        beta=beta,
        c=c,
    )
    return System(elements=[circuit])

Components

class MixedModeChuaCircuitElement(Element):
    """Slow-fast Chua circuit whose oscillations mix small and large amplitudes."""

    def __init__(
        self,
        label: str,
        visual: VisualSpec | None = None,
        *,
        x: Annotated[object, StateSpec()] = 0.01,
        y: Annotated[object, StateSpec()] = 0.0,
        z: Annotated[object, StateSpec()] = 0.0,
        epsilon: Annotated[object, ParameterSpec(domain=RealSet.open(0, oo))] = 0.01,
        a: Annotated[object, ParameterSpec()] = 1.5,
        b: Annotated[object, ParameterSpec()] = -1.0,
        alpha: Annotated[object, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
        R: Annotated[object, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
        beta: Annotated[object, ParameterSpec(domain=RealSet.open(0, oo))] = 0.005,
        c: Annotated[object, ParameterSpec()] = 0.00053,
    ) -> None:
        super().__init__(label, visual=visual)
        self.x = x
        self.y = y
        self.z = z
        self.epsilon = epsilon
        self.a = a
        self.b = b
        self.alpha = alpha
        self.R = R
        self.beta = beta
        self.c = c

    def _get_laws(self, s: SimpleNamespace) -> list[BaseLaw]:
        """Return the rates of the states of the circuit."""
        return [
            RateLaw(s.x, (-s.y + (s.a + s.b * s.x) * s.x**2) / s.epsilon),
            RateLaw(s.y, s.alpha * (s.x - s.R * s.y - s.z)),
            RateLaw(s.z, s.c - s.beta * s.y),
        ]

Equations

circuit.x: diff(x_circuit(t), t) == ((b_circuit*x_circuit(t) + a_circuit)*x_circuit(t)^2 - y_circuit(t))/epsilon_circuit
circuit.y: diff(y_circuit(t), t) == -(R_circuit*y_circuit(t) - x_circuit(t) + z_circuit(t))*alpha_circuit
circuit.z: diff(z_circuit(t), t) == -beta_circuit*y_circuit(t) + c_circuit