Skip to content

Mechanics

Mechanics: pendulums, oscillators, projectiles, collisions, rotating bodies, and potential wells.

BallisticPendulumElement(label, visual=None, *, x=-1.0, x_d1=0.0, theta=0.0, theta_d1=0.0, m=0.01, M=1.0, l=1.0, g=9.81)

Bases: Element

Bullet of mass m at the position x flying towards the block of mass M of a pendulum with the angle theta.

The position and the angle have order 2. The bullet flies freely, the pendulum swings, and the impact transfers the momentum of the bullet to the pendulum.

Parameters:

Name Type Description Default
label str

The label of the bullet and pendulum.

required
visual VisualSpec | None

How the bullet and pendulum are drawn, or None for none.

None
x Annotated[float, StateSpec(order=2)]

The initial position of the bullet, negative before the block.

-1.0
x_d1 float

The initial velocity of the bullet.

0.0
theta Annotated[float, StateSpec(period=2 * pi, order=2)]

The initial angle of the pendulum.

0.0
theta_d1 float

The initial angular velocity of the pendulum.

0.0
m Annotated[float, ParameterSpec(domain=open(0, oo))]

The mass of the bullet.

0.01
M Annotated[float, ParameterSpec(domain=open(0, oo))]

The mass of the block.

1.0
l Annotated[float, ParameterSpec(domain=open(0, oo))]

The length of the rod.

1.0
g Annotated[float, ParameterSpec(domain=open(0, oo))]

The gravitational acceleration.

9.81

Examples:

>>> from skcomplex.models.physics.mechanics import BallisticPendulumElement
>>> BallisticPendulumElement('pendulum', M=2.0).M
2.0
Source code in src/skcomplex/models/physics/mechanics/_ballistic_pendulum.py
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    x: Annotated[float, StateSpec(order=2)] = -1.0,
    x_d1: float = 0.0,
    theta: Annotated[float, StateSpec(period=2 * pi, order=2)] = 0.0,
    theta_d1: float = 0.0,
    m: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 0.01,
    M: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    l: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    g: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 9.81,
) -> None:
    super().__init__(label, visual=visual)
    self.x = x
    self.x_d1 = x_d1
    self.theta = theta
    self.theta_d1 = theta_d1
    self.m = m
    self.M = M
    self.l = l
    self.g = g

BouncingBallElement(label, visual=None, *, h=1.0, h_d1=0.0, g=9.81)

Bases: Element

Ball at the height h above a floor, falling under the gravitational acceleration g.

The height has order 2 and its first derivative h_d1 is the vertical velocity.

Parameters:

Name Type Description Default
label str

The label of the ball.

required
visual VisualSpec | None

How the ball is drawn, or None for none.

None
h Annotated[float, StateSpec(order=2)]

The initial height.

1.0
h_d1 float

The initial vertical velocity.

0.0
g Annotated[float, ParameterSpec()]

The gravitational acceleration.

9.81

Examples:

>>> from skcomplex.models.physics.mechanics import BouncingBallElement
>>> BouncingBallElement('ball', h=2.0).h
2.0
Source code in src/skcomplex/models/physics/mechanics/_bouncing_ball.py
39
40
41
42
43
44
45
46
47
48
49
50
51
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    h: Annotated[float, StateSpec(order=2)] = 1.0,
    h_d1: float = 0.0,
    g: Annotated[float, ParameterSpec()] = 9.81,
) -> None:
    super().__init__(label, visual=visual)
    self.h = h
    self.h_d1 = h_d1
    self.g = g

CartPendulumElement(label, visual=None, *, y=0.0, y_d1=0.0, theta=0.5, theta_d1=0.0, m1=1.0, m2=1.0, l=1.0, g=9.81, F=0.0)

Bases: Element

Cart on a horizontal track at the position y, with a pendulum of angle theta hanging from it.

The position and the angle have order 2. The motion comes from CartPendulumInteraction.

Parameters:

Name Type Description Default
label str

The label of the cart and pendulum.

required
visual VisualSpec | None

How the cart and pendulum are drawn, or None for its class visual.

None
y Annotated[float, StateSpec(order=2)]

The initial position of the cart.

0.0
y_d1 float

The initial velocity of the cart.

0.0
theta Annotated[float, StateSpec(period=2 * pi, order=2)]

The initial angle of the pendulum.

0.5
theta_d1 float

The initial angular velocity of the pendulum.

0.0
m1 Annotated[float, ParameterSpec(domain=open(0, oo))]

The mass of the cart.

1.0
m2 Annotated[float, ParameterSpec(domain=open(0, oo))]

The mass of the bob.

1.0
l Annotated[float, ParameterSpec(domain=open(0, oo))]

The length of the rod.

1.0
g Annotated[float, ParameterSpec(domain=open(0, oo))]

The gravitational acceleration.

9.81
F Annotated[float, ParameterSpec()]

A constant horizontal force on the cart.

0.0

Examples:

>>> from skcomplex.models.physics.mechanics import CartPendulumElement
>>> CartPendulumElement('cart', m1=5.0).m1
5.0
Source code in src/skcomplex/models/physics/mechanics/_cart_pendulum.py
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    y: Annotated[float, StateSpec(order=2)] = 0.0,
    y_d1: float = 0.0,
    theta: Annotated[float, StateSpec(period=2 * pi, order=2)] = 0.5,
    theta_d1: float = 0.0,
    m1: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    m2: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    l: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    g: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 9.81,
    F: Annotated[float, ParameterSpec()] = 0.0,
) -> None:
    super().__init__(label, visual=visual)
    self.y = y
    self.y_d1 = y_d1
    self.theta = theta
    self.theta_d1 = theta_d1
    self.m1 = m1
    self.m2 = m2
    self.l = l
    self.g = g
    self.F = F

CartPendulumInteraction(label, *, target)

Bases: EulerLagrangeInteraction

Motion of a cart and the pendulum that hangs from it, the Euler-Lagrange equations with a force on the cart.

The Lagrangian is (m1 + m2) y'² / 2 + m2 l y' θ' cos θ + m2 l² θ'² / 2 + m2 g l cos θ, and the force F of the target drives the position y of the cart.

Parameters:

Name Type Description Default
label str

The label of the interaction.

required
target Annotated[str, RoleSpec()]

The path of the cart and pendulum.

required

Examples:

>>> from skcomplex.models.physics.mechanics import CartPendulumInteraction
>>> CartPendulumInteraction('motion', target='cart').states
('y', 'theta')
Source code in src/skcomplex/models/physics/mechanics/_cart_pendulum_interaction.py
35
36
def __init__(self, label: str, *, target: Annotated[str, RoleSpec()]) -> None:
    super().__init__(label, target=target, states=('y', 'theta'), forcing=self._get_forces_on_coordinates)

ClockPendulumElement(label, visual=None, *, theta=0.2, theta_d1=0.0, x=0.0, x_d1=0.0, m=1.0, l=1.0, M=1.0, g=9.8, rho=0.1, lam=0.1, theta_c=0.1)

Bases: Element

Pendulum clock, a pendulum of angle theta hanging from a support at the horizontal position x.

The angle and the position have order 2. The motion comes from ClockPendulumInteraction.

Parameters:

Name Type Description Default
label str

The label of the clock.

required
visual VisualSpec | None

How the clock is drawn, or None for none.

None
theta Annotated[float, StateSpec(order=2)]

The initial angle of the pendulum.

0.2
theta_d1 float

The initial angular velocity of the pendulum.

0.0
x Annotated[float, StateSpec(order=2)]

The initial position of the support.

0.0
x_d1 float

The initial velocity of the support.

0.0
m Annotated[float, ParameterSpec(domain=open(0, oo))]

The mass of the bob.

1.0
l Annotated[float, ParameterSpec(domain=open(0, oo))]

The length of the pendulum.

1.0
M Annotated[float, ParameterSpec(domain=open(0, oo))]

The mass of the support.

1.0
g Annotated[float, ParameterSpec(domain=open(0, oo))]

The gravitational acceleration.

9.8
rho Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The damping of the support.

0.1
lam Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The strength of the escapement.

0.1
theta_c Annotated[float, ParameterSpec(domain=open(0, oo))]

The angle below which the escapement pumps energy in.

0.1

Examples:

>>> from skcomplex.models.physics.mechanics import ClockPendulumElement
>>> ClockPendulumElement('clock', lam=0.2).lam
0.2
Source code in src/skcomplex/models/physics/mechanics/_clock_pendulum.py
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    theta: Annotated[float, StateSpec(order=2)] = 0.2,
    theta_d1: float = 0.0,
    x: Annotated[float, StateSpec(order=2)] = 0.0,
    x_d1: float = 0.0,
    m: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    l: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    M: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    g: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 9.8,
    rho: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.1,
    lam: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.1,
    theta_c: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 0.1,
) -> None:
    super().__init__(label, visual=visual)
    self.theta = theta
    self.theta_d1 = theta_d1
    self.x = x
    self.x_d1 = x_d1
    self.m = m
    self.l = l
    self.M = M
    self.g = g
    self.rho = rho
    self.lam = lam
    self.theta_c = theta_c

ClockPendulumInteraction(label, *, target)

Bases: EulerLagrangeInteraction

Motion of a pendulum clock on a moving support, the Euler-Lagrange equations with the escapement.

The Lagrangian is (M + m) x'² / 2 + m l x' θ' cos θ + m l² θ'² / 2 + m g l cos θ. The Rayleigh dissipation λ l sgn(|θ| - θ_c) θ'² + rho x'² damps the support and gives the escapement 2 λ l sgn(|θ| - θ_c) θ', which pumps energy into the pendulum below the angle θ_c and takes it out beyond.

Parameters:

Name Type Description Default
label str

The label of the interaction.

required
target Annotated[str, RoleSpec()]

The path of the clock.

required

Examples:

>>> from skcomplex.models.physics.mechanics import ClockPendulumInteraction
>>> ClockPendulumInteraction('motion', target='clock').states
('theta', 'x')
Source code in src/skcomplex/models/physics/mechanics/_clock_pendulum_interaction.py
36
37
def __init__(self, label: str, *, target: Annotated[str, RoleSpec()]) -> None:
    super().__init__(label, target=target, states=('theta', 'x'), dissipation=self._get_escapement)

DiscBilliardElement(label, visual=None, *, x=0.3, y=0.0, vx=0.4, vy=1.0, radius=1.0)

Bases: Element

Particle moving freely in the plane with the position (x, y) and the velocity (vx, vy), inside a disc.

The position changes at the velocity, which changes only at the reflections off the wall of the disc.

Parameters:

Name Type Description Default
label str

The label of the particle.

required
visual VisualSpec | None

How the particle is drawn, or None for its class visual.

None
x Annotated[float, StateSpec()]

The initial horizontal position.

0.3
y Annotated[float, StateSpec()]

The initial vertical position.

0.0
vx Annotated[float, StateSpec()]

The initial horizontal velocity.

0.4
vy Annotated[float, StateSpec()]

The initial vertical velocity.

1.0
radius Annotated[float, ParameterSpec(domain=open(0, oo))]

The radius of the disc.

1.0

Examples:

>>> from skcomplex.models.physics.mechanics import DiscBilliardElement
>>> DiscBilliardElement('particle', radius=2.0).radius
2.0
Source code in src/skcomplex/models/physics/mechanics/_disc_billiard.py
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    x: Annotated[float, StateSpec()] = 0.3,
    y: Annotated[float, StateSpec()] = 0.0,
    vx: Annotated[float, StateSpec()] = 0.4,
    vy: Annotated[float, StateSpec()] = 1.0,
    radius: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
) -> None:
    super().__init__(label, visual=visual)
    self.x = x
    self.y = y
    self.vx = vx
    self.vy = vy
    self.radius = radius

DoublePendulumElement(label, visual=None, *, theta1=math.pi / 2, theta1_d1=0.0, theta2=math.pi / 2, theta2_d1=0.0, m1=1.0, m2=1.0, l1=1.0, l2=1.0, g=9.81)

Bases: Element

Double pendulum, a bob of mass m2 on a rod of length l2 hanging from a bob of mass m1 on a rod l1.

The angles theta1 and theta2 of the rods from the vertical have order 2. The motion comes from DoublePendulumInteraction.

Parameters:

Name Type Description Default
label str

The label of the double pendulum.

required
visual VisualSpec | None

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

None
theta1 Annotated[float, StateSpec(period=2 * pi, order=2)]

The initial angle of the upper rod.

pi / 2
theta1_d1 float

The initial angular velocity of the upper rod.

0.0
theta2 Annotated[float, StateSpec(period=2 * pi, order=2)]

The initial angle of the lower rod.

pi / 2
theta2_d1 float

The initial angular velocity of the lower rod.

0.0
m1 Annotated[float, ParameterSpec(domain=open(0, oo))]

The mass of the upper bob.

1.0
m2 Annotated[float, ParameterSpec(domain=open(0, oo))]

The mass of the lower bob.

1.0
l1 Annotated[float, ParameterSpec(domain=open(0, oo))]

The length of the upper rod.

1.0
l2 Annotated[float, ParameterSpec(domain=open(0, oo))]

The length of the lower rod.

1.0
g Annotated[float, ParameterSpec(domain=open(0, oo))]

The gravitational acceleration.

9.81

Examples:

>>> from skcomplex.models.physics.mechanics import DoublePendulumElement
>>> DoublePendulumElement('pendulum', m2=2.0).m2
2.0
Source code in src/skcomplex/models/physics/mechanics/_double_pendulum.py
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    theta1: Annotated[float, StateSpec(period=2 * pi, order=2)] = math.pi / 2,
    theta1_d1: float = 0.0,
    theta2: Annotated[float, StateSpec(period=2 * pi, order=2)] = math.pi / 2,
    theta2_d1: float = 0.0,
    m1: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    m2: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    l1: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    l2: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    g: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 9.81,
) -> None:
    super().__init__(label, visual=visual)
    self.theta1 = theta1
    self.theta1_d1 = theta1_d1
    self.theta2 = theta2
    self.theta2_d1 = theta2_d1
    self.m1 = m1
    self.m2 = m2
    self.l1 = l1
    self.l2 = l2
    self.g = g

DoublePendulumInteraction(label, *, target)

Bases: EulerLagrangeInteraction

Motion of a double pendulum, the Euler-Lagrange equations of its kinetic energy and gravity.

The Lagrangian is (m1 + m2) l1² θ1'² / 2 + m2 l2² θ2'² / 2 + m2 l1 l2 θ1' θ2' cos(θ1 - θ2) + (m1 + m2) g l1 cos θ1 + m2 g l2 cos θ2, with the masses, lengths, and gravity of the target.

Parameters:

Name Type Description Default
label str

The label of the interaction.

required
target Annotated[str, RoleSpec()]

The path of the double pendulum.

required

Examples:

>>> from skcomplex.models.physics.mechanics import DoublePendulumInteraction
>>> DoublePendulumInteraction('motion', target='pendulum').states
('theta1', 'theta2')
Source code in src/skcomplex/models/physics/mechanics/_double_pendulum_interaction.py
35
36
def __init__(self, label: str, *, target: Annotated[str, RoleSpec()]) -> None:
    super().__init__(label, target=target, states=('theta1', 'theta2'))

ElasticCollisionInteraction(label, *, left, right, diameter=0.1)

Bases: Interaction

Elastic collision of two equal masses in a row, which exchange their velocities when they touch.

The collision is the event at which the gap x_right - x_left - diameter between the balls falls to zero, and equal masses swap their velocities in an elastic collision.

Parameters:

Name Type Description Default
label str

The label of the interaction.

required
left Annotated[str, RoleSpec()]

The path of the left ball, with a position x of order 2.

required
right Annotated[str, RoleSpec()]

The path of the right ball, with a position x of order 2.

required
diameter Annotated[float, ParameterSpec(domain=open(0, oo))]

The diameter of the balls.

0.1

Examples:

>>> from skcomplex.models.physics.mechanics import ElasticCollisionInteraction
>>> ElasticCollisionInteraction('collision', left='a', right='b', diameter=0.2).diameter
0.2
Source code in src/skcomplex/models/physics/mechanics/_elastic_collision.py
45
46
47
48
49
50
51
52
53
54
55
56
def __init__(
    self,
    label: str,
    *,
    left: Annotated[str, RoleSpec()],
    right: Annotated[str, RoleSpec()],
    diameter: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 0.1,
) -> None:
    super().__init__(label)
    self.left = left
    self.right = right
    self.diameter = diameter

ExponentialWellInteraction(label, *, target, k=2.0)

Bases: PotentialEnergyInteraction

Potential k x² e^-x on a particle in one dimension, a well at x = 0 behind a barrier at x = 2.

Parameters:

Name Type Description Default
label str

The label of the interaction.

required
target Annotated[str, RoleSpec()]

The path of the particle.

required
k Annotated[float, ParameterSpec(domain=open(0, oo))]

The strength of the potential.

2.0

Examples:

>>> from skcomplex.models.physics.mechanics import ExponentialWellInteraction
>>> ExponentialWellInteraction('well', target='particle', k=2.0).k
2.0
Source code in src/skcomplex/models/physics/mechanics/_wells.py
38
39
40
41
42
43
44
45
46
def __init__(
    self,
    label: str,
    *,
    target: Annotated[str, RoleSpec()],
    k: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 2.0,
) -> None:
    super().__init__(label, target=target)
    self.k = k

GyroscopeElement(label, visual=None, *, alpha=0.0, alpha_d1=0.0, beta=0.3, beta_d1=0.0, gamma=0.0, gamma_d1=0.0, A=1.0, C=2.0, A1=0.5, B1=0.5, C1=0.5, A2=1.0, torque=1.0)

Bases: Element

Gyroscope in a gimbal: a rotor in an inner ring that turns in an outer ring about a fixed axis.

The angle alpha of the outer ring, beta of the inner ring, and gamma of the rotor have order 2. The motion comes from GyroscopeInteraction.

Parameters:

Name Type Description Default
label str

The label of the gyroscope.

required
visual VisualSpec | None

How the gyroscope is drawn, or None for none.

None
alpha Annotated[float, StateSpec(period=2 * pi, order=2)]

The initial angle of the outer ring.

0.0
alpha_d1 float

The initial angular velocity of the outer ring.

0.0
beta Annotated[float, StateSpec(order=2)]

The initial angle of the inner ring.

0.3
beta_d1 float

The initial angular velocity of the inner ring.

0.0
gamma Annotated[float, StateSpec(period=2 * pi, order=2)]

The initial angle of the rotor.

0.0
gamma_d1 float

The initial angular velocity of the rotor.

0.0
A Annotated[float, ParameterSpec(domain=open(0, oo))]

The equatorial moment of inertia of the rotor.

1.0
C Annotated[float, ParameterSpec(domain=open(0, oo))]

The polar moment of inertia of the rotor.

2.0
A1 Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The first principal moment of inertia of the inner ring.

0.5
B1 Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The second principal moment of inertia of the inner ring.

0.5
C1 Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The third principal moment of inertia of the inner ring.

0.5
A2 Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The moment of inertia of the outer ring about its axis.

1.0
torque Annotated[float, ParameterSpec()]

The torque on the axis of the inner ring.

1.0

Examples:

>>> from skcomplex.models.physics.mechanics import GyroscopeElement
>>> GyroscopeElement('gyroscope', gamma_d1=10.0).gamma_d1
10.0
Source code in src/skcomplex/models/physics/mechanics/_gyroscope.py
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    alpha: Annotated[float, StateSpec(period=2 * pi, order=2)] = 0.0,
    alpha_d1: float = 0.0,
    beta: Annotated[float, StateSpec(order=2)] = 0.3,
    beta_d1: float = 0.0,
    gamma: Annotated[float, StateSpec(period=2 * pi, order=2)] = 0.0,
    gamma_d1: float = 0.0,
    A: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    C: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 2.0,
    A1: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.5,
    B1: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.5,
    C1: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.5,
    A2: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
    torque: Annotated[float, ParameterSpec()] = 1.0,
) -> None:
    super().__init__(label, visual=visual)
    self.alpha = alpha
    self.alpha_d1 = alpha_d1
    self.beta = beta
    self.beta_d1 = beta_d1
    self.gamma = gamma
    self.gamma_d1 = gamma_d1
    self.A = A
    self.C = C
    self.A1 = A1
    self.B1 = B1
    self.C1 = C1
    self.A2 = A2
    self.torque = torque

GyroscopeInteraction(label, *, target)

Bases: EulerLagrangeInteraction

Motion of a gimballed gyroscope, the Euler-Lagrange equations of its kinetic energy with a torque on beta.

The Lagrangian is the kinetic energy [A2 + (A + A1) cos² β + (C + C1) sin² β] alpha'² / 2 + (A + B1) β'² / 2 + C alpha' gamma' sin β + C gamma'² / 2, and the torque of the target acts on the axis of the inner ring.

Parameters:

Name Type Description Default
label str

The label of the interaction.

required
target Annotated[str, RoleSpec()]

The path of the gyroscope.

required

Examples:

>>> from skcomplex.models.physics.mechanics import GyroscopeInteraction
>>> GyroscopeInteraction('motion', target='gyroscope').states
('alpha', 'beta', 'gamma')
Source code in src/skcomplex/models/physics/mechanics/_gyroscope_interaction.py
35
36
def __init__(self, label: str, *, target: Annotated[str, RoleSpec()]) -> None:
    super().__init__(label, target=target, states=('alpha', 'beta', 'gamma'), forcing=self._get_torques)

HangingBallElement(label, visual=None, *, x=0.0, x_d1=0.0, rest=0.0, g=9.81, length=1.0)

Bases: Element

Ball hanging on a string of length length, which swings with small amplitude about its rest position.

The horizontal position x has order 2, and its first derivative x_d1 is the velocity. For small amplitudes the ball moves in the potential g (x - rest)² / (2 length).

Parameters:

Name Type Description Default
label str

The label of the ball.

required
visual VisualSpec | None

How the ball is drawn, or None for none.

None
x Annotated[float, StateSpec(order=2)]

The initial horizontal position.

0.0
x_d1 float

The initial horizontal velocity.

0.0
rest Annotated[float, ParameterSpec()]

The rest position.

0.0
g Annotated[float, ParameterSpec()]

The gravitational acceleration.

9.81
length Annotated[float, ParameterSpec(domain=open(0, oo))]

The length of the string.

1.0

Examples:

>>> from skcomplex.models.physics.mechanics import HangingBallElement
>>> HangingBallElement('ball', rest=0.5).rest
0.5
Source code in src/skcomplex/models/physics/mechanics/_hanging_ball.py
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    x: Annotated[float, StateSpec(order=2)] = 0.0,
    x_d1: float = 0.0,
    rest: Annotated[float, ParameterSpec()] = 0.0,
    g: Annotated[float, ParameterSpec()] = 9.81,
    length: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
) -> None:
    super().__init__(label, visual=visual)
    self.x = x
    self.x_d1 = x_d1
    self.rest = rest
    self.g = g
    self.length = length

HarmonicForceInteraction(label, *, target, F=0.0, omega=1.0, state='x', divisor='m')

Bases: RateInteraction

Periodic force F cos(omega t) on a target, the second derivative F cos(omega t) / m of its position.

It is the rate F cos(omega t) with the mass as divisor, applied to every component of the position.

Parameters:

Name Type Description Default
label str

The label of the interaction.

required
target Annotated[str, RoleSpec()]

The path of the element or group the force acts on.

required
F Annotated[float, ParameterSpec()]

The amplitude of the force.

0.0
omega Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The angular frequency of the force.

1.0
state str

The name of the state of order 2 that holds the position of the target.

'x'
divisor str

The name of the state or parameter that holds the mass of the target.

'm'

Examples:

>>> from skcomplex.models.physics.mechanics import HarmonicForceInteraction
>>> HarmonicForceInteraction('drive', target='mass', F=0.5, omega=2.0).omega
2.0
Source code in src/skcomplex/models/physics/mechanics/_springs.py
50
51
52
53
54
55
56
57
58
59
60
61
62
def __init__(
    self,
    label: str,
    *,
    target: Annotated[str, RoleSpec()],
    F: Annotated[float, ParameterSpec()] = 0.0,
    omega: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
    state: str = 'x',
    divisor: str = 'm',
) -> None:
    super().__init__(label, target=target, state=state, divisor=divisor)
    self.F = F
    self.omega = omega

MagnetInteraction(label, *, target, center=0.0, gamma=1.0, d=0.3)

Bases: PotentialEnergyInteraction

Magnet below the plane of a pendulum bob, which attracts it with energy -gamma / sqrt(|x - center|**2 + d**2).

Parameters:

Name Type Description Default
label str

The label of the interaction.

required
target Annotated[str, RoleSpec()]

The path of the bob.

required
center Annotated[ArrayLike, ParameterSpec(vector=True)]

The position of the magnet in the plane.

0.0
gamma Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The strength of the magnet.

1.0
d Annotated[float, ParameterSpec(domain=open(0, oo))]

The height of the bob above the plane of the magnets.

0.3

Examples:

>>> from skcomplex.models.physics.mechanics import MagnetInteraction
>>> MagnetInteraction('magnet', target='bob', center=[1.0, 0.0]).gamma
1.0
Source code in src/skcomplex/models/physics/mechanics/_magnet.py
44
45
46
47
48
49
50
51
52
53
54
55
56
def __init__(
    self,
    label: str,
    *,
    target: Annotated[str, RoleSpec()],
    center: Annotated[npt.ArrayLike, ParameterSpec(vector=True)] = 0.0,
    gamma: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
    d: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 0.3,
) -> None:
    super().__init__(label, target=target)
    self.center = center
    self.gamma = gamma
    self.d = d

ParametricPendulumElement(label, visual=None, *, theta=0.1, theta_d1=0.0, g=9.81, l=1.0, epsilon=0.2, Omega=2 * math.sqrt(9.81), c=0.0)

Bases: Element

Pendulum whose suspension point oscillates vertically, with the angle theta of order 2.

The motion comes from ParametricPendulumInteraction.

Parameters:

Name Type Description Default
label str

The label of the pendulum.

required
visual VisualSpec | None

How the pendulum is drawn, or None for none.

None
theta Annotated[float, StateSpec(period=2 * pi, order=2)]

The initial angle.

0.1
theta_d1 float

The initial angular velocity.

0.0
g Annotated[float, ParameterSpec(domain=open(0, oo))]

The gravitational acceleration.

9.81
l Annotated[float, ParameterSpec(domain=open(0, oo))]

The length of the rod.

1.0
epsilon Annotated[float, ParameterSpec()]

The amplitude of the acceleration of the suspension point, relative to g.

0.2
Omega Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The angular frequency of the suspension point.

2 * sqrt(9.81)
c Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The friction coefficient.

0.0

Examples:

>>> from skcomplex.models.physics.mechanics import ParametricPendulumElement
>>> ParametricPendulumElement('swing', epsilon=0.5).epsilon
0.5
Source code in src/skcomplex/models/physics/mechanics/_parametric_pendulum.py
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    theta: Annotated[float, StateSpec(period=2 * pi, order=2)] = 0.1,
    theta_d1: float = 0.0,
    g: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 9.81,
    l: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    epsilon: Annotated[float, ParameterSpec()] = 0.2,
    Omega: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 2 * math.sqrt(9.81),
    c: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
) -> None:
    super().__init__(label, visual=visual)
    self.theta = theta
    self.theta_d1 = theta_d1
    self.g = g
    self.l = l
    self.epsilon = epsilon
    self.Omega = Omega
    self.c = c

ParametricPendulumInteraction(label, *, target)

Bases: EulerLagrangeInteraction

Motion of a parametric pendulum, the Euler-Lagrange equation of θ'² / 2 + (g / l)(1 + ε cos Ω t) cos θ.

The Rayleigh dissipation c θ'² / 2 adds the friction, so the angle follows θ'' + c θ' + (g / l)(1 + ε cos Ω t) sin θ = 0, with the parameters of the target.

Parameters:

Name Type Description Default
label str

The label of the interaction.

required
target Annotated[str, RoleSpec()]

The path of the pendulum.

required

Examples:

>>> from skcomplex.models.physics.mechanics import ParametricPendulumInteraction
>>> ParametricPendulumInteraction('swing', target='pendulum').states
'theta'
Source code in src/skcomplex/models/physics/mechanics/_parametric_pendulum_interaction.py
35
36
def __init__(self, label: str, *, target: Annotated[str, RoleSpec()]) -> None:
    super().__init__(label, target=target, states='theta', dissipation=self._get_friction)

PendulumElement(label, visual=None, *, theta=0.5, theta_d1=0.0, m=1.0, g=9.81, l=1.0)

Bases: 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.

Parameters:

Name Type Description Default
label str

The label of the pendulum.

required
visual VisualSpec | None

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

None
theta Annotated[float, StateSpec(period=2 * pi, order=2)]

The initial angle from the vertical, in radians.

0.5
theta_d1 float

The initial angular velocity.

0.0
m Annotated[float, ParameterSpec(domain=open(0, oo))]

The mass of the bob.

1.0
g Annotated[float, ParameterSpec(domain=open(0, oo))]

The gravitational acceleration.

9.81
l Annotated[float, ParameterSpec(domain=open(0, oo))]

The length of the rod.

1.0

Examples:

>>> from skcomplex.models.physics.mechanics import PendulumElement
>>> PendulumElement('bob', theta=0.2, l=2.0).l
2.0
Source code in src/skcomplex/models/physics/mechanics/_pendulum.py
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
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

PendulumInteraction(label, *, target)

Bases: 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.

Parameters:

Name Type Description Default
label str

The label of the interaction.

required
target Annotated[str, RoleSpec()]

The path of the pendulum.

required

Examples:

>>> from skcomplex.models.physics.mechanics import PendulumInteraction
>>> PendulumInteraction('swing', target='bob').states
'theta'
Source code in src/skcomplex/models/physics/mechanics/_pendulum_interaction.py
37
38
def __init__(self, label: str, *, target: Annotated[str, RoleSpec()]) -> None:
    super().__init__(label, target=target, states='theta')

PendulumSpringInteraction(label, *, first, second, k=1.0)

Bases: EulerLagrangeInteraction

Spring between the bobs of two pendulums, which stretches with the horizontal distance of the bobs.

The spring of stiffness k adds the Lagrangian -k (l2 sin θ2 - l1 sin θ1)² / 2 to the Euler-Lagrange equations of the two angles, so its torque on every pendulum divides by the moment of inertia m l² of the pendulum's own Lagrangian.

Parameters:

Name Type Description Default
label str

The label of the interaction.

required
first Annotated[str, RoleSpec()]

The path of the first pendulum, with the angle theta, the mass m, and the length l.

required
second Annotated[str, RoleSpec()]

The path of the second pendulum.

required
k Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The stiffness of the spring.

1.0

Examples:

>>> from skcomplex.models.physics.mechanics import PendulumSpringInteraction
>>> PendulumSpringInteraction('spring', first='left', second='right', k=2.0).k
2.0
Source code in src/skcomplex/models/physics/mechanics/_pendulum_spring.py
42
43
44
45
46
47
48
49
50
51
52
53
def __init__(
    self,
    label: str,
    *,
    first: Annotated[str, RoleSpec()],
    second: Annotated[str, RoleSpec()],
    k: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
) -> None:
    super().__init__(label, target=first, states='theta')
    self.first = first
    self.second = second
    self.k = k

RotatingBallElement(label, visual=None, *, theta=0.0, m=1.0, r=1.0, beta=math.pi / 2, H=1.0)

Bases: Element

Ball of mass m on a rod at the distance r from the joint, turning about a vertical axis by the angle theta.

The rod is fixed at the angle beta to the axis, so the ball is at the distance r sin β from it, and H is its angular momentum about the axis.

Parameters:

Name Type Description Default
label str

The label of the ball.

required
visual VisualSpec | None

How the ball is drawn, or None for none.

None
theta Annotated[float, StateSpec(period=2 * pi)]

The initial angle about the axis.

0.0
m Annotated[float, ParameterSpec(domain=open(0, oo))]

The mass of the ball.

1.0
r Annotated[float, ParameterSpec(domain=open(0, oo))]

The distance of the ball from the joint of the rod.

1.0
beta Annotated[float, ParameterSpec(domain=open(0, oo))]

The angle of the rod to the axis, in radians.

pi / 2
H Annotated[float, ParameterSpec()]

The angular momentum of the ball about the axis.

1.0

Examples:

>>> from skcomplex.models.physics.mechanics import RotatingBallElement
>>> RotatingBallElement('ball', H=2.0).H
2.0
Source code in src/skcomplex/models/physics/mechanics/_rotating_ball.py
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    theta: Annotated[float, StateSpec(period=2 * pi)] = 0.0,
    m: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    r: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    beta: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = math.pi / 2,
    H: Annotated[float, ParameterSpec()] = 1.0,
) -> None:
    super().__init__(label, visual=visual)
    self.theta = theta
    self.m = m
    self.r = r
    self.beta = beta
    self.H = H

RotatingSwingElement(label, visual=None, *, phi=0.0, phi_d1=0.0, theta=0.1, theta_d1=0.0, J=10.0, m=1.0, n_chairs=8.0, a=1.0, l=1.0, g=9.81, torque=0.0, friction=0.0)

Bases: Element

Rotating swing: chairs on chains hanging from the rods of a stand that turns about a vertical axis.

The angle phi of the stand and the angle theta of the chains from the vertical have order 2. The motion comes from RotatingSwingInteraction.

Parameters:

Name Type Description Default
label str

The label of the swing.

required
visual VisualSpec | None

How the swing is drawn, or None for none.

None
phi Annotated[float, StateSpec(period=2 * pi, order=2)]

The initial angle of the stand.

0.0
phi_d1 float

The initial angular velocity of the stand.

0.0
theta Annotated[float, StateSpec(order=2)]

The initial angle of the chains.

0.1
theta_d1 float

The initial angular velocity of the chains.

0.0
J Annotated[float, ParameterSpec(domain=open(0, oo))]

The moment of inertia of the stand.

10.0
m Annotated[float, ParameterSpec(domain=open(0, oo))]

The mass of every chair.

1.0
n_chairs Annotated[float, ParameterSpec(domain=open(0, oo))]

The number of chairs.

8.0
a Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The length of the rods.

1.0
l Annotated[float, ParameterSpec(domain=open(0, oo))]

The length of the chains.

1.0
g Annotated[float, ParameterSpec(domain=open(0, oo))]

The gravitational acceleration.

9.81
torque Annotated[float, ParameterSpec()]

The torque that drives the stand.

0.0
friction Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The friction coefficient of the stand's axis.

0.0

Examples:

>>> from skcomplex.models.physics.mechanics import RotatingSwingElement
>>> RotatingSwingElement('swing', n_chairs=4.0).n_chairs
4.0
Source code in src/skcomplex/models/physics/mechanics/_rotating_swing.py
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
def __init__(
    self,
    label: str,
    visual: VisualSpec | None = None,
    *,
    phi: Annotated[float, StateSpec(period=2 * pi, order=2)] = 0.0,
    phi_d1: float = 0.0,
    theta: Annotated[float, StateSpec(order=2)] = 0.1,
    theta_d1: float = 0.0,
    J: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 10.0,
    m: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    n_chairs: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 8.0,
    a: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
    l: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 1.0,
    g: Annotated[float, ParameterSpec(domain=RealSet.open(0, oo))] = 9.81,
    torque: Annotated[float, ParameterSpec()] = 0.0,
    friction: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 0.0,
) -> None:
    super().__init__(label, visual=visual)
    self.phi = phi
    self.phi_d1 = phi_d1
    self.theta = theta
    self.theta_d1 = theta_d1
    self.J = J
    self.m = m
    self.n_chairs = n_chairs
    self.a = a
    self.l = l
    self.g = g
    self.torque = torque
    self.friction = friction

RotatingSwingInteraction(label, *, target)

Bases: EulerLagrangeInteraction

Motion of a rotating swing, the Euler-Lagrange equations with a torque and friction on the stand.

The Lagrangian is [J + n m (a + l sin θ)²] φ'² / 2 + n m l² θ'² / 2 + n m g l cos θ, the torque drives the angle φ of the stand, and the Rayleigh dissipation c φ'² / 2 adds the friction of its axis.

Parameters:

Name Type Description Default
label str

The label of the interaction.

required
target Annotated[str, RoleSpec()]

The path of the swing.

required

Examples:

>>> from skcomplex.models.physics.mechanics import RotatingSwingInteraction
>>> RotatingSwingInteraction('motion', target='swing').states
('phi', 'theta')
Source code in src/skcomplex/models/physics/mechanics/_rotating_swing_interaction.py
35
36
37
38
39
40
41
42
def __init__(self, label: str, *, target: Annotated[str, RoleSpec()]) -> None:
    super().__init__(
        label,
        target=target,
        states=('phi', 'theta'),
        dissipation=self._get_friction,
        forcing=self._get_torques,
    )

SupportSpringInteraction(label, *, first, second, k=1.0)

Bases: EulerLagrangeInteraction

Spring of stiffness k between the supports of two pendulum clocks.

The spring adds the Lagrangian -k (x2 - x1)² / 2 to the Euler-Lagrange equations of the supports, so every clock answers with the accelerations of its support and its pendulum through the mass matrix of its own Lagrangian.

Parameters:

Name Type Description Default
label str

The label of the interaction.

required
first Annotated[str, RoleSpec()]

The path of the first clock.

required
second Annotated[str, RoleSpec()]

The path of the second clock.

required
k Annotated[float, ParameterSpec(domain=closed_open(0, oo))]

The stiffness of the spring.

1.0

Examples:

>>> from skcomplex.models.physics.mechanics import SupportSpringInteraction
>>> SupportSpringInteraction('spring', first='clock_0', second='clock_1', k=2.0).k
2.0
Source code in src/skcomplex/models/physics/mechanics/_support_spring.py
41
42
43
44
45
46
47
48
49
50
51
52
def __init__(
    self,
    label: str,
    *,
    first: Annotated[str, RoleSpec()],
    second: Annotated[str, RoleSpec()],
    k: Annotated[float, ParameterSpec(domain=RealSet.closed_open(0, oo))] = 1.0,
) -> None:
    super().__init__(label, target=first, states='x')
    self.first = first
    self.second = second
    self.k = k

ballistic_pendulum(speed=200.0, distance=1.0, *, m=0.01, M=1.0, l=1.0, g=9.81)

Return a bullet that sticks in the block of a pendulum at rest, which then swings up.

The bullet flies horizontally towards the block, and the impact is an event. Momentum is conserved in the perfectly inelastic collision, so the pendulum starts with the speed m v / (m + M) and rises to the height (m v / (m + M))² / (2 g). The bullet position stops at the block after the impact.

Parameters:

Name Type Description Default
speed float

The speed of the bullet.

200.0
distance float

The initial distance of the bullet from the block.

1.0
m float

The mass of the bullet.

0.01
M float

The mass of the block.

1.0
l float

The length of the rod.

1.0
g float

The gravitational acceleration.

9.81

Returns:

Name Type Description
system System

The model, with the bullet position pendulum.x and the angle pendulum.theta, both of order 2.

Examples:

>>> from skcomplex.models.physics.mechanics import ballistic_pendulum
>>> ballistic_pendulum(distance=2.0).elements[0].x
-2.0
Source code in src/skcomplex/models/physics/mechanics/_pendulums.py
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
427
428
429
430
431
432
433
434
435
436
437
438
439
440
441
442
443
444
445
446
447
448
449
450
451
452
453
454
455
456
457
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'events', 'impacts', 'momentum', 'hybrid'),
    references=('Ryzhov et al. 2021',),
    defaults={'until': 3.0, 'dt': 0.01},
)
def ballistic_pendulum(
    speed: float = 200.0,
    distance: float = 1.0,
    *,
    m: float = 0.01,
    M: float = 1.0,
    l: float = 1.0,
    g: float = 9.81,
) -> System:
    """Return a bullet that sticks in the block of a pendulum at rest, which then swings up.

    The bullet flies horizontally towards the block, and the impact is an event. Momentum is conserved in the
    perfectly inelastic collision, so the pendulum starts with the speed `m v / (m + M)` and rises to the height
    `(m v / (m + M))² / (2 g)`. The bullet position stops at the block after the impact.

    Args:
        speed:
            The speed of the bullet.

        distance:
            The initial distance of the bullet from the block.

        m:
            The mass of the bullet.

        M:
            The mass of the block.

        l:
            The length of the rod.

        g:
            The gravitational acceleration.

    Returns:
        system:
            The model, with the bullet position `pendulum.x` and the angle `pendulum.theta`, both of order 2.

    Examples:
        >>> from skcomplex.models.physics.mechanics import ballistic_pendulum
        >>> ballistic_pendulum(distance=2.0).elements[0].x
        -2.0
    """
    swing = EulerLagrangeInteraction(
        'swing',
        target='pendulum',
        states='theta',
        lagrangian=lambda s: s.target.l * derivative(s.target.theta, 1) ** 2 / 2 + s.target.g * cos(s.target.theta),
    )
    impact = CrossingInteraction(
        'impact',
        target='pendulum',
        guard=lambda s: s.target.x,
        reset={
            'x_d1': lambda s: 0,
            'theta_d1': lambda s: s.target.m * s.target.x_d1 / ((s.target.m + s.target.M) * s.target.l),
        },
        direction=1,
    )
    return System(
        elements=[BallisticPendulumElement('pendulum', x=-distance, x_d1=speed, m=m, M=M, l=l, g=g)],
        interactions=[swing, impact],
    )

bouncing_ball(h=1.0, v=0.0, *, g=9.81, e=0.8)

Return a ball that falls under gravity and bounces off the floor, losing speed by a coefficient of restitution.

Every impact is an event that reverses the velocity and scales it by e. The impacts accumulate at the Zeno time sqrt(2 h / g) (1 + e) / (1 - e) for a ball dropped from rest, so runs should end before it.

Parameters:

Name Type Description Default
h float

The initial height.

1.0
v float

The initial vertical velocity.

0.0
g float

The gravitational acceleration.

9.81
e float

The coefficient of restitution.

0.8

Returns:

Name Type Description
system System

The model, with the height ball.h of order 2 and its velocity ball.h_d1.

Examples:

>>> from skcomplex.models.physics.mechanics import bouncing_ball
>>> from skcomplex.engine.runtime import Simulation
>>> trajectory = Simulation(bouncing_ball(), dt=0.01).check().prepare().run(until=1.0).trajectory_
>>> round(trajectory.events[0].time, 6)
0.451524
Source code in src/skcomplex/models/physics/mechanics/_mechanics.py
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'events', 'impacts'),
    references=('Strogatz 2015', 'Hairer, Nørsett, and Wanner 1993'),
    defaults={'until': 3.0, 'dt': 0.01},
)
def bouncing_ball(h: float = 1.0, v: float = 0.0, *, g: float = 9.81, e: float = 0.8) -> System:
    """Return a ball that falls under gravity and bounces off the floor, losing speed by a coefficient of restitution.

    Every impact is an event that reverses the velocity and scales it by `e`. The impacts accumulate at the Zeno time
    `sqrt(2 h / g) (1 + e) / (1 - e)` for a ball dropped from rest, so runs should end before it.

    Args:
        h:
            The initial height.

        v:
            The initial vertical velocity.

        g:
            The gravitational acceleration.

        e:
            The coefficient of restitution.

    Returns:
        system:
            The model, with the height `ball.h` of order 2 and its velocity `ball.h_d1`.

    Examples:
        >>> from skcomplex.models.physics.mechanics import bouncing_ball
        >>> from skcomplex.engine.runtime import Simulation
        >>> trajectory = Simulation(bouncing_ball(), dt=0.01).check().prepare().run(until=1.0).trajectory_
        >>> round(trajectory.events[0].time, 6)
        0.451524
    """
    return System(
        elements=[BouncingBallElement('ball', h=h, h_d1=v, g=g)],
        interactions=[
            RateInteraction('gravity', target='ball', state='h', rate=lambda s: -s.target.g),
            ReflectionInteraction(
                'impact',
                target='ball',
                guard=lambda s: s.target.h,
                positions=('h',),
                velocities=('h_d1',),
                direction=-1,
                restitution=e,
            ),
        ],
    )

cart_pendulum(theta=0.5, y=0.0, *, m1=1.0, m2=1.0, l=1.0, g=9.81, F=0.0)

Return a pendulum hanging from a cart that slides without friction along a horizontal track.

The Lagrange equations (m1 + m2) y'' + m2 l (θ'' cos θ - θ'² sin θ) = F and l θ'' + y'' cos θ + g sin θ = 0 couple the position y of the cart and the angle θ of the pendulum. Without a force the horizontal momentum is conserved, so the cart swings against the pendulum.

Parameters:

Name Type Description Default
theta float

The initial angle of the pendulum.

0.5
y float

The initial position of the cart.

0.0
m1 float

The mass of the cart.

1.0
m2 float

The mass of the bob.

1.0
l float

The length of the rod.

1.0
g float

The gravitational acceleration.

9.81
F float

A constant horizontal force on the cart.

0.0

Returns:

Name Type Description
system System

The model, with the cart position cart.y and the angle cart.theta, both of order 2.

Examples:

>>> from skcomplex.models.physics.mechanics import cart_pendulum
>>> cart_pendulum(m1=5.0).elements[0].m1
5.0
Source code in src/skcomplex/models/physics/mechanics/_pendulums.py
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode', 'second order', 'lagrangian', 'momentum'),
    references=('Ryzhov et al. 2021',),
    defaults={'until': 10.0, 'dt': 0.02},
)
def cart_pendulum(
    theta: float = 0.5,
    y: float = 0.0,
    *,
    m1: float = 1.0,
    m2: float = 1.0,
    l: float = 1.0,
    g: float = 9.81,
    F: float = 0.0,
) -> System:
    """Return a pendulum hanging from a cart that slides without friction along a horizontal track.

    The Lagrange equations `(m1 + m2) y'' + m2 l (θ'' cos θ - θ'² sin θ) = F` and `l θ'' + y'' cos θ + g sin θ = 0`
    couple the position `y` of the cart and the angle `θ` of the pendulum. Without a force the horizontal momentum is
    conserved, so the cart swings against the pendulum.

    Args:
        theta:
            The initial angle of the pendulum.

        y:
            The initial position of the cart.

        m1:
            The mass of the cart.

        m2:
            The mass of the bob.

        l:
            The length of the rod.

        g:
            The gravitational acceleration.

        F:
            A constant horizontal force on the cart.

    Returns:
        system:
            The model, with the cart position `cart.y` and the angle `cart.theta`, both of order 2.

    Examples:
        >>> from skcomplex.models.physics.mechanics import cart_pendulum
        >>> cart_pendulum(m1=5.0).elements[0].m1
        5.0
    """
    return System(
        elements=[CartPendulumElement('cart', y=y, theta=theta, m1=m1, m2=m2, l=l, g=g, F=F)],
        interactions=[CartPendulumInteraction('motion', target='cart')],
    )

circular_billiard(x=0.3, y=0.0, vx=0.4, vy=1.0, *, radius=1.0)

Return a particle moving freely inside a circle and reflecting specularly off its wall.

The angle of incidence is the same at every bounce, so the circular billiard is integrable.

Parameters:

Name Type Description Default
x float

The initial horizontal position.

0.3
y float

The initial vertical position.

0.0
vx float

The initial horizontal velocity.

0.4
vy float

The initial vertical velocity.

1.0
radius float

The radius of the table.

1.0

Returns:

Name Type Description
system System

The model, with the position particle.x, particle.y and the velocity particle.vx, particle.vy.

Examples:

>>> from skcomplex.models.physics.mechanics import circular_billiard
>>> circular_billiard(radius=2.0).elements[0].radius
2.0
Source code in src/skcomplex/models/physics/mechanics/_mechanics.py
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
387
388
389
390
391
392
393
394
395
396
397
398
399
400
401
402
403
404
405
406
407
408
409
410
411
412
413
414
415
416
417
418
419
420
421
422
423
424
425
426
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'events', 'billiard', 'integrable'),
    references=('Tabachnikov 2005', 'Datseris and Parlitz 2022'),
    defaults={'until': 20.0, 'dt': 0.05},
)
def circular_billiard(
    x: float = 0.3,
    y: float = 0.0,
    vx: float = 0.4,
    vy: float = 1.0,
    *,
    radius: float = 1.0,
) -> System:
    """Return a particle moving freely inside a circle and reflecting specularly off its wall.

    The angle of incidence is the same at every bounce, so the circular billiard is integrable.

    Args:
        x:
            The initial horizontal position.

        y:
            The initial vertical position.

        vx:
            The initial horizontal velocity.

        vy:
            The initial vertical velocity.

        radius:
            The radius of the table.

    Returns:
        system:
            The model, with the position `particle.x`, `particle.y` and the velocity `particle.vx`, `particle.vy`.

    Examples:
        >>> from skcomplex.models.physics.mechanics import circular_billiard
        >>> circular_billiard(radius=2.0).elements[0].radius
        2.0
    """
    wall = ReflectionInteraction(
        'bounce',
        target='particle',
        guard=lambda s: s.target.x**2 + s.target.y**2 - s.target.radius**2,
        positions=('x', 'y'),
        velocities=('vx', 'vy'),
        direction=1,
    )
    return System(
        elements=[DiscBilliardElement('particle', x=x, y=y, vx=vx, vy=vy, radius=radius)],
        interactions=[wall],
    )

coupled_pendulums(theta1=0.2, theta2=0.0, *, k=0.5, m1=1.0, m2=1.0, l1=1.0, l2=1.0, g=9.81)

Return two pendulums whose bobs are joined by a spring, which pass their energy back and forth.

The spring stretches with the horizontal distance l2 sin θ2 - l1 sin θ1 of the bobs and adds the torque k l1 (l2 sin θ2 - l1 sin θ1) cos θ1 to the first pendulum and the opposite one to the second. Identical pendulums have the normal modes sqrt(g / l) in phase and sqrt(g / l + 2 k / m) in antiphase.

Parameters:

Name Type Description Default
theta1 float

The initial angle of the first pendulum.

0.2
theta2 float

The initial angle of the second pendulum.

0.0
k float

The stiffness of the spring.

0.5
m1 float

The mass of the first bob.

1.0
m2 float

The mass of the second bob.

1.0
l1 float

The length of the first pendulum.

1.0
l2 float

The length of the second pendulum.

1.0
g float

The gravitational acceleration.

9.81

Returns:

Name Type Description
system System

The model, with the pendulums left and right and their angles left.theta and right.theta.

Examples:

>>> from skcomplex.models.physics.mechanics import coupled_pendulums
>>> coupled_pendulums(k=2.0).interactions[0].k
2.0
Source code in src/skcomplex/models/physics/mechanics/_pendulums.py
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode', 'second order', 'oscillation', 'normal modes'),
    references=('Ryzhov et al. 2021',),
    defaults={'until': 40.0, 'dt': 0.05},
)
def coupled_pendulums(
    theta1: float = 0.2,
    theta2: float = 0.0,
    *,
    k: float = 0.5,
    m1: float = 1.0,
    m2: float = 1.0,
    l1: float = 1.0,
    l2: float = 1.0,
    g: float = 9.81,
) -> System:
    """Return two pendulums whose bobs are joined by a spring, which pass their energy back and forth.

    The spring stretches with the horizontal distance `l2 sin θ2 - l1 sin θ1` of the bobs and adds the torque
    `k l1 (l2 sin θ2 - l1 sin θ1) cos θ1` to the first pendulum and the opposite one to the second. Identical pendulums
    have the normal modes `sqrt(g / l)` in phase and `sqrt(g / l + 2 k / m)` in antiphase.

    Args:
        theta1:
            The initial angle of the first pendulum.

        theta2:
            The initial angle of the second pendulum.

        k:
            The stiffness of the spring.

        m1:
            The mass of the first bob.

        m2:
            The mass of the second bob.

        l1:
            The length of the first pendulum.

        l2:
            The length of the second pendulum.

        g:
            The gravitational acceleration.

    Returns:
        system:
            The model, with the pendulums `left` and `right` and their angles `left.theta` and `right.theta`.

    Examples:
        >>> from skcomplex.models.physics.mechanics import coupled_pendulums
        >>> coupled_pendulums(k=2.0).interactions[0].k
        2.0
    """
    return System(
        elements=[
            PendulumElement('left', theta=theta1, m=m1, l=l1, g=g),
            PendulumElement('right', theta=theta2, m=m2, l=l2, g=g),
        ],
        interactions=[
            PendulumSpringInteraction('spring', first='left', second='right', k=k),
            PendulumInteraction('left_swing', target='left'),
            PendulumInteraction('right_swing', target='right'),
        ],
    )

double_pendulum(theta1=math.pi / 2, theta2=math.pi / 2, omega1=0.0, omega2=0.0, *, m1=1.0, m2=1.0, l1=1.0, l2=1.0, g=9.81)

Return a double pendulum, a pendulum hanging from the bob of another, written as its two Lagrange equations.

The equations (m1 + m2) l1 θ1'' + m2 l2 θ2'' cos(θ1 - θ2) + m2 l2 θ2'² sin(θ1 - θ2) + (m1 + m2) g sin θ1 = 0 and l2 θ2'' + l1 θ1'' cos(θ1 - θ2) - l1 θ1'² sin(θ1 - θ2) + g sin θ2 = 0 are solved jointly for the angular accelerations. The motion is chaotic at large amplitudes.

Parameters:

Name Type Description Default
theta1 float

The initial angle of the upper rod from the vertical, in radians.

pi / 2
theta2 float

The initial angle of the lower rod from the vertical, in radians.

pi / 2
omega1 float

The initial angular velocity of the upper rod.

0.0
omega2 float

The initial angular velocity of the lower rod.

0.0
m1 float

The mass of the upper bob.

1.0
m2 float

The mass of the lower bob.

1.0
l1 float

The length of the upper rod.

1.0
l2 float

The length of the lower rod.

1.0
g float

The gravitational acceleration.

9.81

Returns:

Name Type Description
system System

The model, with the angles pendulum.theta1 and pendulum.theta2 of order 2.

Examples:

>>> from skcomplex.models.physics.mechanics import double_pendulum
>>> double_pendulum(m2=2.0).elements[0].m2
2.0
Source code in src/skcomplex/models/physics/mechanics/_pendulums.py
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode', 'second order', 'chaos', 'lagrangian'),
    references=('Ryzhov et al. 2021', 'Argun et al. 2021'),
    defaults={'until': 20.0, 'dt': 0.02},
)
def double_pendulum(
    theta1: float = math.pi / 2,
    theta2: float = math.pi / 2,
    omega1: float = 0.0,
    omega2: float = 0.0,
    *,
    m1: float = 1.0,
    m2: float = 1.0,
    l1: float = 1.0,
    l2: float = 1.0,
    g: float = 9.81,
) -> System:
    """Return a double pendulum, a pendulum hanging from the bob of another, written as its two Lagrange equations.

    The equations `(m1 + m2) l1 θ1'' + m2 l2 θ2'' cos(θ1 - θ2) + m2 l2 θ2'² sin(θ1 - θ2) + (m1 + m2) g sin θ1 = 0`
    and `l2 θ2'' + l1 θ1'' cos(θ1 - θ2) - l1 θ1'² sin(θ1 - θ2) + g sin θ2 = 0` are solved jointly for the angular
    accelerations. The motion is chaotic at large amplitudes.

    Args:
        theta1:
            The initial angle of the upper rod from the vertical, in radians.

        theta2:
            The initial angle of the lower rod from the vertical, in radians.

        omega1:
            The initial angular velocity of the upper rod.

        omega2:
            The initial angular velocity of the lower rod.

        m1:
            The mass of the upper bob.

        m2:
            The mass of the lower bob.

        l1:
            The length of the upper rod.

        l2:
            The length of the lower rod.

        g:
            The gravitational acceleration.

    Returns:
        system:
            The model, with the angles `pendulum.theta1` and `pendulum.theta2` of order 2.

    Examples:
        >>> from skcomplex.models.physics.mechanics import double_pendulum
        >>> double_pendulum(m2=2.0).elements[0].m2
        2.0
    """
    return System(
        elements=[
            DoublePendulumElement(
                'pendulum',
                theta1=theta1,
                theta2=theta2,
                theta1_d1=omega1,
                theta2_d1=omega2,
                m1=m1,
                m2=m2,
                l1=l1,
                l2=l2,
                g=g,
            ),
        ],
        interactions=[DoublePendulumInteraction('motion', target='pendulum')],
    )

driven_damped_oscillator(x=1.0, v=0.0, *, m=1.0, k=1.0, c=0.1, F=0.0, omega=1.0)

Return a mass on a spring with damping and a periodic force, m x'' = -k x - c x' + F cos(ω t).

Parameters:

Name Type Description Default
x float

The initial displacement.

1.0
v float

The initial velocity.

0.0
m float

The mass.

1.0
k float

The stiffness.

1.0
c float

The damping coefficient.

0.1
F float

The amplitude of the driving force.

0.0
omega float

The angular frequency of the driving force.

1.0

Returns:

Name Type Description
system System

The model, with the mass mass in one dimension, whose displacement mass.x has order 2.

Examples:

>>> from skcomplex.models.physics.mechanics import driven_damped_oscillator
>>> driven_damped_oscillator(F=0.5).interactions[2].F
0.5
Source code in src/skcomplex/models/physics/mechanics/_mechanics.py
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode', 'second order', 'oscillation', 'forcing'),
    references=('Hooke 1678',),
    defaults={'until': 50.0, 'dt': 0.05},
)
def driven_damped_oscillator(
    x: float = 1.0,
    v: float = 0.0,
    *,
    m: float = 1.0,
    k: float = 1.0,
    c: float = 0.1,
    F: float = 0.0,
    omega: float = 1.0,
) -> System:
    """Return a mass on a spring with damping and a periodic force, `m x'' = -k x - c x' + F cos(ω t)`.

    Args:
        x:
            The initial displacement.

        v:
            The initial velocity.

        m:
            The mass.

        k:
            The stiffness.

        c:
            The damping coefficient.

        F:
            The amplitude of the driving force.

        omega:
            The angular frequency of the driving force.

    Returns:
        system:
            The model, with the mass `mass` in one dimension, whose displacement `mass.x` has order 2.

    Examples:
        >>> from skcomplex.models.physics.mechanics import driven_damped_oscillator
        >>> driven_damped_oscillator(F=0.5).interactions[2].F
        0.5
    """
    return System(
        elements=[MassPointElement('mass', x=x, x_d1=v, m=m)],
        interactions=[
            HarmonicPotentialInteraction('spring', target='mass', k=k),
            DampingInteraction('damping', target='mass', gamma=c),
            HarmonicForceInteraction('drive', target='mass', F=F, omega=omega),
        ],
        space=ContinuousSpace(n_dim=1),
    )

gimballed_gyroscope(spin=50.0, beta=0.3, *, A=1.0, C=2.0, A1=0.5, B1=0.5, C1=0.5, A2=1.0, torque=1.0)

Return a gyroscope in a gimbal: a rotor in an inner ring that turns in an outer ring about a fixed axis.

The outer ring turns by alpha about the fixed axis, the inner ring by beta in the outer ring, and the rotor spins by gamma in the inner ring. The equations are the Lagrange equations of the kinetic energy [A2 + (A + A1) cos² beta + (C + C1) sin² beta] alpha'² / 2 + (A + B1) beta'² / 2 + C alpha' gamma' sin beta + C gamma'² / 2, with a torque on the axis of the inner ring, which makes the fast rotor precess about the fixed axis. The momenta conjugate to alpha and gamma are conserved.

Parameters:

Name Type Description Default
spin float

The initial angular velocity of the rotor in the inner ring.

50.0
beta float

The initial angle of the inner ring.

0.3
A float

The equatorial moment of inertia of the rotor.

1.0
C float

The polar moment of inertia of the rotor.

2.0
A1 float

The first principal moment of inertia of the inner ring.

0.5
B1 float

The second principal moment of inertia of the inner ring.

0.5
C1 float

The third principal moment of inertia of the inner ring.

0.5
A2 float

The moment of inertia of the outer ring about its axis.

1.0
torque float

The torque on the axis of the inner ring.

1.0

Returns:

Name Type Description
system System

The model, with the angles gyroscope.alpha, gyroscope.beta, and gyroscope.gamma, all of order 2.

Examples:

>>> from skcomplex.models.physics.mechanics import gimballed_gyroscope
>>> gimballed_gyroscope(spin=10.0).elements[0].gamma_d1
10.0
Source code in src/skcomplex/models/physics/mechanics/_rotations.py
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode', 'second order', 'rotation', 'gyroscope', 'lagrangian'),
    references=('Ryzhov et al. 2021',),
    defaults={'until': 10.0, 'dt': 0.01},
)
def gimballed_gyroscope(
    spin: float = 50.0,
    beta: float = 0.3,
    *,
    A: float = 1.0,
    C: float = 2.0,
    A1: float = 0.5,
    B1: float = 0.5,
    C1: float = 0.5,
    A2: float = 1.0,
    torque: float = 1.0,
) -> System:
    """Return a gyroscope in a gimbal: a rotor in an inner ring that turns in an outer ring about a fixed axis.

    The outer ring turns by `alpha` about the fixed axis, the inner ring by `beta` in the outer ring, and the rotor
    spins by `gamma` in the inner ring. The equations are the Lagrange equations of the kinetic energy
    `[A2 + (A + A1) cos² beta + (C + C1) sin² beta] alpha'² / 2 + (A + B1) beta'² / 2 + C alpha' gamma' sin beta
    + C gamma'² / 2`, with a torque on the axis of the inner ring, which makes the fast rotor precess about the fixed
    axis. The momenta conjugate to `alpha` and `gamma` are conserved.

    Args:
        spin:
            The initial angular velocity of the rotor in the inner ring.

        beta:
            The initial angle of the inner ring.

        A:
            The equatorial moment of inertia of the rotor.

        C:
            The polar moment of inertia of the rotor.

        A1:
            The first principal moment of inertia of the inner ring.

        B1:
            The second principal moment of inertia of the inner ring.

        C1:
            The third principal moment of inertia of the inner ring.

        A2:
            The moment of inertia of the outer ring about its axis.

        torque:
            The torque on the axis of the inner ring.

    Returns:
        system:
            The model, with the angles `gyroscope.alpha`, `gyroscope.beta`, and `gyroscope.gamma`, all of order 2.

    Examples:
        >>> from skcomplex.models.physics.mechanics import gimballed_gyroscope
        >>> gimballed_gyroscope(spin=10.0).elements[0].gamma_d1
        10.0
    """
    return System(
        elements=[
            GyroscopeElement(
                'gyroscope',
                beta=beta,
                gamma_d1=spin,
                A=A,
                C=C,
                A1=A1,
                B1=B1,
                C1=C1,
                A2=A2,
                torque=torque,
            ),
        ],
        interactions=[GyroscopeInteraction('motion', target='gyroscope')],
    )

harmonic_oscillator(x=1.0, v=0.0, *, k=1.0, mass=1.0, damping=0.0)

Return a particle on a spring, optionally damped.

Parameters:

Name Type Description Default
x float

The initial displacement.

1.0
v float

The initial velocity.

0.0
k float

The stiffness.

1.0
mass float

The mass.

1.0
damping float

The linear damping coefficient.

0.0

Returns:

Name Type Description
system System

The model, with the particle mass in one dimension.

Examples:

>>> from skcomplex.models.physics.mechanics import harmonic_oscillator
>>> harmonic_oscillator(k=4.0).interactions[0].k
4.0
Source code in src/skcomplex/models/physics/mechanics/_mechanics.py
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode', 'oscillation'),
    references=('Hooke 1678',),
    defaults={'until': 20.0, 'dt': 0.05},
)
def harmonic_oscillator(
    x: float = 1.0,
    v: float = 0.0,
    *,
    k: float = 1.0,
    mass: float = 1.0,
    damping: float = 0.0,
) -> System:
    """Return a particle on a spring, optionally damped.

    Args:
        x:
            The initial displacement.

        v:
            The initial velocity.

        k:
            The stiffness.

        mass:
            The mass.

        damping:
            The linear damping coefficient.

    Returns:
        system:
            The model, with the particle `mass` in one dimension.

    Examples:
        >>> from skcomplex.models.physics.mechanics import harmonic_oscillator
        >>> harmonic_oscillator(k=4.0).interactions[0].k
        4.0
    """
    interactions: list[BaseInteraction] = [HarmonicPotentialInteraction('spring', target='mass', k=k)]
    if damping:
        interactions.append(DampingInteraction('damping', target='mass', gamma=damping))
    return System(
        elements=[MassPointElement('mass', x=x, x_d1=v, m=mass)],
        interactions=interactions,
        space=ContinuousSpace(n_dim=1),
    )

huygens_clocks(n_clocks=2, *, theta=0.2, k=1.0, m=1.0, l=1.0, M=1.0, g=9.8, rho=0.1, lam=0.1, theta_c=0.1)

Return pendulum clocks on supports joined by springs, which synchronize through the motion of the supports.

Every pendulum hangs from a support of mass M that moves horizontally with damping rho, and neighbouring supports are joined by springs of stiffness k. A pendulum follows m l θ'' + f(θ) θ' + m g sin θ = -m x'' cos θ and its support (M + m) x'' + 2 rho x' + m l θ'' cos θ = m l θ'² sin θ + spring forces. The escapement f(θ) = 2 λ sgn(|θ| - θ_c) is -2 λ for |θ| < θ_c, pumping energy in, and 2 λ beyond, so every clock has a limit cycle. With λ = 0 the clocks are plain pendulums. Two identical clocks synchronize in antiphase. The first clock starts at theta, away from the switch of the escapement at θ_c, and the others at rest.

Parameters:

Name Type Description Default
n_clocks int

The number of clocks in the row.

2
theta float

The initial angle of the first clock. The others start at rest.

0.2
k float

The stiffness of the springs between neighbouring supports.

1.0
m float

The mass of every bob.

1.0
l float

The length of every pendulum.

1.0
M float

The mass of every support.

1.0
g float

The gravitational acceleration.

9.8
rho float

The damping of the supports.

0.1
lam float

The strength of the escapement.

0.1
theta_c float

The angle below which the escapement pumps energy in.

0.1

Returns:

Name Type Description
system System

The model, with the angles clock_<i>.theta and the support positions clock_<i>.x, all of order 2.

Examples:

>>> from skcomplex.models.physics.mechanics import huygens_clocks
>>> len(huygens_clocks(3).interactions)
5
Source code in src/skcomplex/models/physics/mechanics/_pendulums.py
460
461
462
463
464
465
466
467
468
469
470
471
472
473
474
475
476
477
478
479
480
481
482
483
484
485
486
487
488
489
490
491
492
493
494
495
496
497
498
499
500
501
502
503
504
505
506
507
508
509
510
511
512
513
514
515
516
517
518
519
520
521
522
523
524
525
526
527
528
529
530
531
532
533
534
535
536
537
538
539
540
541
542
543
544
545
546
547
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode', 'second order', 'synchronization', 'limit cycle'),
    references=('Huygens 1673', 'Dilão 2023', 'Dilão 2023'),
    defaults={'until': 20.0, 'dt': 0.05},
)
def huygens_clocks(
    n_clocks: int = 2,
    *,
    theta: float = 0.2,
    k: float = 1.0,
    m: float = 1.0,
    l: float = 1.0,
    M: float = 1.0,
    g: float = 9.8,
    rho: float = 0.1,
    lam: float = 0.1,
    theta_c: float = 0.1,
) -> System:
    """Return pendulum clocks on supports joined by springs, which synchronize through the motion of the supports.

    Every pendulum hangs from a support of mass `M` that moves horizontally with damping `rho`, and neighbouring
    supports are joined by springs of stiffness `k`. A pendulum follows `m l θ'' + f(θ) θ' + m g sin θ = -m x'' cos θ`
    and its support `(M + m) x'' + 2 rho x' + m l θ'' cos θ = m l θ'² sin θ + spring forces`. The escapement
    `f(θ) = 2 λ sgn(|θ| - θ_c)` is `-2 λ` for `|θ| < θ_c`, pumping energy in, and `2 λ` beyond, so every clock has a
    limit cycle. With `λ = 0` the clocks are plain pendulums. Two identical clocks synchronize in antiphase. The
    first clock starts at `theta`, away from the switch of the escapement at `θ_c`, and the others at rest.

    Args:
        n_clocks:
            The number of clocks in the row.

        theta:
            The initial angle of the first clock. The others start at rest.

        k:
            The stiffness of the springs between neighbouring supports.

        m:
            The mass of every bob.

        l:
            The length of every pendulum.

        M:
            The mass of every support.

        g:
            The gravitational acceleration.

        rho:
            The damping of the supports.

        lam:
            The strength of the escapement.

        theta_c:
            The angle below which the escapement pumps energy in.

    Returns:
        system:
            The model, with the angles `clock_<i>.theta` and the support positions `clock_<i>.x`, all of order 2.

    Examples:
        >>> from skcomplex.models.physics.mechanics import huygens_clocks
        >>> len(huygens_clocks(3).interactions)
        5
    """
    clocks = [
        ClockPendulumElement(
            f'clock_{index}',
            theta=theta if index == 0 else 0.0,
            m=m,
            l=l,
            M=M,
            g=g,
            rho=rho,
            lam=lam,
            theta_c=theta_c,
        )
        for index in range(n_clocks)
    ]
    springs = [
        SupportSpringInteraction(f'spring_{index}', first=f'clock_{index}', second=f'clock_{index + 1}', k=k)
        for index in range(n_clocks - 1)
    ]
    motions = [ClockPendulumInteraction(f'motion_{index}', target=f'clock_{index}') for index in range(n_clocks)]
    return System(elements=clocks, interactions=[*springs, *motions])

magnetic_pendulum(x=1.0, y=1.5, *, n_magnets=3, omega=1.0, q=0.2, gamma=1.0, d=0.3)

Return a magnetic bob on a string swinging above magnets on a circle, which comes to rest at one of them.

The bob follows x'' = -ω² x - q x' - Σ gamma (x - x_i) / d_i³ with d_i = sqrt(|x - x_i|² + d²). The basins of the magnets have fractal boundaries.

Parameters:

Name Type Description Default
x float

The initial horizontal position of the bob.

1.0
y float

The initial vertical position of the bob.

1.5
n_magnets int

The number of magnets, evenly spaced on the unit circle.

3
omega float

The natural frequency of the pendulum.

1.0
q float

The friction coefficient.

0.2
gamma float

The strength of every magnet.

1.0
d float

The height of the bob above the plane of the magnets.

0.3

Returns:

Name Type Description
system System

The model, with the bob bob and its position bob.x in two dimensions.

Examples:

>>> from skcomplex.models.physics.mechanics import magnetic_pendulum
>>> len(magnetic_pendulum(n_magnets=4).interactions)
6
Source code in src/skcomplex/models/physics/mechanics/_pendulums.py
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
370
371
372
373
374
375
376
377
378
379
380
381
382
383
384
385
386
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode', 'multistability', 'fractal basins'),
    references=('Datseris and Parlitz 2022',),
    defaults={'until': 50.0, 'dt': 0.05},
)
def magnetic_pendulum(
    x: float = 1.0,
    y: float = 1.5,
    *,
    n_magnets: int = 3,
    omega: float = 1.0,
    q: float = 0.2,
    gamma: float = 1.0,
    d: float = 0.3,
) -> System:
    """Return a magnetic bob on a string swinging above magnets on a circle, which comes to rest at one of them.

    The bob follows `x'' = -ω² x - q x' - Σ gamma (x - x_i) / d_i³` with `d_i = sqrt(|x - x_i|² + d²)`. The basins of
    the magnets have fractal boundaries.

    Args:
        x:
            The initial horizontal position of the bob.

        y:
            The initial vertical position of the bob.

        n_magnets:
            The number of magnets, evenly spaced on the unit circle.

        omega:
            The natural frequency of the pendulum.

        q:
            The friction coefficient.

        gamma:
            The strength of every magnet.

        d:
            The height of the bob above the plane of the magnets.

    Returns:
        system:
            The model, with the bob `bob` and its position `bob.x` in two dimensions.

    Examples:
        >>> from skcomplex.models.physics.mechanics import magnetic_pendulum
        >>> len(magnetic_pendulum(n_magnets=4).interactions)
        6
    """
    angles = 2 * np.pi * np.arange(n_magnets) / n_magnets
    magnets = [
        MagnetInteraction(
            f'magnet_{index}',
            target='bob',
            center=[round(float(np.cos(angle)), 12), round(float(np.sin(angle)), 12)],
            gamma=gamma,
            d=d,
        )
        for index, angle in enumerate(angles)
    ]
    return System(
        elements=[MassPointElement('bob', x=[x, y], x_d1=[0.0, 0.0])],
        interactions=[
            HarmonicPotentialInteraction('string', target='bob', k=omega**2),
            DampingInteraction('friction', target='bob', gamma=q),
            *magnets,
        ],
        space=ContinuousSpace(n_dim=2),
    )

newtons_cradle(n_balls=5, *, speed=0.5, diameter=0.1, gap=0.001, length=1.0, g=9.81)

Return Newton's cradle: pendulum balls in a row that pass momentum along by elastic collisions.

The balls swing with small amplitude about rest positions that leave a small gap between neighbours, so collisions are separate events. The first ball starts with a speed towards the others, and equal masses exchange velocities at every collision.

Parameters:

Name Type Description Default
n_balls int

The number of balls.

5
speed float

The initial speed of the first ball.

0.5
diameter float

The diameter of the balls.

0.1
gap float

The gap between neighbouring balls at rest.

0.001
length float

The length of the pendulums.

1.0
g float

The gravitational acceleration.

9.81

Returns:

Name Type Description
system System

The model, with the horizontal positions ball_<i>.x and velocities ball_<i>.x_d1.

Examples:

>>> from skcomplex.models.physics.mechanics import newtons_cradle
>>> len(newtons_cradle(n_balls=3).interactions)
5
Source code in src/skcomplex/models/physics/mechanics/_mechanics.py
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
329
330
331
332
333
334
335
336
337
338
339
340
341
342
343
344
345
346
347
348
349
350
351
352
353
354
355
356
357
358
359
360
361
362
363
364
365
366
367
368
369
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'events', 'impacts', 'momentum'),
    references=('Hinch and Saint-Jean 1999',),
    defaults={'until': 5.0, 'dt': 0.01},
)
def newtons_cradle(
    n_balls: int = 5,
    *,
    speed: float = 0.5,
    diameter: float = 0.1,
    gap: float = 0.001,
    length: float = 1.0,
    g: float = 9.81,
) -> System:
    """Return Newton's cradle: pendulum balls in a row that pass momentum along by elastic collisions.

    The balls swing with small amplitude about rest positions that leave a small gap between neighbours, so collisions
    are separate events. The first ball starts with a speed towards the others, and equal masses exchange velocities
    at every collision.

    Args:
        n_balls:
            The number of balls.

        speed:
            The initial speed of the first ball.

        diameter:
            The diameter of the balls.

        gap:
            The gap between neighbouring balls at rest.

        length:
            The length of the pendulums.

        g:
            The gravitational acceleration.

    Returns:
        system:
            The model, with the horizontal positions `ball_<i>.x` and velocities `ball_<i>.x_d1`.

    Examples:
        >>> from skcomplex.models.physics.mechanics import newtons_cradle
        >>> len(newtons_cradle(n_balls=3).interactions)
        5
    """
    spacing = diameter + gap
    elements = [
        HangingBallElement(
            f'ball_{index}',
            x=index * spacing,
            x_d1=speed if index == 0 else 0.0,
            rest=index * spacing,
            g=g,
            length=length,
        )
        for index in range(n_balls)
    ]
    swings: list[BaseInteraction] = [
        GradientInteraction(
            f'swing_{index}',
            target=f'ball_{index}',
            potential=lambda s: s.target.g * (s.target.x - s.target.rest) ** 2 / (2 * s.target.length),
        )
        for index in range(n_balls)
    ]
    collisions: list[BaseInteraction] = [
        ElasticCollisionInteraction(
            f'collision_{index}',
            left=f'ball_{index}',
            right=f'ball_{index + 1}',
            diameter=diameter,
        )
        for index in range(n_balls - 1)
    ]
    return System(elements=elements, interactions=[*swings, *collisions])

parametric_pendulum(theta=0.1, omega=0.0, *, g=9.81, l=1.0, epsilon=0.2, Omega=2 * math.sqrt(9.81), c=0.0)

Return a pendulum whose suspension point oscillates vertically, θ'' + c θ' + (g / l)(1 + ε cos Ω t) sin θ = 0.

The vertical acceleration of the suspension point modulates gravity. Small oscillations grow when Ω is close to twice the natural frequency sqrt(g / l), the parametric resonance of a swing.

Parameters:

Name Type Description Default
theta float

The initial angle.

0.1
omega float

The initial angular velocity.

0.0
g float

The gravitational acceleration.

9.81
l float

The length of the rod.

1.0
epsilon float

The amplitude of the acceleration of the suspension point, relative to g.

0.2
Omega float

The angular frequency of the suspension point.

2 * sqrt(9.81)
c float

The friction coefficient.

0.0

Returns:

Name Type Description
system System

The model, with the angle swing.theta of order 2.

Examples:

>>> from skcomplex.models.physics.mechanics import parametric_pendulum
>>> parametric_pendulum(epsilon=0.5).elements[0].epsilon
0.5
Source code in src/skcomplex/models/physics/mechanics/_pendulums.py
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode', 'second order', 'forcing', 'parametric resonance'),
    references=('Dilão 2023',),
    defaults={'until': 30.0, 'dt': 0.02},
)
def parametric_pendulum(
    theta: float = 0.1,
    omega: float = 0.0,
    *,
    g: float = 9.81,
    l: float = 1.0,
    epsilon: float = 0.2,
    Omega: float = 2 * math.sqrt(9.81),
    c: float = 0.0,
) -> System:
    """Return a pendulum whose suspension point oscillates vertically, `θ'' + c θ' + (g / l)(1 + ε cos Ω t) sin θ = 0`.

    The vertical acceleration of the suspension point modulates gravity. Small oscillations grow when `Ω` is close to
    twice the natural frequency `sqrt(g / l)`, the parametric resonance of a swing.

    Args:
        theta:
            The initial angle.

        omega:
            The initial angular velocity.

        g:
            The gravitational acceleration.

        l:
            The length of the rod.

        epsilon:
            The amplitude of the acceleration of the suspension point, relative to `g`.

        Omega:
            The angular frequency of the suspension point.

        c:
            The friction coefficient.

    Returns:
        system:
            The model, with the angle `swing.theta` of order 2.

    Examples:
        >>> from skcomplex.models.physics.mechanics import parametric_pendulum
        >>> parametric_pendulum(epsilon=0.5).elements[0].epsilon
        0.5
    """
    return System(
        elements=[
            ParametricPendulumElement(
                'swing',
                theta=theta,
                theta_d1=omega,
                g=g,
                l=l,
                epsilon=epsilon,
                Omega=Omega,
                c=c,
            ),
        ],
        interactions=[ParametricPendulumInteraction('motion', target='swing')],
    )

particle_in_a_well(x=1.0, v=0.5, *, k=2.0, m=1.0)

Return a particle moving without friction in the potential U(x) = k x² e^-x.

The potential has a stable equilibrium at x = 0 and an unstable one at x = 2, the top of a barrier of height 4 k e^-2. A particle in the well with less energy than the barrier oscillates, and one with more escapes to large x.

Parameters:

Name Type Description Default
x float

The initial position.

1.0
v float

The initial velocity.

0.5
k float

The strength of the potential.

2.0
m float

The mass of the particle.

1.0

Returns:

Name Type Description
system System

The model, with the particle particle and its position particle.x in one dimension.

Examples:

>>> from skcomplex.models.physics.mechanics import particle_in_a_well
>>> particle_in_a_well(k=3.0).interactions[0].k
3.0
Source code in src/skcomplex/models/physics/mechanics/_wells.py
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode', 'potential', 'hamiltonian'),
    references=('Starkovich 2021',),
    defaults={'until': 20.0, 'dt': 0.05},
)
def particle_in_a_well(x: float = 1.0, v: float = 0.5, *, k: float = 2.0, m: float = 1.0) -> System:
    """Return a particle moving without friction in the potential `U(x) = k x² e^-x`.

    The potential has a stable equilibrium at `x = 0` and an unstable one at `x = 2`, the top of a barrier of height
    `4 k e^-2`. A particle in the well with less energy than the barrier oscillates, and one with more escapes to
    large `x`.

    Args:
        x:
            The initial position.

        v:
            The initial velocity.

        k:
            The strength of the potential.

        m:
            The mass of the particle.

    Returns:
        system:
            The model, with the particle `particle` and its position `particle.x` in one dimension.

    Examples:
        >>> from skcomplex.models.physics.mechanics import particle_in_a_well
        >>> particle_in_a_well(k=3.0).interactions[0].k
        3.0
    """
    return System(
        elements=[MassPointElement('particle', x=x, x_d1=v, m=m)],
        interactions=[ExponentialWellInteraction('well', target='particle', k=k)],
        space=ContinuousSpace(n_dim=1),
    )

pendulum(theta=0.5, omega=0.0, *, m=1.0, g=9.81, l=1.0)

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

Parameters:

Name Type Description Default
theta float

The initial angle, in radians.

0.5
omega float

The initial angular velocity.

0.0
m float

The mass of the bob.

1.0
g float

The gravitational acceleration.

9.81
l float

The length of the rod.

1.0

Returns:

Name Type Description
system 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
Source code in src/skcomplex/models/physics/mechanics/_mechanics.py
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
@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')],
    )

projectile(speed=10.0, angle=45.0, *, mass=1.0, drag=0.0)

Return a particle launched under uniform gravity, optionally slowed by linear drag.

Parameters:

Name Type Description Default
speed float

The launch speed.

10.0
angle float

The launch angle above the horizon, in degrees.

45.0
mass float

The mass.

1.0
drag float

The linear drag coefficient.

0.0

Returns:

Name Type Description
system System

The model, with the particle ball in two dimensions.

Examples:

>>> from skcomplex.models.physics.mechanics import projectile
>>> projectile(speed=5.0, angle=90.0).elements[0].x_d1[1]
5.0
Source code in src/skcomplex/models/physics/mechanics/_mechanics.py
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode'),
    references=('Newton 1687',),
    defaults={'until': 2.0, 'dt': 0.01},
)
def projectile(speed: float = 10.0, angle: float = 45.0, *, mass: float = 1.0, drag: float = 0.0) -> System:
    """Return a particle launched under uniform gravity, optionally slowed by linear drag.

    Args:
        speed:
            The launch speed.

        angle:
            The launch angle above the horizon, in degrees.

        mass:
            The mass.

        drag:
            The linear drag coefficient.

    Returns:
        system:
            The model, with the particle `ball` in two dimensions.

    Examples:
        >>> from skcomplex.models.physics.mechanics import projectile
        >>> projectile(speed=5.0, angle=90.0).elements[0].x_d1[1]
        5.0
    """
    radians = np.radians(angle)
    velocity = [round(float(speed * np.cos(radians)), 12), round(float(speed * np.sin(radians)), 12)]
    interactions: list[BaseInteraction] = [
        ConstantForceInteraction('gravity', target='ball', F=[0.0, -GRAVITY * mass]),
    ]
    if drag:
        interactions.append(DampingInteraction('drag', target='ball', gamma=drag))
    return System(
        elements=[MassPointElement('ball', x=[0.0, 0.0], x_d1=velocity, m=mass)],
        interactions=interactions,
        space=ContinuousSpace(n_dim=2),
    )

rotating_balls_on_rods(omega=1.0, *, r=1.0, beta=math.pi / 2, m=1.0)

Return two balls on rods that turn about a vertical axis and conserve their angular momentum.

Each ball sits at the distance r along a rod fixed at the angle β to the axis, so its distance from the axis is r sin β. Its angular momentum H = m ω r² sin² β is fixed by the initial angular velocity, and the balls turn at θ' = H / (m r² sin² β), so changing r during a run changes the angular velocity by 1 / r².

Parameters:

Name Type Description Default
omega float

The initial angular velocity.

1.0
r float

The distance of the balls from the joint of the rods.

1.0
beta float

The angle of the rods to the vertical axis, in radians.

pi / 2
m float

The mass of every ball.

1.0

Returns:

Name Type Description
system System

The model, with the angles right.theta and left.theta of the two balls.

Examples:

>>> from skcomplex.models.physics.mechanics import rotating_balls_on_rods
>>> rotating_balls_on_rods(omega=2.0).elements[0].H
2.0
Source code in src/skcomplex/models/physics/mechanics/_rotations.py
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode', 'rotation', 'angular momentum'),
    references=('Ryzhov et al. 2021',),
    defaults={'until': 10.0, 'dt': 0.02},
)
def rotating_balls_on_rods(
    omega: float = 1.0,
    *,
    r: float = 1.0,
    beta: float = math.pi / 2,
    m: float = 1.0,
) -> System:
    """Return two balls on rods that turn about a vertical axis and conserve their angular momentum.

    Each ball sits at the distance `r` along a rod fixed at the angle `β` to the axis, so its distance from the axis is
    `r sin β`. Its angular momentum `H = m ω r² sin² β` is fixed by the initial angular velocity, and the balls turn at
    `θ' = H / (m r² sin² β)`, so changing `r` during a run changes the angular velocity by `1 / r²`.

    Args:
        omega:
            The initial angular velocity.

        r:
            The distance of the balls from the joint of the rods.

        beta:
            The angle of the rods to the vertical axis, in radians.

        m:
            The mass of every ball.

    Returns:
        system:
            The model, with the angles `right.theta` and `left.theta` of the two balls.

    Examples:
        >>> from skcomplex.models.physics.mechanics import rotating_balls_on_rods
        >>> rotating_balls_on_rods(omega=2.0).elements[0].H
        2.0
    """
    momentum = m * omega * r**2 * math.sin(beta) ** 2
    return System(
        elements=[
            RotatingBallElement('right', theta=0.0, m=m, r=r, beta=beta, H=momentum),
            RotatingBallElement('left', theta=math.pi, m=m, r=r, beta=beta, H=momentum),
        ],
        interactions=[
            RateInteraction(
                f'{ball}_turn',
                target=ball,
                state='theta',
                rate=lambda s: s.target.H / (s.target.m * s.target.r**2 * sin(s.target.beta) ** 2),
            )
            for ball in ('right', 'left')
        ],
    )

rotating_swing(omega=2.0, theta=0.1, *, n_chairs=8, J=10.0, m=1.0, a=1.0, l=1.0, g=9.81, torque=0.0, friction=0.0)

Return a rotating swing: chairs on chains hanging from the rods of a stand that turns about a vertical axis.

The stand of inertia J turns by the angle φ, and the chairs, all deflected by the same angle θ, follow the Lagrange equations [J + n m (a + l sin θ)²] φ'' + 2 n m l φ' θ' (a + l sin θ) cos θ = M - c φ' and l θ'' - φ'² (a + l sin θ) cos θ + g sin θ = 0. At a steady speed ω the chairs settle at the angle where g tan θ = ω² (a + l sin θ). Without torque and friction the angular momentum of the swing is conserved.

Parameters:

Name Type Description Default
omega float

The initial angular velocity of the stand.

2.0
theta float

The initial angle of the chains from the vertical.

0.1
n_chairs int

The number of chairs.

8
J float

The moment of inertia of the stand.

10.0
m float

The mass of every chair.

1.0
a float

The length of the rods.

1.0
l float

The length of the chains.

1.0
g float

The gravitational acceleration.

9.81
torque float

The torque that drives the stand.

0.0
friction float

The friction coefficient of the stand's axis.

0.0

Returns:

Name Type Description
system System

The model, with the angle of the stand swing.phi and of the chains swing.theta, both of order 2.

Examples:

>>> from skcomplex.models.physics.mechanics import rotating_swing
>>> rotating_swing(omega=1.0).elements[0].phi_d1
1.0
Source code in src/skcomplex/models/physics/mechanics/_rotations.py
 20
 21
 22
 23
 24
 25
 26
 27
 28
 29
 30
 31
 32
 33
 34
 35
 36
 37
 38
 39
 40
 41
 42
 43
 44
 45
 46
 47
 48
 49
 50
 51
 52
 53
 54
 55
 56
 57
 58
 59
 60
 61
 62
 63
 64
 65
 66
 67
 68
 69
 70
 71
 72
 73
 74
 75
 76
 77
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode', 'second order', 'rotation', 'lagrangian'),
    references=('Ryzhov et al. 2021',),
    defaults={'until': 20.0, 'dt': 0.05},
)
def rotating_swing(
    omega: float = 2.0,
    theta: float = 0.1,
    *,
    n_chairs: int = 8,
    J: float = 10.0,
    m: float = 1.0,
    a: float = 1.0,
    l: float = 1.0,
    g: float = 9.81,
    torque: float = 0.0,
    friction: float = 0.0,
) -> System:
    """Return a rotating swing: chairs on chains hanging from the rods of a stand that turns about a vertical axis.

    The stand of inertia `J` turns by the angle `φ`, and the chairs, all deflected by the same angle `θ`, follow the
    Lagrange equations `[J + n m (a + l sin θ)²] φ'' + 2 n m l φ' θ' (a + l sin θ) cos θ = M - c φ'` and
    `l θ'' - φ'² (a + l sin θ) cos θ + g sin θ = 0`. At a steady speed `ω` the chairs settle at the angle where
    `g tan θ = ω² (a + l sin θ)`. Without torque and friction the angular momentum of the swing is conserved.

    Args:
        omega:
            The initial angular velocity of the stand.

        theta:
            The initial angle of the chains from the vertical.

        n_chairs:
            The number of chairs.

        J:
            The moment of inertia of the stand.

        m:
            The mass of every chair.

        a:
            The length of the rods.

        l:
            The length of the chains.

        g:
            The gravitational acceleration.

        torque:
            The torque that drives the stand.

        friction:
            The friction coefficient of the stand's axis.

    Returns:
        system:
            The model, with the angle of the stand `swing.phi` and of the chains `swing.theta`, both of order 2.

    Examples:
        >>> from skcomplex.models.physics.mechanics import rotating_swing
        >>> rotating_swing(omega=1.0).elements[0].phi_d1
        1.0
    """
    return System(
        elements=[
            RotatingSwingElement(
                'swing',
                phi_d1=omega,
                theta=theta,
                n_chairs=n_chairs,
                J=J,
                m=m,
                a=a,
                l=l,
                g=g,
                torque=torque,
                friction=friction,
            ),
        ],
        interactions=[RotatingSwingInteraction('motion', target='swing')],
    )

two_mass_spring_damper(x1=1.2, x2=2.0, *, m1=1.0, m2=1.0, k1=1.0, k2=1.0, k=1.0, b1=0.0, b2=0.0, l1=1.0, l2=1.0, F1=0.0, F2=0.0, omega1=1.0, omega2=1.0)

Return two masses between two walls, tied to the walls and to each other by springs, with dampers and forces.

The first mass rests at x01 = l1 and the second at x02 = l1 + l2. They follow m1 x1'' = -k1 (x1 - x01) + k ((x2 - x1) - (x02 - x01)) - b1 x1' + F1 cos(ω1 t) and m2 x2'' = -k ((x2 - x1) - (x02 - x01)) - k2 (x2 - x02) - b2 x2' + F2 cos(ω2 t).

Parameters:

Name Type Description Default
x1 float

The initial position of the first mass.

1.2
x2 float

The initial position of the second mass.

2.0
m1 float

The first mass.

1.0
m2 float

The second mass.

1.0
k1 float

The stiffness of the spring between the left wall and the first mass.

1.0
k2 float

The stiffness of the spring between the second mass and the right wall.

1.0
k float

The stiffness of the spring between the masses.

1.0
b1 float

The damping coefficient of the first mass.

0.0
b2 float

The damping coefficient of the second mass.

0.0
l1 float

The rest length of the first spring.

1.0
l2 float

The rest length of the middle spring.

1.0
F1 float

The amplitude of the periodic force on the first mass.

0.0
F2 float

The amplitude of the periodic force on the second mass.

0.0
omega1 float

The angular frequency of the force on the first mass.

1.0
omega2 float

The angular frequency of the force on the second mass.

1.0

Returns:

Name Type Description
system System

The model, with the masses first and second and their positions first.x and second.x.

Examples:

>>> from skcomplex.models.physics.mechanics import two_mass_spring_damper
>>> two_mass_spring_damper(l1=2.0).interactions[0].center
2.0
Source code in src/skcomplex/models/physics/mechanics/_springs.py
 78
 79
 80
 81
 82
 83
 84
 85
 86
 87
 88
 89
 90
 91
 92
 93
 94
 95
 96
 97
 98
 99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
@register_model(
    field='physics/mechanics',
    tags=('mechanics', 'ode', 'oscillation', 'normal modes', 'forcing'),
    references=('Ryzhov et al. 2021',),
    defaults={'until': 30.0, 'dt': 0.05},
)
def two_mass_spring_damper(
    x1: float = 1.2,
    x2: float = 2.0,
    *,
    m1: float = 1.0,
    m2: float = 1.0,
    k1: float = 1.0,
    k2: float = 1.0,
    k: float = 1.0,
    b1: float = 0.0,
    b2: float = 0.0,
    l1: float = 1.0,
    l2: float = 1.0,
    F1: float = 0.0,
    F2: float = 0.0,
    omega1: float = 1.0,
    omega2: float = 1.0,
) -> System:
    """Return two masses between two walls, tied to the walls and to each other by springs, with dampers and forces.

    The first mass rests at `x01 = l1` and the second at `x02 = l1 + l2`. They follow
    `m1 x1'' = -k1 (x1 - x01) + k ((x2 - x1) - (x02 - x01)) - b1 x1' + F1 cos(ω1 t)` and
    `m2 x2'' = -k ((x2 - x1) - (x02 - x01)) - k2 (x2 - x02) - b2 x2' + F2 cos(ω2 t)`.

    Args:
        x1:
            The initial position of the first mass.

        x2:
            The initial position of the second mass.

        m1:
            The first mass.

        m2:
            The second mass.

        k1:
            The stiffness of the spring between the left wall and the first mass.

        k2:
            The stiffness of the spring between the second mass and the right wall.

        k:
            The stiffness of the spring between the masses.

        b1:
            The damping coefficient of the first mass.

        b2:
            The damping coefficient of the second mass.

        l1:
            The rest length of the first spring.

        l2:
            The rest length of the middle spring.

        F1:
            The amplitude of the periodic force on the first mass.

        F2:
            The amplitude of the periodic force on the second mass.

        omega1:
            The angular frequency of the force on the first mass.

        omega2:
            The angular frequency of the force on the second mass.

    Returns:
        system:
            The model, with the masses `first` and `second` and their positions `first.x` and `second.x`.

    Examples:
        >>> from skcomplex.models.physics.mechanics import two_mass_spring_damper
        >>> two_mass_spring_damper(l1=2.0).interactions[0].center
        2.0
    """
    rest1, rest2 = l1, l1 + l2
    return System(
        elements=[MassPointElement('first', x=x1, x_d1=0.0, m=m1), MassPointElement('second', x=x2, x_d1=0.0, m=m2)],
        interactions=[
            HarmonicPotentialInteraction('left_spring', target='first', k=k1, center=rest1),
            HarmonicPotentialInteraction('right_spring', target='second', k=k2, center=rest2),
            SpringInteraction('middle_spring', first='first', second='second', k=k, length=l2),
            DampingInteraction('first_damper', target='first', gamma=b1),
            DampingInteraction('second_damper', target='second', gamma=b2),
            HarmonicForceInteraction('first_force', target='first', F=F1, omega=omega1),
            HarmonicForceInteraction('second_force', target='second', F=F2, omega=omega2),
        ],
        space=ContinuousSpace(n_dim=1),
    )