Pendulum

A pendulum, whose angle follows Newton's second law m l θ'' = -m g sin θ.

• Name: pendulum
• Field: physics/mechanics
• Formalism: ordinary differential equations
• References: Newton 1687

Factory

@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode', 'second order', 'oscillation'),
    references=('Newton 1687',),
    defaults={'until': 10.0, 'dt': 0.01},
)
def pendulum(
    theta: float = 0.5,
    omega: float = 0.0,
    *,
    m: float = 1.0,
    g: float = 9.81,
    l: float = 1.0,
) -> System:
    """Return a pendulum, whose angle follows Newton's second law `m l θ'' = -m g sin θ`.

    Args:
        theta:
            The initial angle, in radians.

        omega:
            The initial angular velocity.

        m:
            The mass of the bob.

        g:
            The gravitational acceleration.

        l:
            The length of the rod.

    Returns:
        system:
            The model, with the pendulum `bob`, whose angle has order 2 and angular velocity `bob.theta_d1`.

    Examples:
        >>> from skcomplex.models.physics.mechanics import pendulum
        >>> pendulum(theta=0.1).elements[0].theta_d1
        0.0
    """
    return System(
        elements=[PendulumElement('bob', theta=theta, theta_d1=omega, m=m, g=g, l=l)],
        interactions=[PendulumInteraction('swing', target='bob')],
    )

Components

class PendulumElement(Element):
    """Pendulum, a bob of mass `m` on a rigid rod of length `l` that swings by the angle `theta` under gravity `g`.

    The angle has order 2 and its first derivative `theta_d1` is the angular velocity. The motion comes from an
    interaction such as `PendulumInteraction`.

    Args:
        label:
            The label of the pendulum.

        visual:
            How the pendulum is drawn, or `None` for its class visual.

        theta:
            The initial angle from the vertical, in radians.

        theta_d1:
            The initial angular velocity.

        m:
            The mass of the bob.

        g:
            The gravitational acceleration.

        l:
            The length of the rod.

    Examples:
        >>> from skcomplex.models.physics.mechanics import PendulumElement
        >>> PendulumElement('bob', theta=0.2, l=2.0).l
        2.0
    """

    visual = VisualSpec('bodies', {'positions': _pendulum_bob, 'pivot': (0.0, 0.0), 'rods': True})

    def __init__(
        self,
        label: str,
        visual: VisualSpec | None = None,
        *,
        theta: Annotated[float, StateSpec(period=2 * pi, order=2)] = 0.5,
        theta_d1: float = 0.0,
        m: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
        g: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 9.81,
        l: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    ) -> None:
        super().__init__(label, visual=visual)
        self.theta = theta
        self.theta_d1 = theta_d1
        self.m = m
        self.g = g
        self.l = l

class PendulumInteraction(EulerLagrangeInteraction):
    """Motion of a pendulum, the Euler-Lagrange equation of `L = m l² θ'² / 2 + m g l cos θ`.

    The equation is Newton's law `m l² θ'' = -m g l sin θ` for the angle `θ` of the target, whose mass `m`, length
    `l`, and gravity `g` are parameters of the target.

    Args:
        label:
            The label of the interaction.

        target:
            The path of the pendulum.

    Examples:
        >>> from skcomplex.models.physics.mechanics import PendulumInteraction
        >>> PendulumInteraction('swing', target='bob').states
        'theta'
    """

    def __init__(self, label: str, *, target: Annotated[str, RoleSpec()]) -> None:
        super().__init__(label, target=target, states='theta')

    def _get_lagrangian(self, namespace: SimpleNamespace) -> Expression:
        """Return the Lagrangian `m l² θ'² / 2 + m g l cos θ` of the pendulum.

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

        Returns:
            lagrangian:
                The Lagrangian.
        """
        bob = namespace.target
        return bob.m * bob.l**2 * derivative(bob.theta, 1) ** 2 / 2 + bob.m * bob.g * bob.l * cos(bob.theta)

Equations

bob.theta: diff(theta_bob(t), t, t) == -g_bob*sin(theta_bob(t))/l_bob