Newmark-beta-Verfahren

Newmark-beta-Verfahren sind Methoden zur impliziten numerischen Integration von Differenzialgleichungen. Die Verfahren gehören zu den Einschrittverfahren, da zur Berechnung der Werte zur Zeit t_n+1 nur die Werte des vorangegangenen Zeitschritts zur Zeit t_n benötigt werden. Dabei werden zwei Parameter β (beta) und γ eingeführt, mit denen die Stabilität und die Genauigkeit des Verfahrens gesteuert werden. Die Verfahrensklasse ist in der numerischen Analyse der Dynamik von Festkörpern wie in der Finite-Elemente-Methode weit verbreitet. Benannt ist sie nach Nathan M. Newmark, der sie 1959 für die Anwendung in der Strukturdynamik entwickelte.^[1]

Herleitung

Annahme linearer oder konstanter Beschleunigung

Im Zeitintervall ${\displaystyle [t_{0},t_{e}]}$, in dem eine Lösung ${\displaystyle x(t)}$ einer Differenzial­gleichung zweiter Ordnung in der Zeit gesucht wird, sei eine streng monoton steigende Folge von Zeitpunkten ${\displaystyle \{t_{0},t_{1},\dotsc ,t_{e}\}}$ vorgegeben, zu denen die Lösung ${\displaystyle x(t)}$ berechnet werden soll. Der Wert der Variable ${\displaystyle x_{n}}$, ihre Rate ${\displaystyle {\dot {x}}_{n}}$ und Beschleunigung ${\displaystyle {\ddot {x}}_{n}}$ seien zur Zeit t_n bekannt. Die Beschleunigung wird im Intervall ${\displaystyle t\in [t_{n},t_{n+1}]}$ linear interpoliert, siehe Bild:

${\displaystyle {\ddot {x}}^{h}(t)={\frac {t_{n+1}-t}{t_{n+1}-t_{n}}}{\ddot {x}}_{n}+{\frac {t-t_{n}}{t_{n+1}-t_{n}}}{\ddot {x}}_{n+1}}$		(I)

worin ${\displaystyle x^{h}}$ eine Näherungslösung der gesuchten Funktion ${\displaystyle x}$ bezeichnet. Integration über die Zeit liefert mit ${\displaystyle \Delta t=t-t_{n}}$:

${\displaystyle {\dot {x}}^{h}(t)={\dot {x}}_{n}+\int _{t_{n}}^{t}{\ddot {x}}^{h}(\tau )\,\mathrm {d} \tau ={\dot {x}}_{n}+\Delta t\,[(1-\gamma ){\ddot {x}}_{n}+\gamma \,{\ddot {x}}^{h}(t)]}$		(II)

${\displaystyle x^{h}(t)=x_{n}+\int _{t_{n}}^{t}{\dot {x}}^{h}(\tau )\mathrm {d} \tau =x_{n}+\Delta t\,{\dot {x}}_{n}+\Delta t^{2}\left[\left({\frac {1}{2}}-\beta \right){\ddot {x}}_{n}+\beta \,{\ddot {x}}^{h}(t)\right]}$		(III)

Mit

β = 1/6 und γ = 1/2

sind diese Formeln exakt und liefern das lineare Beschleunigungsverfahren. Unter der Voraussetzung, dass die Extremwerte der Beschleunigung im Intervall [t_n,t_n+1] an den Grenzen des Intervalls auftreten, stellen die Integrale in Gleichungen (II) und (III) eine abgebrochene Taylorreihe mit Restglied dar, wobei mit 0≤β≤1 und 0≤γ≤1 andere Approximationen ${\displaystyle x^{h}}$ gefunden werden. So können auch andere Werte für die Konstanten β und γ motiviert werden.

Start der Berechnung

Der Newmark-Algorithmus startet zur Zeit t=t₀ mit n=0. Zumeist wird angenommen, dass für t≤t₀ die Beschleunigungen verschwinden. Mit dieser Annahme ist der Algorithmus unter Vorgabe der Anfangswerte x₀ und Anfangsgeschwindigkeit ẋ₀ selbststartend, d. h. die Anfangsbeschleunigungen brauchen nicht in einem ersten Schritt berechnet zu werden.

Aktualisierung der Variablen

Mit dem Newmark-Algorithmus werden aus gegebenen Werten ${\displaystyle x_{n},{\dot {x}}_{n}}$ und ${\displaystyle {\ddot {x}}_{n}}$ zur Zeit t_n die entsprechenden Werte zur Zeit t_n+1 berechnet. Die im Intervall ${\displaystyle [t_{n},t_{n+1}]}$ liegenden Werte können mit den #Gleichungen (I) bis (III) interpoliert werden. Mit ${\displaystyle t=t_{n+1}}$ und ${\displaystyle {\ddot {x}}^{h}(t_{n+1})={\ddot {x}}_{n+1}}$ bekommt man aus Gleichungen (II) und (III):

${\displaystyle {\dot {x}}^{h}(t_{n+1})=:{\dot {x}}_{n+1}={\dot {x}}_{n}+\Delta t[(1-\gamma )\,{\ddot {x}}_{n}+\gamma \,{\ddot {x}}_{n+1}]}$		(IV)

${\displaystyle x^{h}(t_{n+1})=:x_{n+1}=x_{n}+\Delta t{\dot {x}}_{n}+\Delta t^{2}\left[\left({\frac {1}{2}}-\beta \right){\ddot {x}}_{n}+\beta \,{\ddot {x}}_{n+1}\right]}$		(V)

Die beiden Gleichungen (IV) und (V) enthalten drei Unbekannte ${\displaystyle x_{n+1},{\dot {x}}_{n+1}}$ und ${\displaystyle {\ddot {x}}_{n+1}}$. Die dritte zum Abschluss benötigte Gleichung liefert die zu lösende Differenzial­gleichung. Bei ${\displaystyle \beta \neq 0}$ kann auch ${\displaystyle x_{n+1}}$ als primäre Unbekannte gewählt werden:

${\displaystyle {\ddot {x}}_{n+1}={\frac {1}{\beta \Delta t^{2}}}(x_{n+1}-x_{n})-{\frac {1}{\beta \Delta t}}{\dot {x}}_{n}-{\frac {{\frac {1}{2}}-\beta }{\beta }}{\ddot {x}}_{n}}$

${\displaystyle {\dot {x}}_{n+1}={\frac {\gamma }{\beta \Delta t}}(x_{n+1}-x_{n})+\left(1-{\frac {\gamma }{\beta }}\right){\dot {x}}_{n}+\Delta t{\frac {\beta -{\frac {1}{2}}\gamma }{\beta }}{\ddot {x}}_{n}}$.

Sind einmal die Werte ${\displaystyle x_{n+1},{\dot {x}}_{n+1}}$ und ${\displaystyle {\ddot {x}}_{n+1}}$ berechnet, wird der Zähler ${\displaystyle n}$ inkrementiert und die Berechnung fortgesetzt, bis das Ende des interessierenden Zeitintervalls erreicht ist.

Spezialfälle

Konstante Durchschnittsbeschleunigungsverfahren

Die ursprüngliche Form des Newmark-Verfahrens entspricht einer konstanten mittleren Beschleunigung

${\displaystyle {\ddot {x}}^{h}(t)={\frac {1}{2}}({\ddot {x}}_{n+1}+{\ddot {x}}_{n})}$

siehe Bild in der #Herleitung. Vergleich mit den obigen #Gleichungen (IV) und (V) führt auf

${\displaystyle \beta ={\frac {1}{4}}}$ und ${\displaystyle \gamma ={\frac {1}{2}}}$

Gleichung	Folgerung
${\displaystyle {\dot {x}}_{n+1}={\dot {x}}_{n}+\int _{t_{n}}^{t_{n+1}}{\ddot {x}}^{h}(\tau )\,\mathrm {d} \tau ={\dot {x}}_{n}+{\frac {\Delta t}{2}}({\ddot {x}}_{n}+{\ddot {x}}_{n+1})}$	${\displaystyle \gamma ={\frac {1}{2}}}$
${\displaystyle {\begin{aligned}x_{n+1}=&x_{n}+\int _{t_{n}}^{t_{n+1}}{\dot {x}}^{h}\,\mathrm {d} \tau \\=&x_{n}+\int _{t_{n}}^{t_{n+1}}\left({\dot {x}}_{n}+{\frac {\tau -t_{n}}{2}}({\ddot {x}}_{n}+{\ddot {x}}_{n+1})\right)\mathrm {d} \tau \\=&x_{n}+\Delta t\,{\dot {x}}_{n}+{\frac {\Delta t^{2}}{4}}({\ddot {x}}_{n}+{\ddot {x}}_{n+1})\end{aligned}}}$	${\displaystyle \beta ={\frac {1}{4}}}$

Zentrale Differenzenquotienten

→ Hauptartikel: Differenzenquotient

Die zentralen Differenzenquotienten

${\displaystyle {\dot {x}}_{n}={\frac {1}{2\Delta t}}(x_{n+1}-x_{n-1})}$		(VI)

${\displaystyle {\ddot {x}}_{n}={\frac {1}{\Delta t^{2}}}(x_{n+1}-2\,x_{n}+x_{n-1})}$		(VII)

entsprechen den obigen #Gleichungen (IV) und (V) mit

${\displaystyle \beta =0}$ und ${\displaystyle \gamma ={\frac {1}{2}}}$.

Gleichung	Folgerung
${\displaystyle {\begin{array}{lcl}{\ddot {x}}_{n}&=&{\frac {1}{\Delta t^{2}}}(x_{n+1}-2\,x_{n}+x_{n-1})\\[1ex]\rightarrow x_{n-1}&=&\Delta t^{2}{\ddot {x}}_{n}-x_{n+1}+2\,x_{n}\\[1ex]{\dot {x}}_{n}&=&{\frac {1}{2\Delta t}}(x_{n+1}-x_{n-1})\\[1ex]\rightarrow \Delta t{\dot {x}}_{n}&=&{\frac {1}{2}}(x_{n+1}-x_{n-1})=x_{n+1}-{\frac {\Delta t^{2}}{2}}{\ddot {x}}_{n}-x_{n}\\[1ex]\rightarrow x_{n+1}&=&x_{n}+\Delta t{\dot {x}}_{n}+{\frac {\Delta t^{2}}{2}}{\ddot {x}}_{n}\end{array}}}$	${\displaystyle \beta =0}$
${\displaystyle {\begin{array}{lcl}{\frac {\Delta t}{2}}{\ddot {x}}_{n}&=&{\frac {1}{2\Delta t}}(x_{n+1}-2\,x_{n}+x_{n-1})\\[1ex]{\frac {\Delta t}{2}}{\ddot {x}}_{n+1}&=&{\frac {1}{2\Delta t}}(x_{n+2}-2\,x_{n+1}+x_{n})\\[2ex]{\dot {x}}_{n}&=&{\frac {1}{2\Delta t}}(x_{n+1}-x_{n-1})\\[1ex]{\dot {x}}_{n+1}&=&{\frac {1}{2\Delta t}}(x_{n+2}-x_{n})\\[1ex]\rightarrow {\frac {\Delta t}{2}}({\ddot {x}}_{n}+{\ddot {x}}_{n+1})&=&{\frac {1}{2\Delta t}}(x_{n+2}-x_{n}-x_{n+1}+x_{n-1})={\dot {x}}_{n+1}-{\dot {x}}_{n}\\[1ex]\rightarrow {\dot {x}}_{n+1}&=&{\dot {x}}_{n}+{\frac {\Delta t}{2}}({\ddot {x}}_{n}+{\ddot {x}}_{n+1})\end{array}}}$	${\displaystyle \gamma ={\frac {1}{2}}}$

Explizite Zeitintegration

→ Hauptartikel: Explizites Euler-Verfahren

Das explizite Zeitintegrationsverfahren gehört nicht zur Familie der (impliziten!) Newmark-beta Algorithmen und wird hier nur zu Vergleichszwecken angegeben. Die obigen Gleichungen (VI) und (VII) für die zentralen Differenzenquotienten sind äquivalent zu

${\displaystyle {\begin{aligned}{\dot {x}}_{n+1/2}=&{\frac {1}{\Delta t}}(x_{n+1}-x_{n})\\{\ddot {x}}_{n}=&{\frac {1}{\Delta t}}({\dot {x}}_{n+1/2}-{\dot {x}}_{n-1/2})\end{aligned}}}$.

Hier fällt auf, dass die Geschwindigkeiten immer in der Mitte der Zeitintervalle berechnet werden. Mit der Annahme

${\displaystyle {\begin{array}{lcl}{\dot {x}}_{n+1}&\approx &{\dot {x}}_{n+1/2}={\dot {x}}_{n-1/2}+\Delta t\,{\ddot {x}}_{n}\\x_{n+1}&=&x_{n}+\Delta t\,{\dot {x}}_{n+1/2}=x_{n}+\Delta t\,({\dot {x}}_{n-1/2}+\Delta t\,{\ddot {x}}_{n})\end{array}}}$

können die Werte ${\displaystyle x_{n+1}}$ und die Geschwindigkeiten ${\displaystyle {\dot {x}}_{n+1}}$ zum Zeitpunkt t_n+1 auf bereits bekannte Ergebnisse ${\displaystyle x_{n},{\dot {x}}_{n-1/2},{\ddot {x}}_{n}}$ zurückgeführt werden und die Differenzial­gleichung liefert die Bestimmungsgleichung für die nunmehr einzige Unbekannte ${\displaystyle {\ddot {x}}_{n+1}}$.

Beispiel

Zeitintegration mit Algorithmen der Newmark Familie

Ein einfacher Schwinger gehorche in Abwesenheit einer Erregung der homogenen Differenzial­gleichung

${\displaystyle {\ddot {x}}(t)+2\,D\,{\dot {x}}(t)+x(t)=0}$.

mit Dämpfungsgrad D. Die analytische Lösung dieser Gleichung lautet

${\displaystyle {\begin{aligned}x(t)=\exp(-D\,t)[&A\sin(\omega t)+B\cos(\omega t)]\\{\dot {x}}(t)=\exp(-D\,t)[&(A\omega -B\,D)\cos(\omega t)\\&-(A\,D+B\omega )\sin(\omega t)]\end{aligned}}}$

mit ${\displaystyle \omega ={\sqrt {1-D^{2}}}}$, der Exponentialfunktion exp sowie dem Sinus und Cosinus sin bzw. cos. Mit den Parametern aus der Tabelle

Parameter	A	B	D	ω
Wert	13	84	13⁄85	84⁄85

resultieren die Anfangsauslenkung und -geschwindigkeit

${\displaystyle x(t=0)=x_{0}=84,\quad {\dot {x}}(t=0)={\dot {x}}_{0}=0}$

sowie eine Anfangsbeschleunigung

${\displaystyle {\ddot {x}}(t=0)={\ddot {x}}_{0}=-x_{0}-2\,D\,{\dot {x}}_{0}=-84}$

Die Differenzialgleichung liefert eine Gleichung für die aktuellen Unbekannten zur Zeit t_n+1:

${\displaystyle {\ddot {x}}_{n+1}+2\,D\,{\dot {x}}_{n+1}+x_{n+1}=0}$

Zusammen mit obigen Gleichungen (IV) und (V) resultiert ein geschlossenes System von drei Gleichungen für drei Unbekannte, das von

${\displaystyle {\begin{aligned}{\ddot {x}}_{n+1}=&{\frac {[2\,D\,(\gamma -1)+(\beta -1/2)\,\Delta t]\,{\ddot {x}}_{n}\,\Delta t-(2\,D+\Delta t)\,{\dot {x}}_{n}-x_{n}}{1+2\,D\,\gamma \,\Delta t+\beta \,\Delta t^{2}}}\\{\dot {x}}_{n+1}=&{\dot {x}}_{n}+[(1-\gamma )\,{\ddot {x}}_{n}+\gamma \,{\ddot {x}}_{n+1}]\Delta t\\x_{n+1}=&-{\ddot {x}}_{n+1}-2\,D\,{\dot {x}}_{n+1}\end{aligned}}}$

gelöst wird. Die Zeitintegration mit dem Newmark-Verfahren im Intervall t∈[0,10π] und 𝚫t=0,5 liefert die Verläufe im Bild. Die mittlere Abweichung über N Zeitschritte

${\displaystyle e={\sqrt {\frac {\sum _{n=1}^{N}{(x_{n}-x(t_{n}))}^{2}}{N}}}}$

als Maß für die Genauigkeit enthält die Tabelle mit vier signifikanten Stellen:

Verfahren	β	γ	Mittlere Abweichung ${\displaystyle e}$
			𝚫t=0,5, N=64	𝚫t=0,05, N=630,	𝚫t=0,005, N=6285
Lineare Beschleunigung	1/6	1/2	0,8780	8,839e-03	8,850e-05
Zentrale Differenzen	0	1/2	1,262	1,239e-02	1,241e-04
Konstante mittlere Beschleunigung	1/4	1/2	1,859	1,861e-02	1,863e-04
Explizite Integration	–	–	8,478	7,644e-01	7,577e-02

In diesem Beispiel ergibt eine zehntel so große Zeitschrittweite eine etwa hundertfache Genauigkeit bei Benutzung der Newmark-beta-Verfahren und eine etwa zehnfache Genauigkeit der expliziten Integration.

Literatur

• R. Gasch, K. Knothe, R. Liebich: Strukturdynamik. Springer Verlag, 2012, ISBN 978-3-540-88977-9.
• T. Belytschko, T.J.R. Hughes (Hrsg.): Computational methods for transient analysis. North-Holland, 1986, ISBN 978-0-444-86479-6.