Interactive Plotly visualizers: SO(3) momentum sphere and orbit#

Unlike the static Matplotlib plots elsewhere in physicskit.classical.visualizers, these are draggable/zoomable in a live Python session or a browser – useful for directly inspecting how close a rigid-body trajectory passes to the unstable intermediate-axis pole, or for zooming into a precessing orbit’s perihelion. The two systems are governed by the equations already introduced elsewhere in the gallery: EulerTop’s torque-free Euler equations, \(I_i\dot\omega_i = (I_j - I_k)\,\omega_j\omega_k\) (cyclic in \(i,j,k\)) for the body-frame angular velocity (see The intermediate axis theorem (Dzhanibekov effect)), and KeplerSystem’s separable Hamiltonian \(H = |p|^2/(2\mu) + V(r)\) with the post-Newtonian perturbation \(c_\mathrm{pn}/r^3\) in \(V(r)\) that drives the perihelion precession (see Kepler orbits and perihelion precession).

from physicskit.classical.core.base_system import SimulationResult
from physicskit.classical.systems.newtonian import KeplerSystem
from physicskit.classical.systems.rotations import EulerTop
from physicskit.classical.visualizers.interactive import interactive_orbit, interactive_so3_momentum_sphere

# The integrations below use a fine dt for numerical accuracy, but an
# interactive plot of every single step (potentially hundreds of
# thousands of points) is both unnecessarily large and sluggish to pan
# and zoom -- downsample with a stride to a few hundred points, plenty
# to see the geometry clearly.


def _stride(result: SimulationResult, stride: int) -> SimulationResult:
    idx = list(range(0, len(result.t), stride))
    if idx[-1] != len(result.t) - 1:
        idx.append(len(result.t) - 1)
    q = result.q[idx] if result.q is not None else None
    p = result.p[idx] if result.p is not None else None
    return SimulationResult(t=result.t[idx], y=result.y[idx], q=q, p=p)

SO(3) momentum sphere for a tumbling (intermediate-axis) Euler top#

top = EulerTop([0.01, 1.0, 0.01], I1=1.0, I2=2.0, I3=3.0)
result_top = _stride(top.integrate((0, 100), dt=1e-1, method="implicit_midpoint"), stride=200)
fig_sphere = interactive_so3_momentum_sphere(result_top.y[:, :3], I1=1.0, I2=2.0, I3=3.0, title="EulerTop: tumbling about the intermediate axis")

A precessing Kepler orbit with its LRL vector rotating#

kepler = KeplerSystem.from_orbital_elements(a=1.0, e=0.3, c_pn=0.02)
result_orbit = _stride(kepler.integrate((0, 300), dt=1e-3, method="yoshida4"), stride=300)
fig_orbit = interactive_orbit(result_orbit, system=kepler, title="Precessing Kepler orbit (c_pn=0.02)")

Total running time of the script: (0 minutes 0.434 seconds)

Gallery generated by Sphinx-Gallery