Allgemeine Lösung für gedämpfte Schwingungen

Mehr realistische Modelle von Schwingungen beziehen Dämpfung ein. Oft ist die Dämpfung proportional zur Geschwindigkeit \(\dot{x} = dx/dt\). Es wird daher ein Reibungsterm \(2 \cdot \gamma \cdot m \cdot \dot{x}\) eingefügt, so dass wir folgende Bewegungsgleichung bekommen:

\[ m \cdot \ddot{x} = - k_R \cdot x - 2 \cdot \gamma \cdot m \cdot \dot{x} \]

Hier ist \(m\) wieder die Masse und \(2 \cdot \gamma\) der Reibungskoeffizient mit der Einheit \(1/s\) (wie eine Frequenz). (Warum der Faktor 2? - Das muss man nicht so machen, die Lösungen lassen sich dann aber eleganter ausdrücken, siehe unten …)

Das negative Vorzeichen resultiert daraus, dass die Dämpfung - ebenso wie die Rückstellkraft - die Bewegung behindert.

Diese Gleichung kann man mit dem Ansatz:

\[ x(t) = \exp{(\lambda \cdot t)} \]

lösen, wobei \(\lambda\) eine noch unbekannte (komplexe) Konstante ist, die durch Einsetzen in die Bewegungsgleichung gefunden werden soll. Wir brauchen daher wieder die zeitlichen Ableitungen des Ansatzes:

\[ \begin{align*} x(t) &= \exp{(\lambda \cdot t)} \\[6pt] \dot{x}(t) &= \lambda \cdot \exp{(\lambda \cdot t)} \\[6pt] \ddot{x}(t) &= \lambda^2 \cdot \exp{(\lambda \cdot t)} \end{align*} \]

Wenn wir dies in die obige Bewegungsgleichung einsetzen, erhalten wir nach Herauskürzen von \(\exp{(\lambda \cdot t)}\):

\[ \begin{align*} m \cdot \lambda^2 &= - k_R - 2 \cdot \gamma \cdot m \cdot \lambda \\[6pt] \lambda^2 &+ 2 \cdot \gamma \cdot \lambda + k_R/m = 0 \\[6pt] \lambda^2 &+ 2 \cdot \gamma \cdot \lambda + \omega^2 = 0 \end{align*} \]

Im letzten Schritt haben wir die Abkürzung \(\omega^2 = k_R/m\) benutzt, denn das war ja die Kreisfrequenz der ungedämpften Schwingung.

Aus der quadratischen Gleichung ergeben sich für \(\lambda\) nun zwei Lösungen:

\[ \lambda_{1/2} = -\gamma \; \pm \sqrt{\gamma^2 - \omega^2} \]

(Nun kommt die oben etwas künstlich eingeführte 2 nicht mehr vor …).

Wir haben nun zwei verschiedene Werte für \(\lambda\) und damit zwei Lösungsfunktionen \(\exp{(\lambda_1 \cdot t)}\) und \(\exp{(\lambda_2 \cdot t)}\).

Die (komplexe) Lösung der Bewegungsgleichung ist eine Linearkombination dieser beiden Lösungsfunktionen:

\[ x(t) = C_1 \cdot \exp{(\lambda_1 \cdot t)} + C_2 \cdot \exp{(\lambda_2 \cdot t)} \]


Schwache Dämpfung

Sogenannte schwache Dämpfung liegt vor, wenn \(\gamma < \omega\) ist. In diesem Fall sind \(\lambda_1\) und \(\lambda_2\) komplexe Zahlen.

Unter Beachtung von \(i^2=-1\) können wir die Wurzel so umformen: \(\sqrt{\gamma^2 - \omega^2} = \sqrt{i^2 \cdot (\omega^2 - \gamma^2)} = i \cdot \sqrt{\gamma^2 - \omega^2}\). Der Wurzelausdruck ist nun reell und größer Null.

Wir erhalten dann für die beiden \(\lambda\):

\[ \lambda_{1/2} = -\gamma \; \pm i \cdot \sqrt{\omega^2 - \gamma^2} \] und damit für die Lösung:

\[ x(t) = \exp{(-\gamma \cdot t)} \cdot \Bigg[ C_1 \cdot \exp{ \left( i \cdot \sqrt{\omega^2 - \gamma^2} \cdot t \right)} + C_2 \cdot \exp{ \left(- i \cdot \sqrt{\omega^2 - \gamma^2} \cdot t \right)}\Bigg] \hspace{1.5cm} (1) \]

Werden statt \(C_1\) und \(C_2\) neue Konstanten \(C_3 = C_1 + C_2\) und \(C_4 = i \cdot (C_1-C_2)\) eingeführt, erhält man eine reelle Lösung:

\[ x(t) = \exp{(-\gamma \cdot t)} \cdot \Bigg[ C_3 \cdot \cos{ \left( \sqrt{\omega^2 - \gamma^2} \cdot t \right)} + C_4 \cdot \sin{ \left( \sqrt{\omega^2 - \gamma^2} \cdot t \right)}\Bigg] \hspace{2.8cm} (2) \]

Wie die Umrechnung gemacht wird, kann man bei Bedarf in Anhang nachlesen.

Schließlich kann man mit einer weiteren Umrechnung auch die folgende Form erhalten (mit der Phase \(\delta\)):

\[ x(t) = A \cdot \exp{(-\gamma \cdot t)} \cdot \cos{ \left( \sqrt{\omega^2 - \gamma^2} \cdot t - \delta\right)} = A \cdot \exp{(-\gamma \cdot t)} \cdot \cos{ \left( \alpha \cdot t - \delta\right)} \hspace{1.5cm} (3) \]

Wir benutzen im Folgenden oft die Abkürzung \(\alpha =\sqrt{\omega^2 - \gamma^2}\), um Schreibarbeit einzusparen.

Auch die zu Gleichung (3) führende Umrechnung gibt es im Anhang.

Formel (3) ist vielleicht am intuitivsten. Die Größe \(A\) ist die maximale Auslenkung, denn der Kosinus bewegt sich im Intervall \(-1 \dots +1\), während der exponentielle Term \(\le 1\) ist.
Der Kosinus-Term beschreibt eine Schwingung (mit etwas kleinerer Frequenz als \(\omega\)), während der Term \(\exp{\left(-\gamma \cdot t \right)}\) dafür sorgt, dass die Auslenkung mit der Zeit immer kleiner wird.

Die Amplitude \(A\) und die Phase \(\delta\) in Gleichung (3) müssen aus den Anfangsbedingungen bestimmt werden. Mit \(x(0)=x_0\) und \(v(0)=v_0\) mit bekannten \(x_0\) und \(v_0\) erhalten wir:

\[ \begin{align*} \delta &= \arctan{\Big[ \frac{v_0 + \gamma \cdot x_0}{\alpha \cdot x_0} \Big]} \\[6pt] A &= \sqrt{x_0^2 + \frac{\left(v_0 + \gamma \cdot x_0 \right)^2}{\alpha^2}} \end{align*} \]

Die Herleitung dieser Ausdrücke wird im Anhang gezeigt.

Ausgehend von \(x(t)\) können wir auch die Geschwindigkeit berechnen. Diese ergibt sich aus der zeitlichen Ableitung von Gleichung (3):

\[ \begin{align*} x(t) &= A \cdot \exp{(-\gamma \cdot t)} \cdot \cos{ \left( \alpha \cdot t - \delta\right)} \\[6pt] \frac{dx}{dt} = v(t) &= -A \cdot \exp{(-\gamma \cdot t)} \cdot \Bigg[ \gamma \cdot \cos{ \left( \alpha \cdot t - \delta\right)} + \alpha \cdot \sin{ \left( \alpha \cdot t - \delta\right)} \Bigg] \end{align*} \]


xt <- function(t, gamma, omega, x0, v0) {  # extension
  alpha = sqrt(omega^2-gamma^2)
  delta = atan((v0+gamma*x0)/(alpha*x0))
  A = sqrt(x0^2 + (v0 + gamma*x0)^2/alpha^2)
  A*exp(-gamma*t)*cos(alpha*t - delta) 
} 

vt <- function(t, gamma, omega, x0, v0) {   # velocity  
  alpha = sqrt(omega^2-gamma^2)
  delta = atan((v0+gamma*x0)/(alpha*x0))
  A = sqrt(x0^2 + (v0 + gamma*x0)^2/alpha^2)
  -A*exp(-gamma*t)*(gamma*cos(alpha*t - delta) + alpha * sin(alpha*t - delta))
}

x0 = 1                # initial extension
v0 = 0                # initial velocity 
T = 1                 # period
omega = 2*pi/T        # intrinsic frequency 
gamma = c(0.5, 0.7)   # damping coefficient; gamma < omega ==> weak attenuation

mtxt = "Extension for damped oscillation"
curve(xt(x, gamma[1], omega, x0, v0), 
      from = 0, to = 8, n = 300, col = "red",  lwd = 1.5, lty = 1, 
      ylab = "x(t)", xlab = "t(s)", main = mtxt, font.main = 1)
curve(xt(x, gamma[2], omega, x0, v0), 
      from = 0, to = 8, n = 300, col = "blue", lwd = 1.5, lty = 3, add = TRUE)
labels <- c(parse(text = sprintf("gamma == %f", gamma[1])),
            parse(text = sprintf("gamma == %f", gamma[2])))
legend("topright", labels, lwd=2, lty=c(1,3), col=c("red", "blue"))

mtxt = "Velocity for damped oscillation"
curve(vt(x, gamma[1], omega, x0, v0), 
      from = 0, to = 8, n = 300, col = "red",  lwd = 1.5, lty = 1, ylab = "v(t)", xlab = "t(s)", main = mtxt, font.main = 1)
curve(vt(x, gamma[2], omega, x0, v0), 
      from = 0, to = 8, n = 300, col = "blue", lwd = 1.5, lty = 3, add = TRUE)
labels <- c(parse(text = sprintf("gamma == %f", gamma[1])),
            parse(text = sprintf("gamma == %f", gamma[2])))
legend("topright", labels, lwd=2, lty=c(1,3), col=c("red", "blue"))

Wir kontrollieren noch, ob die Anfangsbedingungen tatsächlich erfüllt sind:

cat(paste("x0 =", xt(0, gamma[1], omega, x0, v0), "  v0 =", round(vt(0, gamma[1], omega, x0, v0), 15))) 
## x0 = 1   v0 = 0
cat(paste("x0 =", xt(0, gamma[2], omega, x0, v0), "  v0 =", round(vt(0, gamma[2], omega, x0, v0), 15))) 
## x0 = 1   v0 = 0

Das stimmt, jedenfalls wenn auf 15 Stellen (bei der Geschwindigkeit) gerundet wird.
(Es ergibt sich für \(\gamma_2\) ein Rechenfehler von \(\approx 10^{-16} \; m/s\) für \(v_0\)).


Potentielle und kinetische Energie bei schwacher Dämpung

Wie bereits hier gezeigt ist die potentielle Energie:

\[ W_{pot} = \frac{k_R}{2}\cdot x^2 \]

(Wir erinnern uns: \(k_R\) war der Proportionalitätsfaktor zwischen Kraft und Auslenkung \(F = -k_R \cdot x\).)

Wir brauchen also nur den berechneten Ausdruck für \(x(t)\) einsetzen:

\[ \begin{align*} W_{pot} &= \frac{k_R}{2} \cdot \Bigg[ A \cdot \exp{(-\gamma \cdot t)} \cdot \cos{ \left( \alpha \cdot t \right)} \Bigg]^2 \\[6pt] &= \frac{k_R}{2} \cdot A^2 \cdot \exp{(- 2 \cdot\gamma \cdot t)} \cdot \cos^2{ \left( \alpha \cdot t \right)} \hspace{1.5cm} (4) \end{align*} \]


Mit Hilfe der oben berechneten Geschwindigkeit erhalten wir auch die kinetische Energie:

\[ W_{kin} = \frac{m}{2} \cdot v^2 = \frac{m}{2} \cdot A^2 \cdot \exp{(-2 \cdot\gamma \cdot t)} \cdot \Bigg[ \gamma \cdot \cos{ \left( \alpha \cdot t - \delta\right)} + \alpha \cdot \sin{ \left( \alpha \cdot t - \delta\right)} \Bigg]^2 \hspace{1.5cm} (5) \]

Im folgenden Code-Snippet benutzen wir einfach die bekannten Ausdrücke für \(x(t)\) und \(v(t)\):

Wpot <- function(t, kR, m, gamma, x0, v0) {
  omega = sqrt(kR/m)
  kR/2*xt(t, gamma, omega, x0, v0)^2 
}

Wkin <- function(t, kR, m, gamma, x0, v0) {
  omega = sqrt(kR/m)
  m/2*vt(t, gamma, omega, x0, v0)^2
}  

Wges <- function(t, kR, m, gamma, x0, v0) { 
  Wpot(t, kR, m, gamma, x0, v0) + Wkin(t, kR, m, gamma, x0, v0) 
} 

x0 = 1                # initial extension                   m
v0 = 0                # initial velocity                    m/s
kR = 4*pi^2           # characteristic constant             kg/s^2
m = 1                 # mass                                kg
gamma = 0.5           # damping coefficient; gamma < omega ==> weak attenuation

mtxt = "Energy for damped  Oscillation"
curve(Wpot(x, kR, m, gamma, x0, v0), from = 0, to = 3, n = 300, col = "red",   lwd = 1.5, lty = 1, 
      ylab = "Energy (J)", xlab = "t(s)", main = mtxt, font.main = 1)
curve(Wkin(x, kR, m, gamma, x0, v0), from = 0, to = 3, n = 300, col = "blue",  lwd = 1.5, lty = 1, add = TRUE)
curve(Wges(x, kR, m, gamma, x0, v0), from = 0, to = 3, n = 300, col = "green", lwd = 1.5, lty = 3, add = TRUE)
omega = sqrt(kR/m)             # for subtitle only
alpha = sqrt(omega^2-gamma^2)  # for subtitle only 
mtext(substitute(T == a ~~ omega == b ~~ gamma == c ~~ alpha == d, 
      list(a=signif(T,3), b=signif(omega,3), c=signif(gamma,3), d=signif(alpha,3))))
labels = c(expression('W'[pot]), expression('W'[kin]), expression('W'[ges])   )
legend("topright",legend=labels, bty="n", col = c("red", "blue", "green"), lty= c(1,1,3), lwd = 2)

Es tritt ein Verlust von potentieller und kinetischer Energie ein.
Dieser ist besonders groß, wo die Geschwindigkeit groß ist, denn dort sind ja die Reibungsverluste groß.

Der Energie”verlust” ist das Integral über die Reibungskraft (Es handelt sich natürlich nicht um einen wirklichen Verlust, sondern um eine Umwandlung in Reibungswärme):

\[ W_{R} = -\int F_R(x) \;dx = 2 \cdot \gamma \cdot m \cdot \int v(x) \; dx \] Allerdings haben wir \(v(x)\) nicht. Wir können aber \(dx = v \cdot dt\) ausnutzen, so dass wir für das Integral schreiben können:

\[ W_{R} = 2 \cdot \gamma \cdot m \cdot \int v \cdot \; v \cdot dt = 2 \cdot \gamma \cdot m \cdot \int_0^t v(t)^2 \cdot \; dt \]

Die Ausführung der Integration überlassen wir allerdings gern der Software.
Wir berechnen hier das Integral für 300 verschiedene obere Grenzen, um die Zeitabhängigkeit darzustellen:

x0 = 1                             # initial extension                   m
v0 = 0                             # initial velocity                    m/s
gamma = 0.5                        # attenuation  parameter              1/s
kR = 4*pi^2                        # characteristic constant             kg/s^2
m = 1                              # mass                                kg
T = 2*pi*sqrt(m/kR)                # period                              s
omega = 2*pi/T                     # frequency of undamped oscillation   1/s 
alpha = sqrt(omega^2 - gamma^2)    # effective frequency                 1/s

to = seq(0, 3, length = 300)       # upper limits for integration
Wr = numeric()                     # vector to store the friction energy
i = 0
for (time in to) {
  i = i + 1
  integrand <- function(t, gamma = gamma, omega = omega, x0 = x0, v0 = v0) {  
    2*gamma*m*vt(t, gamma, omega, x0, v0)^2
  }
  res = integrate(integrand, lower = 0, upper = time, gamma = gamma, omega = omega, x0 = x0, v0 = v0)
  Wr[i] = res$value
}

mtxt = "Friction Energy and Sum of Kinetic and Potential Energy"
plot(to, Wr,  col = "red",   lwd = 1.5, lty = 1, type = "l", 
     ylab = "Energy (J)", xlab = "t(s)", main = mtxt, font.main = 1)
curve(Wges(x, kR, m, gamma, x0, v0), 
      from = 0, to = 3, n = 300, col = "green", lwd = 2, lty = 3, add = TRUE)
labels = c(expression(paste('W'[pot], "+", 'W'[kin])), expression('W'[friction]))
legend("right", labels, col = c("green", "red"), lwd = 2, bty="n", lty = c(3,1))

Die Summe aller beteiligten Energien ist offensichtlich konstant, wie wir das erwarten.


Starke Dämpfung und Aperiodischer Grenzfall

Bei starker Dämpfung (\(\gamma > \omega\)) erhalten wir eine exponentiell abfallende Funktion ohne Schwingungsanteil, denn \(\lambda_1\) und \(\lambda_2\) sind dann reell und positiv:

\[ \lambda_{1/2} = -\gamma \; \pm \sqrt{\gamma^2 - \omega^2} \]

Die Lösung ist daher:

\[ \begin{align*} x(t) &= C_1 \cdot \exp{(\lambda_1 \cdot t)} + C_2 \cdot \exp{(\lambda_2 \cdot t)} \\[6pt] &= C_1 \cdot \exp{ \left( -\gamma \cdot t + \sqrt{\gamma^2 - \omega^2} \cdot t \right)} + C_2 \cdot \exp{ \left( -\gamma \cdot t - \sqrt{\gamma^2 - \omega^2} \cdot t \right)} \\[6pt] &= \exp{ \left( -\gamma \cdot t \right)} \cdot \Big[ C_1 \cdot \exp{(\beta \cdot t)} + C_2 \cdot \exp{(-\beta \cdot t)} \Big] \end{align*} \]

Der Einfachheit halber wurde in der letzten Zeile die Abkürzung \(\beta = \sqrt{\gamma^2 - \omega^2}\) eingeführt. \(\beta\) ist positiv und reell, da \(\gamma > \omega\).

Die Geschwindigkeit ergibt sich durch Ableitung von \(x(t)\):

\[ \begin{align*} x(t) &= \exp{ \left( -\gamma \cdot t \right)} \cdot \Big[ C_1 \cdot \exp{(\beta \cdot t)} + C_2 \cdot \exp{(-\beta \cdot t)} \Big] \\[6pt] \frac{dx}{dt} = v(t) &= \exp{ \left( -\gamma \cdot t \right)} \cdot \Big[ C_1 \cdot (\beta - \gamma) \cdot \exp{(\beta \cdot t)} - C_2 \cdot (\beta + \gamma)\cdot \exp{(-\beta \cdot t)} \Big] \end{align*} \]

Die noch unbekannten Konstanten \(C_1\) und \(C_2\) können nun mit Hilfe der Anfangsbedingungen \(x(0)=x_0\) und \(v(0)=v_0\) bestimmt werden:

\[ \begin{align*} x(0) = x_0 &= C_1 + C_2 \\[6pt] v(0) = v_0 &= \beta \cdot (C_1 - C_2) - \gamma \cdot (C_1 + C_2) \end{align*} \] Damit ergibt sich für die Konstanten:

\[ \begin{align*} C_1 &= \frac{(\beta + \gamma) \cdot x_0 + v_0}{2 \beta} \\[6pt] C_2 &= \frac{(\beta - \gamma) \cdot x_0 - v_0}{2 \beta} \end{align*} \]

Für \(x(0)=1\) und \(v(0) = 0\) ergeben sich folgende Kurven:

x0 = 1                             # amplitude (initial extension)       m
v0 = 0                             # initial velocity                    m/s
gamma = c(10, 15, 20)              # attenuation  parameter              1/s
kR = 4*pi^2                        # characteristic constant             kg/s^2
m = 1                              # mass                                kg
T = 2*pi*sqrt(m/kR)                # period of undamped oscillation      s
omega = 2*pi/T                     # frequency of undamped oscillation   1/s 

xt_strong <- function(t, gamma, omega, x0, v0) { 
  beta = sqrt(gamma^2 - omega^2) 
  C1 = ((beta+gamma)*x0 + v0)/(2*beta)
  C2 = ((beta-gamma)*x0 - v0)/(2*beta)
  exp(-gamma*t)*(( C1*exp(beta*t) + C2*exp(-beta*t)))
} 

mtxt = "Attenuated  Oscillation, strong attenuation"
curve(xt_strong(x, gamma[1], omega, x0, v0), 
      from = 0, to = 2, n = 300, col = "red", lwd = 1.5, lty = 1, 
      ylab = "x(t)", xlab = "t(s)", main = mtxt, font.main = 1)
curve(xt_strong(x, gamma[2], omega, x0, v0), 
      from = 0, to = 2, n = 300, col = "blue",  lwd = 1.5, lty = 1, add = TRUE)
curve(xt_strong(x, gamma[3], omega, x0, v0), 
      from = 0, to = 2, n = 300, col = "green",  lwd = 1.5, lty = 1, add = TRUE)
labels <- c(parse(text = sprintf("gamma == %f", gamma[1])),
            parse(text = sprintf("gamma == %f", gamma[2])),
            parse(text = sprintf("gamma == %f", gamma[3])))
legend("topright", labels, col = c("red", "blue", "green"), lwd = 2, lty = c(1,1,1))

Ein etwas unerwartetes Resultat: mit steigender Dämpfung (\(\gamma\)) fällt die Kurve langsamer, da sich der Körper langsamer in Richtung Ruhelage bewegt.

Der sogenannte aperiodische Grenzfall tritt auf, wenn \(\gamma = \omega\) ist. Durch den oben angewendeten Ansatz \(x(t) = \exp{(\lambda \cdot t)}\) ergibt sich jetzt nur eine Lösung \(x_1(t) = \exp{(\lambda \cdot t)} = \exp{(-\gamma \cdot t)}\). Im allgemeinen brauchen wir jedoch zwei (linear unabhängige) Lösungen für eine Differentialgleichung 2. Ordnung. Durch Einsetzen in die Ausgangsgleichung können wir uns überzeugen, dass auch die Funktion \(x_2(t) = t \cdot \exp{(-\gamma \cdot t)}\) die Differentialgleichung löst (siehe Anhang).
Um uns zu überzeugen, dass \(x_1(t)\) und \(x_2(t)\) linear unabhängig sind, können wir die Wronski-Determinante berechnen, siehe Anhang.

Die Lösung ist dann wieder eine Linearkombination beider Teillösungen, also

\[ x(t) = C_1 \cdot e^{-\gamma \cdot t} + C_2 \cdot t \cdot e^{-\gamma \cdot t} = e^{-\gamma \cdot t} \cdot \left( C_1 + C_2 \cdot t \right) \]

Die Geschwindigkeit \(\dot{x}(t)\) ergibt sich wieder aus der Zeitableitung von \(x(t)\):

\[ \dot{x}(t) = -\gamma \cdot e^{-\gamma \cdot t} \cdot \left( C_1 + C_2 \cdot t \right) + e^{-\gamma \cdot t} \cdot C_2 = e^{-\gamma \cdot t} \cdot \left( C_2 - \gamma \cdot C_1 - \gamma \cdot C_2 \cdot t\right) \]

Die beiden Konstanten \(C_1\) und \(C_2\) müssen wieder aus den Anfangsbedingungen \(x(0)=x_0\) und \(v(0)=v_0\) bestimmt werden:

\[ \begin{align*} x_0 = x(0) &= C_1 \\[6pt] v_0 = v(0) &= C_2 - \gamma \cdot C_1 \end{align*} \] Daraus ergibt sich:

\[ \begin{align*} C_1 &= x_0 \\[6pt] C_2 &= v_0 + \gamma \cdot x_0 \end{align*} \]

Die Lösung lautet nun also:

\[ x(t) = e^{-\gamma \cdot t} \cdot \Big[ x_0 + \left( v_0 + \gamma \cdot x_0 \right) \cdot t \Big] \]

Dieses Ergebnis können wir zur obigen Abbildung für die starke Dämpfung hinzufügen:

x0 = 1                             # amplitude (initial extension)       m
v0 = 0                             # initial velocity                    m/s
gamma = c(10, 15, 20)              # attenuation  parameter              1/s
kR = 4*pi^2                        # characteristic constant             kg/s^2
m = 1                              # mass                                kg
T = 2*pi*sqrt(m/kR)                # period of undamped oscillation      s
omega = 2*pi/T                     # frequency of undamped oscillation   1/s 

# xt_strong defined above
xt_ape <- function(t, omega, x0, v0) {
  # gamma = omega
  exp(-omega*t)*(x0 + (v0 + omega*x0)*t)
}  

mtxt = "Strong attenuation and critical damping"
curve(xt_strong(x, gamma[1], omega, x0, v0), 
      from = 0, to = 2, n = 300, col = "red", lwd = 1.5, 
      lty = 1, ylab = "x(t)", xlab = "t(s)", main = mtxt, font.main = 1)
curve(xt_strong(x, gamma[2], omega, x0, v0), 
      from = 0, to = 2, n = 300, col = "blue",  lwd = 1.5, lty = 1, add = TRUE)
curve(xt_strong(x, gamma[3], omega, x0, v0), 
      from = 0, to = 2, n = 300, col = "green",  lwd = 1.5, lty = 1, add = TRUE)
curve(xt_ape(x, omega, x0, v0), 
      from = 0, to = 2, n = 300, col = "black",  lwd = 2, lty = 3, add = TRUE)
labels <- c(parse(text = sprintf("gamma == %f", gamma[1])),
            parse(text = sprintf("gamma == %f", gamma[2])),
            parse(text = sprintf("gamma == %f", gamma[3])),
            parse(text = sprintf("gamma == omega")))
legend("topright", labels, col = c("red", "blue", "green", "black"), lwd = 2, lty = c(1,1,1,3))

Im aperiodischen Grenzfall \(\gamma = \omega\) kehrt der Körper am schnellsten in seine Ruhelage zurück. Dies ist sehr wichtig z. B. beim Bau von Messgeräten.


Die Phasenraumdarstellung für die gedämpfte Schwingung

Die Dynamik der gedämpften Schwingung kann auch im Phasenraum untersucht werden. Ein Phasenraum-Plot zeigt die Geschwindigkeit \(v\) als Funktion der Auslenkung \(x\). Eine analytische Beziehung \(v(x)\) könnte man im Prinzip herleiten, indem die Zeit aus den beiden Gleichungen \(x(t)\) und \(v(t)\) eliminiert wird. In vielen Fällen wird dies jedoch nicht unkompliziert sein. Wir werden daher für eine Reihe von Zeitpunkten \(t\) jeweils die Werte \(x(t)\) und \(v(t)\) berechnen und das Ergebnis plotten. Wir betrachten zunächst den Fall der schwachen Dämpfung mit folgenden, oben für \(x(0) = x_0\) und \(v(0) = v_0\) hergeleiteten Beziehungen:

\[ \begin{align*} x(t) &= A \cdot \exp{(-\gamma \cdot t)} \cdot \cos{ \left( \alpha \cdot t - \delta\right)} \\[6pt] v(t) &= -A \cdot \exp{(-\gamma \cdot t)} \cdot \Bigg[ \gamma \cdot \cos{ \left( \alpha \cdot t - \delta\right)} + \alpha \cdot \sin{ \left( \alpha \cdot t - \delta\right)} \Bigg] \end{align*} \]

Hier war \(\alpha = \sqrt{\omega^2 - \gamma^2}\) die effektive Frequenz des Systems.

x0 = 1                # initial extension                   m
v0 = 0                # initial velocity                    m/s
T = 1                 # period                              s
omega = 2*pi/T        # intrinsic (undamped) frequency      1/s

t = seq(0, 10, length = 501)
x = xt(t, gamma = 0.5, omega, x0, v0)
v = vt(t, gamma = 0.5, omega, x0, v0)

cols = rep("blue", length(t))
cols[1] = "red"
cols[length(t)] = "red"
mtxt = "Phase portrait, damped oscillation"

plot(x, v, pch = 20, type = "b", lwd = 1.5, cex = 0.6, col = cols, main = mtxt, font.main = 1)
mtext(substitute(gamma == a, list(a=0.5)))
abline(h = 0, lty = 3, lwd = 1, col = "darkgrey")
abline(v = 0, lty = 3, lwd = 1, col = "darkgrey")

Der Anfangspunkt bei \(t=0\) ist rot gekennzeichnet (in der Abbildung rechts). Hier ist die Auslenkung 1 und die Geschwindigkeit 0, wie durch die Anfangsbedingungen vorgegeben. Die Auslenkung wird danach immer kleiner, währenddessen die Geschwindigkeit negativ ist, d. h. der Schwinger bewegt sich zurück in Richtung Ruhelage. Sobald \(x=0\) ist, nimmt der Betrag der Geschwindigkeit wieder ab und wird schließlich Null, sobald die maximale Auslenkung auf der linken Seite erreicht ist. Bedingt durch die Dämpfung werden im Laufe der Zeit die maximalen Auslenkungen und die maximalen Geschwindigkeiten immer kleiner. Schließlich kommt das System bei \(x=0\) und \(v=0\) zur Ruhe (roter Punkt in der Mitte).


uwe.menzel@matstat.org