BA Marvin Henke

Herleitung inkompressible,inviskose, wirbelfreie Strömung um einen Zylinder

Betrachtet wird ein sich mit Geschwindigkeit \(\vec{v_0}=v_0 \hat{x}\) durch ein Fluid bewegender Zylinder mit Radius \(R\). Es wird in Zylinderkoordinaten \((r,\varphi,z)\) gerechnet, wobei die \(z\)-Dimension irrelevant für die Rechnung ist. Aufgrund der Annahmen (inkompressibel, inviskos, wirbelfrei) gilt folgendes: \begin{align} \frac{\mathrm{d} \rho}{\mathrm{d} t} = 0\ ,\ \frac{\partial \rho}{\partial t} + \vec{\nabla}\cdot(\rho\vec{v}) = \frac{\mathrm{d} \rho}{\mathrm{d} t} + \rho (\vec{\nabla}\cdot\vec{v}) = 0 \ \Rightarrow\ \vec{\nabla}\cdot\vec{v}=0\\ \vec{\nabla}\times\vec{v}=0 \ \Rightarrow\ \exists \phi : \vec{v} = \vec{\nabla} \phi \end{align} Für das Potential \(\phi\) folgt Aufgrund von \(\vec{\nabla}\cdot\vec{v}=0\), dass \(\Delta\phi = 0\) gilt (Laplace-Gleichung).

Im Unendlichen soll die Geschwindigkeit des Fluids verschwinden, d.h. \( \vec{v} \xrightarrow[]{r \to \infty} \vec{0}\).

Für den Rand des Zylinders soll die Relativgeschwindigkeit des Fluids senkrecht zur Oberflächennormale sein, d.h. \(\forall \varphi : (\vec{v}(r=R,\varphi) - v_0 \hat{x})\cdot \hat{n} = 0 \).

Es ist also eine Lösung der Laplace-Gleichung auf \( \mathbb{R}^2 \setminus B_R(0) \) gesucht, welche im Unendlichen einen verschwindenden Gradienten hat und auf \( \partial B_R(0) \) die Neumann-Randbedingung \( \partial_n \phi = v_0 \hat{n}\cdot\hat{x} \) erfüllt. Die Fundamentallösungen der Laplace-Gleichung inspirieren folgenden Ansatz: \begin{align} \phi(r,\varphi) = \frac{c(\varphi)}{r} \end{align} Für beschränkte \( c(\varphi) \) verschwindet das Potential und die Geschwindigkeit im Unendlichen. Einsetzen in die Laplace-Gleichung liefert die folgende Bedingung an \( c(\varphi) \): \begin{align} c^{\prime\prime}(\varphi) + c(\varphi) = 0 \ \Rightarrow\ c(\varphi) = A \sin (\varphi) + B \cos (\varphi) \end{align}

Durch die Neumann-Bedingung lassen sich die Koeffizienten bestimmen: \(A = 0\) und \(B=-v_0 R^2\)

Damit ist das gesuchte Potential gefunden: \begin{align} \phi(r,\varphi) = -v_0 R^2 \frac{\cos \varphi}{r} \end{align} Es ergibt sich folgende Flussgeschwindigkeit als Gradient des Potentials: \begin{align} \vec{\nabla}\phi = \vec{v}(r,\varphi) = \frac{R^2}{r^2} v_0 \left[ \hat{r} \cos\varphi + \hat{\varphi} \sin\varphi \right] \end{align}

Für die Flussgeschwindigkeit um einen umströmten statischen Zylinder, wird nun die Geschwindigkeit des bewegten Zylinders subtrahiert. Es ergibt sich: \begin{align} \vec{v}(r,\varphi) = \frac{R^2}{r^2} v_0 \left[ \hat{r} \cos\varphi + \hat{\varphi} \sin\varphi \right] - v_0 \hat{x}\\ \vec{v}(r,\varphi) = \frac{R^2}{r^2} v_0 \left[ \hat{r} \cos\varphi + \hat{\varphi} \sin\varphi \right] - v_0 (\hat{r} \cos\varphi - \hat{\varphi} \sin\varphi)\\ \vec{v}(r,\varphi) = v_0 \left[ \hat{r} (\frac{R^2}{r^2} - 1) \cos\varphi + \hat{\varphi} (\frac{R^2}{r^2} + 1) \sin\varphi \right] \end{align}

Für den Druck lässt sich mit der Bernoulli-Gleichung \(\displaystyle \frac{v^2}{2} + \frac{p}{\rho} = \mathrm{const.}\) folgender Ausdruck herleiten: \begin{align} p=\frac{\rho}{2}v_0^2\left(2\frac{R^2}{r^2}\cos2\varphi-\frac{R^4}{r^4}\right) \end{align}

Numerik und Simulationsparameter

Die Simulation läuft in Zylinderkoordinaten, der Radius des Zylinders beträgt \(R = 1\), der Radius des Simulationsgebiets beträgt \(R_{\infty}=10\).

Die externe Einströmgeschwindigkeit wird über die Mach-Zahl \(\mathcal{M}=\frac{v_{\mathrm{ext}}}{c}\) geregelt, die Reynoldszahl \(\mathrm{Re} = \frac{v L}{\nu}\) regelt die kinematische Viskosität \(\nu\).

Die externe Dichte \(\rho_{\mathrm{ext}}\) wird auf \(1\) normiert, genauso die Schallgeschwindigkeit \(c\).

Mit \(p V^{\gamma} = \mathrm{const.}\) und \(c^2 = \left(\frac{\partial p}{\partial\rho}\right)_S\) folgt für den externen Druck \(p_{\mathrm{ext}} = \frac{c^2 \rho_{\mathrm{ext}}}{\gamma}\) , wobei \(\gamma\) der Isentropenexponent \(\gamma = 1 + \frac{2}{f}\) mit der Anzahl an Freiheitsgraden \(f\) ist.

Die Anfangsbedingungen sind folgende: \begin{align} \vec{v}_0 = -v_{\mathrm{ext}} \hat{x}\\ \rho_0 = \rho_{\mathrm{ext}} = 1\\ p_0 = p_{\mathrm{ext}} \end{align}

Für die Winkel sind periodische Randbedingungen definiert, in radiale Richtung werden für die Zylinderoberfläche reflektive Randbedingungen verwendet. In \(z\)-Richtung werden outflow-Randbedingungen benutzt.

In radiale Richtung auf dem äußeren Rand des Simulationsgebiets sind die Randbedingungen benutzerdefiniert. Standardmäßig wird hier \(\rho = \rho_{\mathrm{ext}}\), \(\vec{v} = \vec{v}_0\) und \(p = p_{\mathrm{ext}}\) für die Geisterzellen verwendet. Möglicherweise ist es sinnvoller für \(x < 0\) outflow/inflow Randbedingungen zu verwenden. Um dies zu testen sollte eine Simulation mit sehr großem \(R_{\infty} \approx 50\) gemacht werden und mit den beiden Simulationen für \(R_{\infty} = 10\) verglichen werden.

Bestimmung der Regimeübergänge

Im folgenden sollen die Reynoldszahlen der Übergänge verschiedener Strömungsregime bestimmt werden. Da die Strömungskenngrößen (\(z.B. \rho, \vec{v}\)) aufgrund der Numerik diskret sind, helfen folgende Definitionen für Skalarfelder: \begin{align} t = k\Delta t \ \ \Rightarrow s(\vec{x},t)\ \widehat{=}\ s_i^k \\ \dot{s}_i^{k+\frac{n}{2}}\ \widehat{=}\ \frac{s_i^{k+n} - s_i^k}{n\Delta t} \end{align} Außerdem sei \(\Delta t_k\) die \(k\)-te Zeitschrittweite und \(\Delta V_i\) das Volumen der \(i\)-ten Zelle. Im Folgenden wird über Zeitintervalle und Volumina summiert, dafür gilt \(T = \sum_k \Delta t_k\) und \( V = \sum_i \Delta V_i \).

Methodik zur Bestimmung des Übergangs von L1 in L2

Im Folgenden ist die implementierte Methodik zur Bestimmung des ersten Übergangs von L1 (creeping flow / non-separation regime [0 < RE < 4-5]) zu L2 (closed near-wake regime [4-5 < RE < 30-48]) beschrieben. Der wesentliche Unterschied liegt in der Symmetrie des Flusses bezüglich der y-Achse. Für L1 ist eine hohe Symmetrie zu erwarten, während für L2 aufgrund der sich bildenden Vortices hinter dem Zylinder eine niedrigere Symmetrie zu erwarten ist. Als Maß für die Asymmetrie eines Skalarfelds \(s\) (z.B. \(\rho\)) wird folgende Definition verwendet: \begin{align} A(s) := \frac{1}{2 T} \sum_k \frac{\Delta t_k}{\sum_i \Delta V_i \cdot (s_i^k)^2} \sum_{-\frac{\pi}{2} < \varphi < \frac{\pi}{2}} \Delta V(\varphi) \left[s^k(\varphi) - s^k(\pi - \varphi) \right]^2 \end{align} Diese Definition garantiert, dass für Skalarfelder \(s(x,y) = -s(-x,y)\) die Asymmetrie \(A(s) = 1\) beträgt. Aufgrund der Normierung hat \(A(s)\) immer die Dimension Zahl.

Methodik zur Bestimmung des Übergangs von L2 in L3

Beim Übergang von L2 (closed near-wake regime [4-5 < RE < 30-48]) zu L3 (periodic laminar regime [30-48 < RE < 180-200]) ändert sich das Verhalten des Systems im steady state. In L2 wird nach einiger Zeit ein zeitlich konstanter steady state angenommen, in L3 hingegen verhält sich das System periodisch, d.h. man müsste den Übergang anhand der Zeitableitung der Strömungskenngrößen erkennen können. Dazu wird folgende Größe definiert: \begin{align} \frac{D}{Dt}(s,n) := \frac{1}{V T} \sum_{i, k} \Delta t_k \Delta V_i {(\dot{s}_i^{k+\frac{n}{2}})}^2 \end{align} Prinzipiell handelt es sich hier um den Erwartungswert bezüglich Zeit und Ort der Größe \({(\dot{s}_i^{k+\frac{n}{2}})}^2\). Das \(n\) gibt an wie viele Mittelungen bei der Berechnung des diskreten Differenzenquotienten gemacht werden. Die Dimension von \(\frac{D}{Dt}(s,n)\) ist \(\left[\frac{s^2}{t^2}\right]\) und somit ist \(\frac{D}{Dt}(s,n)\) nicht einheitenlos. Eine mögliche Normierung wäre durch \(\langle s^2 \rangle_{\vec{x},t} \left( \frac{v_{\mathrm{ext}}}{R} \right)^2\)