Der Lagrange-Formalismus

Um die Bewegungsgleichungen für physikalische Objekte zu bestimmen, kann das 2. Newtonsche Gesetz angewendet werden, also \(\vec{F} = m \cdot \vec{a}\) für geradlinige Bewegungen bzw. \(\vec{M} = J \cdot \vec{a_\theta}\) für Drehbewegungen. Jedoch ist es manchmal nicht sehr leicht, die dafür notwendigen Ausdrücke für die Kraft \(\vec{F}\) oder das Drehmoment \(\vec{M}\) zu finden. Oft ist es einfacher, Ausdrücke für die kinetische und potentielle Energie der beteiligten Objekte zu finden. Aus diesem Grunde kann es günstig sein, den Lagrange-Formalismus zur Herleitung der Bewegungsgleichungen heranzuziehen, wie wir das bereits hier getan haben.

Wie bereits erwähnt, kommt es beim Lagrange-Formalismus darauf an, geeignete Koordinaten zur Beschreibung des physikalischen Sachverhaltes zu wählen. Es muss ein minimaler Satz an Koordinaten gefunden werden, der den Bewegungszustand des Systems zu jedem Zeitpunkt eindeutig beschreiben kann. So reicht es z. B. beim Pendel aus, den Auslenkungswinkel \(\phi\) und dessen zeitliche Ableitung \(\dot{\phi}\) zu kennen: \(\phi\) gibt uns genau den Ort des Pendels an, während \(\dot{\phi}\) eindeutig seine Geschwindigkeit angibt. Weitere Variablen zur Beschreibung der Bewegung sind nicht notwendig, daher ist \(\phi\) die (einzige) generalisierte Koordinate für das Pendelproblem.

Kinetische und potentielle Energie des Fadenpendels

Wir nehmen an, dass ein Körper der Masse \(m\) an einem Pendel der Länge \(L\) befestigt ist. Die Masse des Pendelfadens ist vernachlässigbar, ebenso die Reibung (in der Luft oder am Aufhängepunkt). Die Anfangsauslenkung ist \(\phi_0\). Das Fadenpendel (auch als Mathematisches Pendel bezeichnet) kann als Punktmasse betrachtet werden, welche eine definierte Bahngeschwindigkeit \(v\) hat.

Um die Lagrange-Gleichung aufzustellen, brauchen wir die kinetische Energie des Körpers sowie das Potential (die potentielle Energie) in Abhängigkeit von den generalisierten Koordinaten und deren Zeitableitungen, hier also von \(\phi\) und \(\dot{\phi}\). Der Nullpunkt der potentiellen Energie kann beliebig gewählt werden. Hier ist es günstig, den tiefsten Punkt der Pendelbewegung mit \(E_{pot}=0\) zu verknüpfen. Die potentielle Energie kann dann von der generalisierten Koordinate, dem Winkel \(\phi\), abhängig gemacht werden.

Die potentielle Energie ist:

\[ E_{pot} = m \cdot g \cdot \left( L - L \cdot \cos{\phi} \right) = m \cdot g \cdot L \cdot \left( 1 - \cos{\phi} \right) \] Die kinetische Energie ist:

\[ E_{kin} = \frac{m}{2} \cdot v^2 = \frac{m}{2} \cdot \left( \omega \cdot L \right)^2 = \frac{m}{2} \cdot L^2 \cdot \dot{\phi}^2 \] Wir haben hier die Bahngeschwindigkeit durch die Kreisfrequenz \(\omega\) ersetzt: \(v = \omega \cdot L\) und dann \(\omega = \dot{\phi}\) verwendet, siehe Formelsammlung. Damit ist die kinetische Energie nur von der zeitlichen Ableitung der generalisierten Koordinate abhängig, alle anderen Größen sind (bekannte) Konstanten.

Lagrange-Gleichung für das Pendel

Die Kenntnis von kinetischer und potentieller Energie reicht aus, um die Lagrange-Funktion \(L\) zu finden. Es ist

\[ L ( \phi \; , \; \dot{\phi} ) = E_{kin} - E_{pot} = \frac{m}{2} \cdot L^2 \cdot \dot{\phi}^2 - m \cdot g \cdot L \cdot \left( 1 - \cos{\phi} \right) \hspace{2cm} (1) \] Die Lagrange-Gleichung lautet allgemein:

\[ \frac{d}{dt} \frac{\partial L}{\partial \dot{q}} - \frac{\partial L}{\partial q} = 0 \] Dabei ist \(q\) die generalisierte Koordinate, \(\dot{q}\) deren Zeitableitung. Sind mehrere generalisierte Koordinaten notwendig, muss für jede dieser Koordinaten eine Lagrange-Gleichung aufgestellt werden und das System von Differentialgleichungen gelöst werden.

Für das Pendel haben wir \(q \rightarrow \phi\) und \(\dot{q} \rightarrow \dot{\phi}\), also

\[ \frac{d}{dt} \frac{\partial L}{\partial \dot{\phi}} - \frac{\partial L}{\partial \phi} = 0 \]

Wir berechnen zunächst die partiellen Ableitungen:

\[ \frac{\partial L}{\partial \dot{\phi}}= \frac{\partial }{\partial \dot{\phi}} \Bigg[ E_{kin}-E_{pot}\Bigg] = \frac{\partial }{\partial \dot{\phi}} \Bigg[ \frac{m}{2} \cdot L^2 \cdot \dot{\phi}^2 - m \cdot g \cdot L \cdot \left( 1 - \cos{\phi} \right) \Bigg] = m \cdot L^2 \cdot \dot{\phi} \]

\[ \frac{\partial L}{\partial \phi}=\frac{\partial }{\partial \phi} \Bigg[ E_{kin}-E_{pot}\Bigg]=\frac{\partial }{\partial \phi}\Bigg[ \frac{m}{2} \cdot L^2 \cdot \dot{\phi}^2 - m \cdot g \cdot L \cdot \left( 1 - \cos{\phi} \right) \Bigg]= - m \cdot g \cdot L \cdot \sin{\phi} \] Nun noch die Zeitableitung des ersten Ausdruckes:

\[ \frac{d}{dt} \frac{\partial L}{\partial \dot{\phi}} = \frac{d}{dt} \Bigg[ m \cdot L^2 \cdot \dot{\phi} \Bigg] = m \cdot L^2 \cdot \ddot{\phi} \] Es taucht hier die Winkelbeschleunigung \(\ddot{\phi}\) auf. Die Lagrange-Gleichung wird nun:

\[ m \cdot L^2 \cdot \ddot{\phi} + m \cdot g \cdot L \cdot \sin{\phi} = 0 \] Nach Kürzen erhalten wir daraus die Bewegungsgleichung für das Pendel:

\[ \ddot{\phi} + \frac{g}{L} \cdot \sin{\phi} = 0 \hspace{2cm} (2) \]

In diesem einfachen Fall hätten wir auch das 2. Newtonsche Gesetz anwenden können und relativ leicht dieselbe Lösung erhalten.


Approximative Lösung der Bewegungsgleichung

Allerdings lässt sich die Bewegungsgleichung (2) nicht analytisch lösen. Wir können jedoch für kleine Auslenkungswinkel \(\sin{\phi} \approx \phi\) setzen, wir berücksichtigen damit nur das 1. Glied in der Reihenentwicklung von \(\sin{\phi}\). Dann haben wir für die Bewegungsgleichung die sehr viel einfachere Form

\[ \ddot{\phi} + \frac{g}{L} \cdot \phi = 0 \hspace{2cm} (3) \]

Diese Gleichung kann mit einem geeigneten Ansatz gelöst werden. Wir versuchen es mit

\[ \phi(t) = \phi_0 \cdot \cos{(\Omega \cdot t)} \] Hier ist \(\Omega\) eine noch unbekannte Frequenz. Die Größe \(\phi_0\) ist der Auslenkungswinkel bei \(t=0\) (denn für \(t=0\) ergibt der obige Ansatz \(\phi(0) = \phi_0\)). Zweimaliges Differenzieren dieses Ansatzes nach \(t\) ergibt:

\[ \ddot{\phi} = -\phi_0 \cdot \Omega^2 \cdot \cos{(\Omega \cdot t)} \] Wir setzen den Ansatz und diese 2. Ableitung in die Bewegungsgleichung ein:

\[ -\phi_0 \cdot \Omega^2 \cdot \cos{(\Omega \cdot t)} + \frac{g}{L} \cdot \phi_0 \cdot \cos{(\Omega \cdot t)} = 0 \] Nach Kürzen und Umstellen nach \(\Omega\) erhalten wir:

\[ \Omega = \sqrt{\frac{g}{L}} \] Die Lösung der approximativen Bewegungsgleichung (mit \(\sin{\phi} \approx \phi\)) ist also:

\[ \phi = \phi_0 \cdot \cos{\left(\sqrt{\frac{g}{L}} \cdot t \right)} \hspace{2cm} (4) \]

Ein Plot des Auslenkungswinkels sieht also so aus:

g = 9.81  # earth' acceleration
L = 2     # length of pendulum
m = 3     # mass of pendulum
phi0 = 5  # initial angle in degrees

phi <- function(t) {phi0*cos(sqrt(g/L)*t)}  # solution of equation of motion 

mtext = "Angle of pendulum"
curve(phi(x), from = 0, to = 20, n = 300, col = "red", lwd = 1.5, lty = 1, ylab = "angle(degrees)", xlab = "t(s)", main = mtext, font.main = 1)

Der Winkel beginnt bei \(\phi = 5^o\) und schwingt harmonisch zwischen \(+5^o\) und \(-5^o\), so wie wir das erwarten (wir hatten ja Reibungsverluste vernachlässigt).

Wir können auch noch einen Plot für die Bahngeschwindigkeit hinzufügen. Durch Differentiation der obigen Lösung erhalten wir:

\[ \omega = \dot{\phi} = -\phi_0 \cdot \sqrt{\frac{g}{L}} \cdot \sin{\left (\sqrt{\frac{g}{L}} \cdot t \right)} \] Die Bahngeschwindigkeit ist \(v = \omega \cdot L\), also

\[ v = -L \cdot \phi_0 \cdot \sqrt{\frac{g}{L}} \cdot \sin{(\sqrt{\frac{g}{L}} \cdot t)} \]

v_Bahn <- function(t) {-L*phi0*sqrt(g/L)*sin(sqrt(g/L)*t)}  # tangential velocity 

mtext = "Tangential velocity"
curve(v_Bahn(x), from = 0, to = 20, n = 300, col = "blue", lwd = 1.5, lty = 1, ylab = "velocity(m/s)", xlab = "t(s)", main = mtext, font.main = 1)

Die Bahngeschwindigkeit ist Null bei \(t=0\) und nimmt dann sofort negative Werte an, da das Pendel nach links schwingt (Der initiale Auslenkungswinkel war ja positiv, das Pendel befand sich also anfangs auf der rechten Seite.)

Wir könnten jetzt auch die kinetische und die potentielle Energie plotten, sowie deren Summe, die natürlich konstant sein muss. Allerdings wollen wir jetzt darauf verzichten, um eventuellen Enttäuschungen aus dem Wege zu gehen. Wir haben nämlich bisher die Gleichung nur ungefähr und für kleine Auslenkungswinkel gelöst, aufgrund der von uns gemachten Näherung \(\sin{\phi} \approx \phi\). Die zeitliche Ableitung der Gesamtenergie ist:

\[ \frac{dE_{tot}}{dt} = \frac{d}{dt} \Bigg[ E_{kin} + E_{pot} \Bigg] = \frac{d}{dt} \Bigg[ \frac{m}{2} \cdot L^2 \cdot \dot{\phi}^2 + m \cdot g \cdot L \cdot \left( 1 - \cos{\phi} \right) \Bigg] = m \cdot L^2 \cdot \dot{\phi} \cdot \ddot{\phi} + m \cdot g \cdot L \cdot \sin{\phi} \cdot \dot{\phi} = m \cdot L \cdot \dot{\phi} \Big[ L \cdot \ddot{\phi} + g \cdot \sin{\phi} \Big] \] Man beachte hier, dass z. B. beim Bilden der totalen zeitlichen Ableitung von \(\dot{\phi}^2\) die innere Ableitung nicht vergessen werden darf, also \(\frac{d}{dt} \left( \dot{\phi}^2 \right) = 2 \cdot \dot{\phi} \cdot \ddot{\phi}\).

Wir sehen in der obigen Gleichung, dass die Energie nur zeitlich konstant ist, wenn entweder \(\dot{\phi}=0\) oder \(L \cdot \ddot{\phi} + g \cdot \sin{\phi}=0\), also \(\phi\) der exakten Lösung der Lagrange-Gleichung, Gleichung (2), folgt.


Numerische Lösung der exakten Bewegungsgleichung

Eine numerische Lösung der exakten Bewegungsgkeichung (2) kann mit dem Package deSove berechnet werden. Dazu wird die Differentialgleichung 2. Ordnung (2) in eine System von zwei Differentialgleichungen 1. Ordnung umgeschrieben:

\[ \begin{align*} \dot{\phi} &= \omega \\[8pt] \dot{\omega} &= -g/L \cdot \sin{\phi} \hspace{2cm} (5) \end{align*} \]

Die rechte Seite des Differentialgleichungs-Systems muss als Funktion (hier: pendel) definiert werden und im Aufruf von ode mit dem Parameter func assoziiert werden:

library(deSolve)

phi0 =   25                             # initial angle (t=0) 
phi0 = pi/180*phi0                      # angle in radians, as required by the sine function in R
omega0 = 0                              # initial angular velocity (t=0) (1/s)
tmax = 20                               # maximum time to calculate (s) 
m = 3                                   # mass of the pendulum (kg)
g = 9.81                                # earth' acceleration (m/s^2)
L = 2                                   # length of pendulum (m)
params = c(g = g, L = L)                # parameters needed in the function definition   
times = seq(0, tmax, len = 401)         # time grid for calculation of the solution 
y_init = c(phi = phi0, omega = omega0)  # initial conditions for the differential equation 

pendel <- function(t, y, params) {      # function definition, right side of the system of differential equations
  dphi = y[2]                                  # dphi = omega                     omega ==> y[2]
  domega = -params["g"]/params["L"]*sin(y[1])  # domega = -g/L* sin(phi)            phi ==> y[1]
  return(list(c(dphi, domega))) 
}

solution <- ode(y_init, times, func = pendel, params)
head(solution)
##      time       phi      omega
## [1,] 0.00 0.4363323  0.0000000
## [2,] 0.05 0.4337437 -0.1034558
## [3,] 0.10 0.4260064 -0.2057604
## [4,] 0.15 0.4132063 -0.3057682
## [5,] 0.20 0.3954866 -0.4023436
## [6,] 0.25 0.3730474 -0.4943683

Mit Hilfe der erhaltenen Lösung können wir sofort die zeitlichen Verläufe des Auslenkungswinkels \(\phi\) und der Bahngeschwindigkeit \(v\) plotten (Wir plotten den Winkel in Grad):

phi_degrees = 180/pi*solution[,"phi"]   # angle in degrees
v_tan = solution[,"omega"]*L            # tangetial velocity
plot(solution[,"time"], phi_degrees, type = "l", col = "red", lwd = 2, xlab = "time", ylab = "angle", main = "Angle", font.main = 1)

plot(solution[,"time"], v_tan, type = "l", col = "blue", lwd = 2, xlab = "time", ylab = "tangential velocity", main = "Velocity", font.main = 1)

Selbst bei dem oben vereinbarten initialen Auslenkungswinkel von \(25^o\) sieht dies noch sehr nach Sinus- bzw. Kosinusfunktion aus. Sehen wir uns einen Vergleich der numerischen Lösung der exakten Bewegungsgleichung (2) und der analytischen Lösung der approximativen Gleichung (3) für den Startwinkel von \(25^o\) an:

t = seq(0, 20, len = 401)
phi_appr = 25*cos(sqrt(g/L)*t)            # approxiamtive solution for phi in degrees
v_appr = -L*25*sqrt(g/L)*sin(sqrt(g/L)*t) # approxiamtive tangential velocity

plot(solution[,"time"], phi_degrees, type = "l", col = "red", lwd = 2, xlab = "time", ylab = "angle", main = "Angle", font.main = 1)
lines(t, phi_appr, type = "l", col = "blue", lwd = 2)

Es scheint sich eine immer größer werdende Phasenverschiebung anzudeuten (rot = exakte Gleichung, blau = Näherung).

Jetzt berechnen wir auch die potentielle und kinetische Energie (siehe oben):

E_kin = m/2*L^2*solution[,"omega"]^2
E_pot = m*g*L*(1-cos(solution[,"phi"]))
E_tot = E_kin + E_pot

plot(solution[,"time"], E_kin, type = "l", col = "red", lwd = 2, xlab = "time", ylab = "energy", main = "Energy", font.main = 1)
lines(solution[,"time"], E_pot, type = "l", col = "blue", lwd = 2)

(rot = kinetische Energie, blau = potentialle Energie)

Die Gesamtenergie, berechnet mit der exakten Gleichung, sollte nun zeitlich konstant sein:

plot(solution[,"time"], E_tot, type = "l", col = "red", lwd = 2, xlab = "time", ylab = "energy", main = "Total Energy", font.main = 1)

Dies ist hinreichend genau für eine numerische Lösung.


uwe.menzel@matstat.org