v3.1 – Transportprobleme

 

3. Anfangs-Randwertaufgaben

Wir betrachten die Transportgleichung

\frac{{\partial u}}{{\partial t}}+a\frac{{\partial u}}{{\partial x}} = 0,\quad \quad a = const

Dies ist der lineare und skalare Spezialfall der nichtlinearen Erhaltungsgleichung

\frac{{\partial \vec u}}{{\partial t}}+\frac{{\partial \vec F\left( {\vec u} \right)}}{{\partial x}} = 0

mit der z.B. Strömungen von Gasen simuliert werden können.

3.1.1 Reine Anfangswertaufgabe (unendliches Gebiet)

3.1.1.1 Aufgabenstellung

Die Aufgabe

\frac{{\partial u}}{{\partial t}}+a\frac{{\partial u}}{{\partial x}} = 0\quad in\:\:\:\mathbb{R} \times \left( {0,\infty } \right)

u\left( {x,0} \right) = {u_0}\left( x \right)\quad in\:\:\:\mathbb{R}

besitzt die exakte Lösung

u\left( {x,t} \right) = {u_0}\left( {x-at} \right)

num-301-transportproblem-charakteristiken-losung-konstant

Die Geraden x-at = const heißen Charakteristiken. Die Lösung ist entlang der Charakteristiken konstant. Informationen ({u_0}\left( x \right)) pflanzen sich entlang der Charakteristiken fort.

Die Lösung erfüllt ein Maximumprinzip:

\max \limits_{x \in \mathbb{R},t \in \left( {0,\infty } \right)} u\left( {x,t} \right) =  \max \limits_{x \in \mathbb{R}} {u_0}\left( x \right)

Dieses sollte möglichst auch im Diskreten gelten.

3.1.1.2 Diskretisierung

Sei h die Ortsschrittweite und \tau die Zeitschrittweite. Sei weiter durch

{x_j} = jh,\quad j \in \mathbb{Z}

{t^k} = k\tau ,\quad k \in \mathbb{N}

das Orts- bzw. Zeitgitter definiert und bezeichne {u_{h,\tau }} die gesuchte Näherungslösung. Wir vereinbaren

u_j^k: = {u_{h,\tau }}\left( {{x_j},{t^k}} \right)

Bei unendlichem Gebiet kommen nur explizite Verfahren in Betracht:

FTBS: \frac{1}{\tau }\left( {u_j^{k+1}-u_j^k} \right)+\frac{a}{h}\left( {u_j^k-u_{j-1}^k} \right) = 0

FTFS: \frac{1}{\tau }\left( {u_j^{k+1}-u_j^k} \right)+\frac{a}{h}\left( {u_{j+1}^k-u_j^k} \right) = 0

FTCS: \frac{1}{\tau }\left( {u_j^{k+1}-u_j^k} \right)+\frac{a}{{2h}}\left( {u_{j+1}^k-u_{j-1}^k} \right) = 0

Mit \gamma : = \frac{{a\tau }}{h} sind diese Gleichungen äquivalent zu

FTBS: u_j^{k+1} = \left( {1-\gamma } \right)u_j^k+\gamma u_{j-1}^k

FTFS: u_j^{k+1} = \left( {1+\gamma } \right)u_j^k-\gamma u_{j-1}^k

FTCS: u_j^{k+1} = u_j^k+\frac{\gamma }{2}\left( {u_{j-1}^k-u_{j+1}^k} \right)

Dabei ist u_j^0 = {u_0}\left( {{x_j}} \right). Die Abkürzungen stehen für:

  • FTBS: Forward in time, backward in space
  • FTFS: Forward in time, forward in space
  • FTCS: Forward in time, central in space

3.1.1.3 Fehleranalyse

Den Konsistenzfehler erhält man wieder, indem man die exakte Lösung in die Differenzengleichung einsetzt.

Für das FTBS-Verfahren gilt zum Beispiel:

{\ell _{h,\tau }}\left( {x,t} \right) = \frac{{u\left( {x,t+\tau } \right)-u\left( {x,t} \right)}}{\tau }+\frac{a}{h}\left( {u\left( {x,t} \right)-u\left( {x-h,t} \right)} \right)

Die Diskretisierung ist stabil, wenn das diskrete Maximumprinzip für beliebige Anfangsbedingungen {u_0}\left( x \right) erfüllt wird:

\max \limits_{j \in \mathbb{Z}} \left| {u_j^{k+1}} \right| \leq  \max \limits_{j \in \mathbb{Z}} \left| {u_j^k} \right|\quad \quad \quad \quad \left( 1 \right)

Es gilt wie im Fall von Randwertaufgaben: Ist die Diskretisierung stabil und konsistent mit der Ordnung p, dann ist sie auch konvergent mit der Ordnung p.

3.1.1.4 Konsistenz für FTBS

Wie oben erwähnt gilt:

{\ell _{h,\tau }}\left( {x,t} \right) = \frac{{u\left( {x,t+\tau } \right)-u\left( {x,t} \right)}}{\tau }+\frac{a}{h}\left( {u\left( {x,t} \right)-u\left( {x-h,t} \right)} \right)

Mit Hilfe der Taylorentwicklung erhalten wir daraus:

{\ell _{h,\tau }} = \frac{1}{\tau }\left[ {u+\tau {u_t}+\frac{{{\tau ^2}}}{2}{u_{tt}}+O\left( {{\tau ^3}} \right)-u} \right]+\frac{a}{h}\left[ {u-\left( {u-h{u_x}+\frac{{{h^2}}}{2}{u_{xx}}+O\left( {{h^3}} \right)} \right)} \right]

\quad = \left[ {{u_t}+\frac{\tau }{2}{u_{tt}}+O\left( {{\tau ^2}} \right)} \right]+a\left[ {{u_x}-\frac{h}{2}{u_{xx}}+O\left( {{h^2}} \right)} \right]

\quad = \underbrace {{u_t}+a{u_x}}_{ = 0}+\frac{\tau }{2}{u_{tt}}-\frac{{ah}}{2}{u_{xx}}+O\left( {{h^2}+{\tau ^2}} \right)

\quad = \frac{1}{2}\left( {\tau {u_{tt}}-ah{u_{xx}}} \right)+O\left( {{h^2}+{\tau ^2}} \right)

Im Allgemeinen ist die Konsistenz also O\left( {h+\tau } \right). Speziell für \tau = \frac{h}{a} gilt sogar eine bessere Ordnung, denn mit

{u_{tt}} = -{\left( {a{u_x}} \right)_t} = -a{\left( {{u_t}} \right)_x} = {a^2}{u_{xx}}

folgt

\tau {u_{tt}}-ah{u_{xx}} = \tau {a^2}{u_{xx}}-ah{u_{xx}} = \left( {\tau {a^2}-ah} \right){u_{xx}}

3.1.1.5 Stabilität für FTBS

Wir zeigen die Beziehung \left( 1 \right) wie folgt:

\max \limits_j \left| {u_j^{k+1}} \right| =  \max \limits_j \left| {\left( {1-\gamma } \right)u_j^k+\gamma u_{j-1}^k} \right|

\quad \leq \left| {1-\gamma } \right| \max \limits_j \left| {u_j^k} \right|+\left| \gamma \right| \max \limits_j \left| {u_{j-1}^k} \right|

\quad =  \max \limits_j \left| {u_j^k} \right|

wenn 0 \leq \gamma \leq 1, d.h. wenn \alpha > 0 und

\tau \leq \frac{h}{a}

Diese Bedingung heißt CFL-Bedingung nach Courant, Friedrichs und Levy. Sie bedeutet, dass Informationen nur um maximal eine Zelle pro Zeitschritt transportiert werden können.

num-301-charakteristik-cfl-bedingung

Die Gerade x-at = const = {x_j}-a{t^{k+1}} schneidet die Linie t = {t^k} zwischen {x_{j-1}} und {x_j}, wenn

x = {x_j}+a\left( {{t^k}-{t^{k+1}}} \right) = {x_j}-\underbrace {a\tau }_{ \leq h} \geq {x_{j-1}}

Dass die CFL-Bedingung notwendig für die Stabilität des Verfahrens ist, kann man an einem Beispiel verdeutlichen. Sei

u_j^k = 1 für j < {j_0} und u_j^k = 0 für j \geq {j_0}.

Dann ist

u_{j0}^{k+1} = \left( {1-\gamma } \right)\underbrace {u_{j0}^k}_{ = 0}+\gamma \underbrace {u_{j0-1}^k}_{ = 1} = \gamma

Daraus folgt \max \limits_{j \in \mathbb{Z}} \left| {u_j^{k+1}} \right| >  \max \limits_{j \in \mathbb{Z}} \left| {u_j^k} \right|, wenn \gamma > 1 ist.

Stabilität der anderen Verfahren:
FTFS ist stabil, wenn a < 0 und \tau \leq h/\left| a \right|.
FTCS ist nicht stabil.

Lax-Friedrichs-Verfahren:

\frac{1}{\tau }\left( {u_j^{k+1}-\frac{{u_{j-1}^k+u_{j+1}^k}}{2}} \right)+\frac{a}{{2h}}\left( {u_{j+1}^k-u_{j-1}^k} \right) = 0

Mit \gamma : = a\tau /h ist diese Gleichung äquivalent zu

u_j^{k+1} = \frac{1}{2}\left( {\left( {1+\gamma } \right)u_{j-1}^k+\left( {1-\gamma } \right)u_{j+1}^k} \right)

Dieses Verfahren kann auch als Mittelung von FTBS und FTFS interpretiert werden. Das Lax-Friedrichs-Verfahren ist unabhängig vom Vorzeichen von a stabil, wenn \tau \leq h/\left| a \right|.

3.1.1.6 Numerische Diffusion

Die Analyse des Konsistenzfehlers zeigt, dass das FTBS-Schema die Differentialgleichung

{u_t}+a{u_x}-\frac{{a\left( {h-a\tau } \right)}}{2}{u_{xx}} = 0

mit der Ordnung 2 approximiert. Man kann das so interpretieren, dass das Schema einen Diffusionsterm in die Gleichung hineinträgt. Man spricht dann von numerischer Diffusion, die ein „Verschmieren“ der Näherungslösung bewirkt.

Wenn die CFL-Bedingung erfüllt ist, hat der Koeffizient vor {u_{xx}} das richtige Vorzeichen. Ist die CFL-Bedingung nicht erfüllt, ergibt sich eine nicht korrekt gestellte Differentialgleichung.

3.1.1.7 Beispiel Hütchenfunktion

Zur Illustration berechnen wir die numerische Lösung im Fall a = 1 mit den verschiedenen Verfahren und verschiedenen Schrittweiten. Als Anfangsbedingung verwenden wir die Funktion

{u_0}\left( x \right) = \left\{ {\begin{array}{*{20}{c}}{x-1} & {x \in \left[ {1,2} \right]} \\{3-x} & {x \in \left[ {2,3} \right]} \\  0 & {sonst} \\   \end{array} } \right.

num-301-transportproblem-hutchenfunktion-diffusion-beispiel

3.1.1.8 Zusammenfassung

  • Das FTBS-Verfahren konvergiert mit der Ordnung p = 1, der Fehler ist von der Ordnung O\left( {h+\tau } \right), wenn \tau \leq h/a und a > 0
  • Das FTFS-Verfahren konvergiert mit der Ordnung p = 1, der Fehler ist von der Ordnung O\left( {h+\tau } \right), wenn \tau \leq h/\left| a \right| und a < 0
  • Das FTCS-Verfahren hat zwar einen Konsistenzfehler von O\left( {{h^2}+\tau } \right), ist aber instabil.
  • Das Lax-Friedrichs-Verfahren konvergiert mit der Ordnung p = 1, d.h. der Fehler ist von der Ordnung O\left( {h+\tau } \right), wenn \tau \leq h/\left| a \right|. Dem Vorteil der Unabhängigkeit vom Vorzeichen von a steht der Nachteil der größeren numerischen Diffusion gegenüber

3.1.2 Transportgleichung, endliches Gebiet

3.1.2.1 Aufgabenstellung

Wir betrachten die Transportgleichung

{u_t}+a{u_x} = 0\quad in\:\:\:\left( {0,1} \right) \times \left( {0,\infty } \right),\quad a = const > 0

u\left( {x,0} \right) = {u_0}\quad in\:\:\:\left( {0,1} \right)

Welche Randbedingungen müssen gestellt werden?

num-301-randbedingung-endliches-gebiet-transportproblem

Auf den Charakteristiken, die {u_0} transportieren, brauchen wir keine Randbedingung, da wir dort schon die Anfangsbedingung haben.
Auf allen anderen Charakteristiken brauchen wir nur an einer Seite eine Randbedingung (da die Funktion auf der Charakteristik konstant ist). Sinnvollerweise wählt man Randbedingungen am Einflussrand, also hier links.

Die exakte Lösung lautet dann:

u\left( {x,t} \right) = \left\{ {\begin{array}{*{20}{c}}{{u_0}\left( {x-at} \right)} & {x > at} \\{{g_0}\left( {t-\frac{x}{a}} \right)} & {x < at} \\   \end{array} } \right.

Wenn die Lösung stetig sein soll, müssen {u_0} und {g_0} kompatibel sein: {u_0}\left( 0 \right) = {g_0}\left( 0 \right). Wenn a < 0 gilt, dann muss eine Randbedingung am rechten Rand vorgegeben werden.

Das Maximumprinzip wird modifiziert zu

\max \limits_{x \in \left( {0,1} \right),t \in \left( {0,{t_*}} \right)} \left| {u\left( {x,t} \right)} \right| \leq \max \left\{ { \max \limits_{x \in \left( {0,1} \right)} \left| {{u_0}\left( x \right)} \right|, \max \limits_{t \in \left( {0,{t_*}} \right)} \left| {{g_0}\left( t \right)} \right|} \right\}

3.1.2.2 Diskretisierung

Sei h = \frac{1}{N} die Ortsschrittweite. Wie bei den Randwertaufgaben ersetzen wir das Gebiet \Omega = \left( {0,1} \right) durch das Ortsgitter {\Omega _h} = \left\{ {{x_i} = ih,\:\:i = 1, \ldots ;N-1} \right\}. Wie bei den Anfangswertaufgaben definieren wir ein Zeitgitter {\omega _\tau } = \left\{ {{t^k} = k\tau ,\:\:k \in \mathbb{N}} \right\}. Wie im Fall des unendlichen Gebiets suchen wir eine Näherungslösung {u_{h,\tau }} und vereinbaren u_j^k: = {u_{h,\tau }}\left( {{x_j},{t^k}} \right). Konsistenzfehler und Diskretisierungsfehler werden analog zum Fall des unendlichen Gebiets definiert.

Eine Näherungslösung erfüllt das diskrete Maximumprinzip, wenn

\max \limits_{j = 0, \ldots ,N} \left| {u_j^{k+1}} \right| \leq \max \left\{ {{g_0}\left( {{t^{k+1}}} \right), \max \limits_{j = 0, \ldots ,N} \left| {u_j^k} \right|} \right\}

FTBS: Diese Diskretisierung ist wie im Fall des unendlichen Gebiets konsistent mit Ordnung 1 und erfüllt das Maximumprinzip für \tau \leq h/a.

BTBS: Für a > 0 wird diese Diskretisierung durch

\frac{1}{\tau }\left( {u_j^{k+1}-u_j^k} \right)+\frac{a}{h}\left( {u_j^{k+1}-u_{j-1}^{k+1}} \right) = 0,\quad j = 1, \ldots ,N

u_0^{k+1} = {g_0}\left( {{t^{k+1}}} \right)

beschrieben und kann rekursiv berechnet werden:

u_j^{k+1} = {\left( {\frac{1}{\tau }+\frac{a}{h}} \right)^{-1}}\left( {\frac{1}{\tau }u_j^k+\frac{a}{h}u_{j-1}^{k+1}} \right)\quad \quad \quad \quad \left( 2 \right)

Aus \left( 2 \right) folgt

\min \limits_{i = 0, \ldots ,j} \left( {u_i^k,u_0^{k+1}} \right) \leq u_j^{k+1} \leq  \max \limits_{i = 0, \ldots ,j} \left( {u_i^k,u_0^{k+1}} \right)

Das heißt, es gilt das diskrete Maximumprinzip ohne Bedingung an \tau. Das Verfahren ist unbedingt stabil.

Zur Illustration berechnen wir die numerische Lösung mit verschiedenen Schrittweiten:

num-203-btbs-stabilitat-beispiel-transportproblem

3.1.3 Transportgleichung, variabler Koeffizient

3.1.3.1 Randbedingung

Im vorherigen Abschnitt haben wir gesehen, dass bei positiver Transportgeschwindigkeit a die Randbedingung am linken Rand, x = 0, und bei negativer Transportgeschwindigkeit die Randbedingung am rechten Rand, x = 1, vorgegeben werden muss. Diese Eigenschaft bedeutet, übertragen auf den Fall einer variablen Geschwindigkeit a = a\left( {x,t} \right): Eine Randbedingung u\left( {0,t} \right) muss für diejenigen t vorgegeben werden, für die a\left( {0,t} \right) > 0 ist, und eine Randbedingung u\left( {1,t} \right) muss für diejenigen t vorgegeben werden, für die a\left( {1,t} \right) < 0 ist.

Insbesondere gibt es für a\left( {x,t} \right) = 1-2x zwei Einströmränder und für a\left( {x,t} \right) = 2x-1 zwei Ausströmränder:

num-301-transportgleichung-variabler-koeffizient

3.1.3.2 Charakteristiken

Auch im Fall variabler Transportgeschwindigkeit pflanzt sich die Information entlang von Charakteristiken fort. Die Charakteristik x = X\left( {t;\xi ,\tau } \right) durch den Punkt \left( {\xi ,\tau } \right) wird durch die Anfangswertaufgabe

\frac{{dX}}{{dt}} = a\left( {x,t} \right),\quad \quad X\left( {\tau ;\xi ,\tau } \right) = \xi

beschrieben.

3.1.3.3 Upwind-Idee

Als numerisches Verfahren eignet sich eine Kombination aus Vorwärts- und Rückwärtsdifferenz. Das FT-Verfahren lautet dann

u_j^{k+1} = \left\{ {\begin{array}{*{20}{c}}{\left( {1-\gamma _j^k} \right)u_j^k+\gamma _j^ku_{j-1}^k} & {a\left( {{x_j},{t^k}} \right) > 0} \\{\left( {1+\gamma _j^k} \right)u_j^k-\gamma _j^ku_{j-1}^k} & {a\left( {{x_j},{t^k}} \right) < 0} \\   \end{array} } \right.

mit \gamma _j^k = a\left( {{x_j},{t^k}} \right) \cdot \tau /h und der Stabilitätsbedingung \left| {\gamma _j^k} \right| \leq 1. Analog kann das BT-Verfahren als Kombination zweier stabiler Verfahren formuliert werden.

3.1.3.4 Burgers-Gleichung

Die Burgers-Gleichung ist eine nichtlineare Erhaltungsgleichung der Form

\frac{{\partial u}}{{\partial t}}+\frac{{\partial f\left( u \right)}}{{\partial x}} = 0,\quad \quad f\left( u \right) = \frac{1}{2}{u^2}

Als Anfangs-(Randwert-)Aufgabe kann also zum Beispiel

\frac{{\partial u}}{{\partial t}}+u\frac{{\partial u}}{{\partial x}} = 0\quad \quad in\:\:\:\Omega \times \left( {0,T} \right)

u\left( {x,0} \right) = {u_0}\left( x \right)\quad \quad in\:\:\:\Omega

formuliert werden. Je nach Gebiet und Entwicklung der Lösung benötigt man noch Randbedingungen.

Ähnliche Artikel

Kommentar verfassen