|
Programmierpraktikum NPDGL I
|
Unter einem Gitter \({\cal T}_h = \{e_i\}_{i=0}^H\) verstehen wir eine Zerlegung eines Gebietes \( \Omega \subset \mathbb{R}^d \) in \(H\) Gitterzellen \(e_i\), so dass
Hat der Schnitt zweier Gitterzellen positives \(d-1\)-dimensionales Hausdorff-Maß
\[ H^{d-1}(\bar{e_i} \cap \bar{e_j}) > 0, \]
so bezeichnen wir diesen Schnitt \(s_{ij} = \bar{e_i} \cap \bar{e_j}\) als Interface.
Wir wollen skalare Funktionen \( f:\mathbb{R} \to \mathbb{R} \) diskret beschreiben. Allgemein wollen wir dazu die Funktion als Linearkombination von lokalen Gitterfunktionen darstellen, die nur auf Gitterzellen definiert sind, d.h.
\[ f_h(x) = \sum_{i=1}^{H} \sum_{j=i}^{l(i)} F_{i,j} \varphi_{i,j}(x) \]
mit skalaren Größen \(F_{i,j} \in \mathbb{R}\), den sogenannten Freiheitsgraden (engl. degrees of freedom, Abkürzung: Dofs), und lokalen Basisfunktionen \(\varphi_{i,j}: \bar{e_i} \to \mathbb{R}\) für \(j=1\ldots,l(i)\) bei \(i=1,\ldots,H\). Die Funktion \(l:\mathbb{N} \to \mathbb{N}\) gibt die Anzahl der Freiheitsgrade auf einer Gitterzelle an.
In Abhängigkeit von der Wahl des zugrunde liegenden Gitters \({\cal T}_h\) und der Basisfunktionen erhält man somit einen Funktionenraum \(V_h\). Dieser Funktionenraum sollte die Lösungen des Problems möglichst gut approximieren.
Das einfachste Beispiel für einen diskreten Funktionenraum sind Finite-Volumen Funktionen. Hier gibt es jeweils genau eine Basisfunktion für jede Gitterzelle, d.h. \(l(i)=1\) für alle \(i=1,\ldots,H\) und diese ist dort konstant: \(\varphi_{i,1} = 1\) für alle \(i=1,\ldots,H\). Der Funktionenraum ist also unstetig.
Eine Projektion \( {\cal P}_h: L^2(\Omega) \to V_h^{\text{FV}}\) auf diesen Funktionenraum \(V_h^{\text{FV}}\) lässt sich sinnvollerweise durch
realisieren. Die Projektionsfunktion ist dann durch
\[ {\cal P}_h[f] = \sum_{i=1}^H F_{i,1} \varphi_{i,1}(x) \]
gegeben.
Zuerst wollen wir eine skalare Evolutionsgleichung in 1D lösen. Gesucht ist also eine Funktion \(u:\Omega\times[0,T] \to \mathbb{R}\) definiert auf einem Intervall \(\Omega:=[a,b]\), die die Gleichungen
\begin{align} \partial_t u(x,t) + \partial_x f(u(x,t)) &= 0 & \text{in }&\Omega\times[0,T] \\ u(x,0) &= u_0(x) &\text{in }&\Omega \\ u(x,t) &= u_{dir}(x,t) &\text{auf }&\Gamma_{in} \times [0,T] \end{align}
erfüllt. Die sogenannten Dirichlet-Randbedingungen gelten nur auf dem Einflußrand \(\Gamma_{in}:=\{x \in \partial \Omega | f'(u) \mathbf{n} < 0 \} \). Auf dem übrigen Rand \(\partial \Omega \backslash \Gamma_{in}\) wollen wir entweder No-Flow-Bedingungen ( \( \partial_x f(u) = 0 \)) oder Out-Flow-Bedingungen vorschreiben. Mehr dazu später.
Das numerische Verfahren zum Lösen der obigen Evolutionsgleichung mittels eines Finite Volumen Verfahrens in 1D besteht aus den folgenden drei Schritten:
\[ u_h^0 = {\cal P}_h[u_0] \]
(c.f. Diskrete Funktionenräume (Ortsdiskretisierung))\[ {\cal L}_h[u_h] \approx \partial_x f(u_h) \]
(c.f. Diskreter Finite Volumen Operator in 1D)\[ u_h^{k+1} = u_h^{k} - \Delta t^k {\cal L}_h[u_h^k] \]
(c.f. Zeitdiskretisierung mit Euler-Verfahren)In diesem Abschnitt wollen wir einen diskreten Operator definieren, der für 1D Finite Volumen Funktionen \( \partial_{x} f(u_h) \) löst.
Hierzu machen wir auf einer Zelle (einem Intervall) \(e_i = [x_{i,1}, x_{i,2}]\) den Finite Volumen Ansatz, d.h. wir mitteln über diese Zelle und wenden anschließend den Hauptsatz der Differential- und Integralrechnung an.
\[ \frac{1}{|e_i|} \int_{e_i} \partial_{x} f(u(x)) dx = \frac{1}{|e_i|} \left( f(u(x_{i,2})) - f(u(x_{i,1})) \right) \]
Weil Finite Volumen Funktionen unstetig sind, ist der Funktionswert \(f(u(x_{i,\cdot}))\) auf den Kantenelementen (Knotenpunkten) jedoch nicht definiert, weswegen wir in der Anwendung auf diskrete Funktionen \(u_h \in V_h^{\text{FV}}\) numerische Flüsse verwenden, d.h.
\[ \frac{1}{|e_i|} \int_{e_i} \partial_{x} f(u_h(x)) dx \approx \frac{1}{|e_i|} \left( g(U_{i},U_{i+1}) - g(U_{i-1}, U_{i}) \right) \]
mit numerischen Flüssen, die jeweils von den Funktionswerten auf den Gitterzellen, die dem Knotenpunkt \(x_{i,1}\) bzw. \(x_{i,2}\) benachbart sind abhängen.
Beispiele für Flüsse sind
\[ g(u,v):= \frac{1}{2} (f(u) + f(v)) \]
\[ g(u,v):= \frac{1}{2} (f(u) + f(v)) + \frac{\lambda}{2} (u-v), \qquad \lambda = \frac{\Delta x}{\Delta t} \]
\begin{align*} f^{+}(u) &= f(0) + \int_{0}^u \max(f'(s), 0) ds \\ f^{-}(u) &= \hphantom{f(0) +} \int_{0}^u \min(f'(s), 0) ds \end{align*}
\[ g(u,v) := f^{+}(u) + f^{-}(v) \]
Bei einer Zelle am Rand des Gebiets Omega, also \( e_i \cap \partial \Omega \neq \emptyset \) wird eigentlich ein Freiheitsgrad \(U_0\) oder \(U_{H+1}\) benötigt.
Bei der Behandlung der Randwerte unterscheiden wir die folgenden Fälle:
Zur Diskretisierung des Zeitintervalls \([0,T]\) zerlegen wir dieses in \(K\) Zeitschritte der Länge \(\Delta t^k>0, k=1,\ldots,K\), so dass wir \(K+1\) Zeitpunkte \(t^0=0 \leq t^1=t^0+\Delta t^1 \leq \ldots \leq t^k=t^{k-1} + \Delta t^k = T \) erhalten, und aus (1) eine explizite semi-diskrete Differentialgleichung machen können:
\[ \frac{1}{\Delta t^k} \left( u(x,t^k) - u(x,t^{k+1}) \right) + \partial_x f(u(x,t^{k})) = 0 \]
Im zweiten Term könnte alternativ auch der Zeitschritt \(t^{k+1}\) eingesetzt werden. Dann wäre die Gleichung implizit.
Mit der Ortsdiskretisierung \(u_h^k(\cdot) \approx u(\cdot,t^k)\) für alle \(k=0,\ldots,K \) lautet die vollständige explizite Diskretisierung der Gleichung dann:
\[ u_h^{k+1} = u_h^k - \Delta t^k {\cal L}_h[u_h^k] \]
für alle \(k=0, \ldots, K\) wobei \(\Delta t^k\) so klein gewählt werden muss, dass es die CFL Bedingung
\[ \max_{x\in[a,b]}\frac{f'(u_h^k(x)) \Delta t^k}{\Delta x} \leq \frac{1}{2} \]
erfüllt.
1.7.3