Doppelpendelsimulation

Beginner

This tutorial is from open-source community. Access the source code

Einführung

Das Doppelpendel ist ein klassisches Problem in der Physik und Mathematik. Es besteht aus zwei miteinander verbundenen Pendeln, und die Bewegung des zweiten Pendels wird durch die Bewegung des ersten Pendels beeinflusst. In diesem Lab verwenden wir Python und Matplotlib, um die Bewegung eines Doppelpendels zu simulieren.

Tipps für die VM

Nachdem der VM-Start abgeschlossen ist, klicken Sie in der oberen linken Ecke, um zur Registerkarte Notebook zu wechseln und Jupyter Notebook für die Übung zu nutzen.

Manchmal müssen Sie einige Sekunden warten, bis Jupyter Notebook vollständig geladen ist. Die Validierung von Vorgängen kann aufgrund der Einschränkungen von Jupyter Notebook nicht automatisiert werden.

Wenn Sie bei der Lernphase Probleme haben, können Sie Labby gerne fragen. Geben Sie nach der Sitzung Feedback, und wir werden das Problem für Sie prompt beheben.

Importieren der erforderlichen Bibliotheken

Wir beginnen mit dem Importieren der erforderlichen Bibliotheken für unsere Simulation.

from collections import deque
import matplotlib.pyplot as plt
import numpy as np
from numpy import cos, sin
import matplotlib.animation as animation

Parameter festlegen

Als nächstes definieren wir die Parameter für unsere Simulation. Diese Parameter umfassen die Gravitationsbeschleunigung, die Länge und Masse jedes Pendels sowie den Zeitintervall für die Simulation.

G = 9.8  ## Gravitationsbeschleunigung in m/s^2
L1 = 1.0  ## Länge des ersten Pendels in m
L2 = 1.0  ## Länge des zweiten Pendels in m
L = L1 + L2  ## maximale Länge des kombinierten Pendels
M1 = 1.0  ## Masse des ersten Pendels in kg
M2 = 1.0  ## Masse des zweiten Pendels in kg
t_stop = 2.5  ## Anzahl der zu simulierenden Sekunden
history_len = 500  ## Anzahl der anzuzeigenden Trajektoriepunkte

Die Ableitungsfunktion definieren

Wir definieren eine Funktion, die die Ableitungen des Systems zu einem beliebigen Zeitpunkt berechnet. Diese Funktion erhält den aktuellen Zustand des Systems (die Winkel und Winkelgeschwindigkeiten jedes Pendels) und gibt die Ableitungen dieser Werte zurück.

def derivs(t, state):
    dydx = np.zeros_like(state)

    dydx[0] = state[1]

    delta = state[2] - state[0]
    den1 = (M1+M2) * L1 - M2 * L1 * cos(delta) * cos(delta)
    dydx[1] = ((M2 * L1 * state[1] * state[1] * sin(delta) * cos(delta)
                + M2 * G * sin(state[2]) * cos(delta)
                + M2 * L2 * state[3] * state[3] * sin(delta)
                - (M1+M2) * G * sin(state[0]))
               / den1)

    dydx[2] = state[3]

    den2 = (L2/L1) * den1
    dydx[3] = ((- M2 * L2 * state[3] * state[3] * sin(delta) * cos(delta)
                + (M1+M2) * G * sin(state[0]) * cos(delta)
                - (M1+M2) * L1 * state[1] * state[1] * sin(delta)
                - (M1+M2) * G * sin(state[2]))
               / den2)

    return dydx

Anfangsbedingungen festlegen und die gewöhnliche Differentialgleichung (ODE) integrieren

Wir legen nun die Anfangsbedingungen für unsere Simulation fest. Wir setzen die Anfangswinkel und -winkelgeschwindigkeiten jedes Pendels sowie das Zeitintervall für die Simulation. Anschließend integrieren wir die ODE mit der Euler-Methode, um die Position und Geschwindigkeit jedes Pendels in jedem Zeitschritt zu erhalten.

## create a time array from 0..t_stop sampled at 0.02 second steps
dt = 0.01
t = np.arange(0, t_stop, dt)

## th1 and th2 are the initial angles (degrees)
## w10 and w20 are the initial angular velocities (degrees per second)
th1 = 120.0
w1 = 0.0
th2 = -10.0
w2 = 0.0

## initial state
state = np.radians([th1, w1, th2, w2])

## integrate the ODE using Euler's method
y = np.empty((len(t), 4))
y[0] = state
for i in range(1, len(t)):
    y[i] = y[i - 1] + derivs(t[i - 1], y[i - 1]) * dt

Berechne die Position jedes Pendels

Wir verwenden nun die Position und Geschwindigkeit jedes Pendels in jedem Zeitschritt, um die x- und y-Koordinaten jedes Pendels zu berechnen.

x1 = L1*sin(y[:, 0])
y1 = -L1*cos(y[:, 0])

x2 = L2*sin(y[:, 2]) + x1
y2 = -L2*cos(y[:, 2]) + y1

Das Diagramm einrichten

Wir richten nun das Diagramm für unsere Simulation ein. Wir erstellen eine Figur mit einer x- und y-Begrenzung, die der maximalen Länge des Pendels entspricht, setzen das Seitenverhältnis gleich und fügen ein Gitter hinzu.

fig = plt.figure(figsize=(5, 4))
ax = fig.add_subplot(autoscale_on=False, xlim=(-L, L), ylim=(-L, 1.))
ax.set_aspect('equal')
ax.grid()

Die Animationsfunktion definieren

Wir definieren eine Funktion, die die Bewegung des Doppelpendels animiert. Diese Funktion erhält den aktuellen Zeitschritt und verwendet ihn, um die Position jedes Pendels zu aktualisieren. Sie aktualisiert auch die Bahn des zweiten Pendels.

line, = ax.plot([], [], 'o-', lw=2)
trace, = ax.plot([], [], '.-', lw=1, ms=2)
time_template = 'time = %.1fs'
time_text = ax.text(0.05, 0.9, '', transform=ax.transAxes)
history_x, history_y = deque(maxlen=history_len), deque(maxlen=history_len)

def animate(i):
    thisx = [0, x1[i], x2[i]]
    thisy = [0, y1[i], y2[i]]

    if i == 0:
        history_x.clear()
        history_y.clear()

    history_x.appendleft(thisx[2])
    history_y.appendleft(thisy[2])

    line.set_data(thisx, thisy)
    trace.set_data(history_x, history_y)
    time_text.set_text(time_template % (i*dt))
    return line, trace, time_text

Die Animation erstellen

Wir erstellen nun die Animation mit der FuncAnimation-Funktion aus Matplotlib.

ani = animation.FuncAnimation(
    fig, animate, len(y), interval=dt*1000, blit=True)
plt.show()

Zusammenfassung

In diesem Lab verwendet haben wir Python und Matplotlib, um die Bewegung eines Doppelpendels zu simulieren. Wir haben die Parameter für unsere Simulation festgelegt, eine Funktion definiert, um die Ableitungen des Systems zu berechnen, die gewöhnliche Differentialgleichung (ODE) mit der Euler-Methode integriert, die Position jedes Pendels in jedem Zeitschritt berechnet und eine Animation der Bewegung des Doppelpendels erstellt.