BA Emilio Schmidt

(Der Echtzeit-Previewer rendert leider keine Formeln. Am besten stets mit "edit source" statt dem visual "edit" arbeiten und bitte ausgiebig den Knopf "Show preview" benutzen anstatt zig mal abzuspeichern. Diskussionen sollten unter Discussion geführt werden. -- Lothar)

Stokes-Fluss um eine Kugel

Kraftberechnung

\(\newcommand{\et}{ {\large\mathbb 1} }\) \begin{gather} \vec F=\iint\limits_{\partial K} \underline{\underline\sigma}\cdot\text d\vec S\\ \text{und}\quad F_z = \vec e_z\cdot\iint\limits_{\partial K}\underline{\underline\sigma}\cdot\text d\vec S=\iint\limits_{\partial K}\vec e_z\cdot\underline{\underline\sigma}\cdot\text d\vec S \stackrel{\text{Kugel}}{=}\iint\limits_{\partial K}\vec e_z\cdot\underline{\underline\sigma}\cdot\vec e_r\,\text dS \end{gather} Bei Zylindersymmetrie und Verwendung von KK wird daraus \begin{align} F_z &= \int_0^\pi\int_0^{2\pi}(\vec e_r\cos\theta-\vec e_\theta\sin\theta)\cdot\underline{\underline\sigma}(r{=}R,\theta)\cdot\vec e_r\,R^2\sin\theta\,\text d\phi\,\text d\theta\\ &=2\pi R^2\int_0^\pi \Big(\sigma_{rr}(r{=}R,\theta)\cos\theta\sin\theta-\sigma_{\theta r}(r{=}R,\theta)\sin^2\theta\Big)\,\text d\theta \end{align}

Dabei lautet der Spannungstensor (s. auch Wikipedia) \begin{equation} \underline{\underline\sigma} = -\et(p - \zeta\vec\nabla\cdot\vec u) + \mu\Big(\vec\nabla\vec u + (\vec\nabla\vec u)^\text{T}-\et\frac{2}{3}\nabla\cdot\vec u\Big) \end{equation} und bei der Formulierung von \(\vec\nabla\vec u\) und \(\vec\nabla\cdot\vec u\) in KK hilft ebenfalls die Wikipedia. Die Komponenten des Einheitstensors \(\et\) schließlich bilden in jeder Basis die Einheitsmatrix.

Kraftberechnung (Emilio)

Die Kraft, welche auf die umströmte Kugel wirkt, ist gegeben durch folgendes Oberflächenintegral: \(\newcommand{\sigt}{\underline{\underline\sigma}}\) \begin{align} \vec{F}=\oint_{\partial V} \sigt\cdot\text d\vec A \end{align}

Dabei gilt für ein Oberflächenelement \(\text d\vec{A}\) einer Kugel mit Radius \(r=R\), in Kugelkoordinaten: \begin{align} \text d\vec{A}=R^2\sin(\theta)\text d\varphi\text d\theta\vec{e}_r \end{align}

Einsetzen und weitere Berechnung in Kugelkoordinaten liefert: \begin{align} \vec{F}=\int_{0}^{\pi}\int_{0}^{2\pi}\sigt(r=R,\theta,\varphi)\cdot\vec{e}_r R^2\sin(\theta)\text d\varphi\text d\theta \end{align}

Auswerten des Matrix-Vektor-Produktes liefert, dass nur die \(rr\)-, \(\theta r\)- und \(\varphi r\)-Komponenten des Spannungstensors multipliziert mit entsprechenden Basisvektoren übrig bleiben. Es gilt also: \begin{align} \vec{F}=\int_{0}^{\pi}\int_{0}^{2\pi}[\sigma_{rr}\vec{e}_r+\sigma_{\theta r}\vec{e}_{\theta} +\sigma_{\varphi r}\vec{e}_{\varphi}]R^2\sin(\theta)\text d\varphi\text d\theta \end{align}

Mit der Annahme, dass vorliegendes Problem Zylindersymmetrie hat, lässt sich folgern, dass der Spannungstensor und damit auch dessen Komponenten unabhängig von \(\varphi\) ist, weshalb gelten muss: \begin{align} \sigma_{\varphi r}=0 \end{align}

Damit kann nun das Integral über \(\varphi\) berechnet werden: \begin{align} \vec{F}=2\pi R^2\int_{0}^{\pi}[\sigma_{rr}\vec{e}_r+ \sigma_{\theta r}\vec{e}_{\theta}]\sin(\theta)\text d\theta \end{align}

Um nun das Integral über \(\theta\) berechnen zu können, wird die explizite Form des Spannungstensors, für dieses Problem benötigt. Im Allgemeinen ist der Spannungstensor für reale, newtonsche Fluide gegeben durch: \begin{align} \sigt= -\bigg[p-\zeta(\vec{\nabla}\cdot\vec{u})\bigg]\mathbb{1} +\mu\bigg[\vec{\nabla}\vec{u}+(\vec{\nabla}\vec{u})^T -\frac{2}{3}(\vec{\nabla}\cdot\vec{u})\mathbb{1}\bigg] \end{align}

Dabei ist \(\zeta\) die Volumenviskosität, welche jedoch aufgrund (Begründung), vernachlässigt werden kann. Damit folgt: \begin{align} \sigt= -p\mathbb{1} +\mu\bigg[\vec{\nabla}\vec{u}+(\vec{\nabla}\vec{u})^T \frac{2}{3}(\vec{\nabla}\cdot\vec{u})\mathbb{1}\bigg] \end{align}

Dabei ist \(\mathbb{1}\) der Einheitstensor, \(\mu\) die dynamische Viskosität und \(\vec{u}\) die Geschwindigkeit des Fluides. Aufgrund dessen, dass bei der Kraftberechnung in Kugelkoordinaten nur die \(rr\)- und \(\theta r\)-Komponente des Spannungstensors eingehen, werden im Folgenden auch nur diese bestimmt. Die \(rr\)-Komponente ist dabei gegeben durch: \begin{align} \sigma_{rr}=-p+\mu\bigg[2\vec{\nabla}\vec{u}_{rr} -\frac{2}{3}\vec{\nabla}\cdot\vec{u}\bigg] \end{align}

Hierbei ist \(\vec{\nabla}\vec{u}_{rr}\) die \(rr\)-Komponente des Vektorgradienten \(\vec{\nabla}\vec{u}\), welche gegeben ist durch: \begin{align} \vec{\nabla}\vec{u}_{rr}=\frac{\partial u_r}{\partial r} \end{align}

Die Divergenz, des Geschwindigkeitsfeldes, ist bei Zylindersymmetrie gegeben durch: \begin{align} \vec{\nabla}\cdot\vec{u}= \frac{1}{r^2}\frac{\partial r^2 u_r}{\partial r} +\frac{1}{r\sin(\theta)} \frac{\partial \sin(\theta)u_{\theta}}{\partial\theta} \end{align}

Einsetzen in die \(rr\)-Komponente des Vektorgradienten liefert: \begin{align} \sigma_{rr}=-p +\mu\bigg[2\frac{\partial u_r}{\partial r}-\frac{2}{3}\bigg(\frac{1}{r^2}\frac{\partial r^2 u_r}{\partial r} +\frac{1}{r\sin(\theta)}\frac{\partial \sin(\theta)u_{\theta}}{\partial\theta}\bigg)\bigg] \end{align}

Aufgrund dessen, dass bei der Kraftberechnung das Oberflächenintegral über eine Kugel mit dem Radius \(r=R\) bestimmt wird, muss zusätzlich gelten: \begin{align} \sigma_{rr}=-p+\mu\bigg[2\frac{\partial u_r}{\partial r}\Big|_{r=R}-\frac{2}{3}\bigg(\frac{1}{R^2} \frac{\partial r^2 u_r}{\partial r}\Big|_{r=R}+\frac{1}{R\sin(\theta)} \frac{\partial \sin(\theta)u_{\theta}}{\partial\theta}\Big|_{r=R}\bigg)\bigg] \end{align}

Die \(\theta r\)-Komponente ist analog gegeben durch: \begin{align} \sigma_{\theta r}=\mu\bigg[\vec{\nabla}\vec{u}_{\theta r}+\vec{\nabla}\vec{u}_{r\theta} -\frac{2}{3}\vec{\nabla}\cdot\vec{u}\bigg] \end{align}

Einsetzen der jeweiligen Komponenten des Vektorgradienten und der Divergenz und Auswerten an \(r=R\) liefert: \begin{align} \sigma_{\theta r}=\mu\bigg[\frac{\partial u_{\theta}}{\partial r}\Big|_{r=R} +\frac{1}{R^2}\frac{\partial u_r}{\partial\theta}\Big|_{r=R} -\frac{u_{\theta}}{R}\Big|_{r=R}-\frac{2}{3}\bigg(\frac{1}{R^2} \frac{\partial r^2 u_r}{\partial r}\Big|_{r=R}+\frac{1}{R\sin(\theta)} \frac{\partial \sin(\theta)u_{\theta}}{\partial\theta}\Big|_{r=R}\bigg)\bigg] \end{align}

Einsetzen in das Integral zur Kraftberechnung liefert: \begin{align} &\phantom{\leftrightarrow,,} \vec{F}=2\pi R^2\int_{0}^{\pi}\bigg[-p(r=R,\theta) +\mu(r=R,\theta)\bigg[2\frac{\partial u_r}{\partial r}\Big|_{r=R} -\frac{2}{3}\bigg(\frac{1}{R^2} \frac{\partial r^2 u_r}{\partial r}\Big|_{r=R} +\frac{1}{R\sin(\theta)} \frac{\partial \sin(\theta)u_{\theta}} {\partial\theta}\Big|_{r=R}\bigg)\bigg]\vec{e}_r\\ &\phantom{\leftrightarrow,,,,,,} +\mu(r=R,\theta)\bigg[ \frac{\partial u_{\theta}}{\partial r}\Big|_{r=R}+\frac{1}{R^2}\frac{\partial u_r}{\partial\theta}\Big|_{r=R}-\frac{u_{\theta}}{R}\Big|_{r=R} -\frac{2}{3}\bigg(\frac{1}{R^2} \frac{\partial r^2 u_r}{\partial r}\Big|_{r=R} +\frac{1}{R\sin(\theta)} \frac{\partial \sin(\theta)u_{\theta}} {\partial\theta}\Big|_{r=R}\bigg)\bigg]\vec{e}_{\theta}\bigg] \sin(\theta)\text d\theta\\ &\Leftrightarrow \vec{F}=2\pi R^2\bigg[-\int_{0}^{\pi}p(r=R,\theta)\vec{e}_r\sin(\theta)d\theta +\int_{0}^{\pi}\mu(r=R,\theta)\bigg[ 2\frac{\partial u_r}{\partial r}\Big|_{r=R} -\frac{2}{3}\bigg(\frac{1}{R^2} \frac{\partial r^2 u_r}{\partial r}\Big|_{r=R} +\frac{1}{R\sin(\theta)} \frac{\partial \sin(\theta)u_{\theta}} {\partial\theta}\Big|_{r=R}\bigg)\bigg] \vec{e}_r\sin(\theta)\text d\theta\\ &\phantom{\leftrightarrow,,,,,,} +\int_{0}^{\pi}\mu(r=R,\theta)\bigg[ \frac{\partial u_{\theta}}{\partial r}\Big|_{r=R} +\frac{1}{R^2} \frac{\partial u_r}{\partial\theta}\Big|_{r=R} -\frac{u_{\theta}}{R}\Big|_{r=R} -\frac{2}{3}\bigg(\frac{1}{R^2} \frac{\partial r^2 u_r}{\partial r}\Big|_{r=R} +\frac{1}{R\sin(\theta)} \frac{\partial \sin(\theta)u_{\theta}} {\partial\theta}\Big|_{r=R}\bigg)\bigg]\bigg] \vec{e}_{\theta} \sin(\theta)\text d\theta\bigg] \end{align}

Mit der Produktregel folgt: \begin{align} &\phantom{\leftrightarrow,,} \vec{F}=2\pi R^2\bigg[-\int_{0}^{\pi}p(r=R,\theta)\vec{e}_r\sin(\theta)d\theta +\int_{0}^{\pi}\mu(r=R,\theta)\bigg[\frac{4}{3}\frac{\partial u_r}{\partial r}\Big|_{r=R} -\frac{4}{3R}u_r\Big|_{r=R} -\frac{2}{3R}\cot(\theta)u_{\theta}\Big|_{r=R} -\frac{2}{3R} \frac{\partial u_{\theta}}{\partial\theta}\Big|_{r=R}\bigg]\vec{e}_r\sin(\theta) \text d\theta\\ &\phantom{\leftrightarrow,,,,,,} +\int_{0}^{\pi}\mu(r=R,\theta)\bigg[ \frac{\partial u_{\theta}}{\partial r}\Big|_{r=R} -\frac{2}{3}\frac{\partial u_r}{\partial r}\Big|_{r=R} -\frac{4}{3R}u_r\Big|_{r=R}+\frac{1}{R^2} \frac{\partial u_r}{\partial\theta}\Big|_{r=R} -\frac{2}{3R}\cot(\theta)u_{\theta}\Big|_{r=R} -\frac{2}{3R} \frac{\partial u_{\theta}}{\partial\theta}\Big|_{r=R} \bigg]\vec{e}_{\theta} \sin(\theta)\text d\theta\bigg] \end{align}

Um Schreibarbeit zu sparen, werden die Größen \(\alpha\) und \(\beta\) eingeführt, welche gegeben sind durch: \begin{align} &\alpha(r=R,\theta):=\frac{4}{3}\frac{\partial u_r}{\partial r}\Big|_{r=R} -\frac{4}{3R}u_r\Big|_{r=R} -\frac{2}{3R}\cot(\theta)u_{\theta}\Big|_{r=R} -\frac{2}{3R} \frac{\partial u_{\theta}}{\partial\theta}\Big|_{r=R}\\ &\beta(r=R,\theta):=\frac{\partial u_{\theta}}{\partial r}\Big|_{r=R} -\frac{2}{3}\frac{\partial u_r}{\partial r}\Big|_{r=R} -\frac{4}{3R}u_r\Big|_{r=R}+\frac{1}{R^2} \frac{\partial u_r}{\partial\theta}\Big|_{r=R} -\frac{2}{3R}\cot(\theta)u_{\theta}\Big|_{r=R} -\frac{2}{3R} \frac{\partial u_{\theta}}{\partial\theta}\Big|_{r=R} \end{align}

Damit folgt für die Kraft: \begin{align} \vec{F}=2\pi R^2\bigg[ -\int_{0}^{\pi}p(r=R,\theta)\vec{e}_r \sin(\theta)d\theta +\int_{0}^{\pi}\mu(r=R,\theta) \big[\alpha(r=R,\theta)\vec{e}_r +\beta(r=R,\theta)\vec{e}_{\theta}\big] \sin(\theta)\text d\theta\bigg] \end{align}

Da in vorliegenden Problem jedoch nicht die Kraft in radialer bzw. polarer Richtung von Interesse ist, sondern insbesondere die Kraft in \(z\)-Richtung, muss nun die Projektion des Kraftvektors auf die \(z\)-Achse bestimmt werden. Es entsteht jedoch ein Problem, wenn die Projektion des ausintegrierten Kraftvektors berechnet wird. Dies liegt daran, da der Basisvektor \(\vec{e}_z\) in Kugelkoordinaten den Polarwinkel enthält und somit von der Wahl des Koordinatensystems abhängt. Aufgrund dessen, dass der ausintegrierte Kraftvektor ein invarianter Vektor ist, ist die physikalische Bedeutung dieser Kraft unabhängig vom Koordinatensystem.\\ Die, auf ein bestimmtes Flächenelement \(\text d\vec{A}\) der Kugel wirkende Kraft \(\text d\vec{F}\), ist jedoch abhängig von der Wahl des Koordinatensystems. Daher lässt sich auch die Projektion der auf dieses bestimmte Flächenelement wirkenden Kraft, auf die z-Achse bestimmen: \begin{align} \text d F_z=\text d\vec{F}\cdot\vec{e}_z \end{align}

Indem nun die Projektion, über gegebene Kugel mit Radius \(r=R\), durch eine analoge Rechnung, wie oben aufintegriert wird, folgt für die Kraft in \(z\)-Richtung: \begin{align} F_z=2\pi R^2\bigg[ -\int_{0}^{\pi}p(r=R,\theta)\vec{e}_r\cdot\vec{e}_z \sin(\theta)d\theta +\int_{0}^{\pi}\mu(r=R,\theta) \big[\alpha(r=R,\theta)\vec{e}_r\cdot\vec{e}_z +\beta(r=R,\theta)\vec{e}_{\theta}\cdot\vec{e}_z\big] \sin(\theta)\text d\theta\bigg] \end{align}

Der Basisvektor \(\vec{e}_z\) ist dabei in Kugelkoordinaten gegeben durch: \begin{align} \vec{e}_z=\cos(\theta)\vec{e}_r -\sin(\theta)\vec{e}_{\theta} \end{align}

Damit folgt für die z-Komponente der Kraft: \begin{align} F_z=2\pi R^2\bigg[ -\int_{0}^{\pi}p(r=R,\theta)\cos(\theta) \sin(\theta)d\theta +\int_{0}^{\pi}\mu(r=R,\theta) \big[\alpha(r=R,\theta)\cos(\theta) -\beta(r=R,\theta)\sin(\theta)\big] \sin(\theta)\text d\theta\bigg] \end{align}

Eine analoge Rechnung kann für jede andere Komponente durchgeführt werden, wobei \(x\) und \(y\) Null liefern.