Ziel dieses Kapitels ist, die Impulsdiffusion und Dissipation auf dem dreidimensionalen deformierten C-Gitter zu diskretisieren.
Zunächst geht es darum, die Impulsdiffusion im Rahmen der shallow water equations auf dem zweidimensionalen regulären hexagonalen Gitter zu diskretisieren. Im zweidimensionalen Fall kann man laut Glg. (8.50) für die Reibungsbeschleunigung schreiben
\[ \begin{align} \frac{\partial\mathbf{v}_\text{diff}}{\partial t} = \frac{1}{h}\nabla\cdot\left[2K\left(\begin{array}{cc} S_{1, 1} & S_{1, 2}\\ S_{2, 1} & S_{2, 2} \end{array}\right)\right] = \frac{1}{h}\nabla\cdot\left[K\left(\begin{array}{cc} 2\frac{\partial v_1}{\partial x_1} & \frac{\partial v_1}{\partial x_2} + \frac{\partial v_2}{\partial x_1}\\ \frac{\partial v_2}{\partial x_1} + \frac{\partial v_1}{\partial x_2} & 2\frac{\partial v_2}{\partial x_2} \end{array}\right)\right].\tag{29.1}\label{eq:momentum_diff_hex_deriv_0} \end{align} \]
Hierbei wurde die Dichte $\rho$ durch die Schichtdicke $h$ ersetzt, $K$ ist ein allgemeiner massengewichteter Diffusionskoeffizient. Die dreidimensionale Divergenz $\nabla\cdot\mathbf{v}$ verschwindet im inkompressiblen Fall. Als Ansatz für die Diskretisierung verwendet man
\[ \begin{align} \frac{\partial\mathbf{v}^T_\text{diff}}{\partial t} = \left(\delta_1, \delta_2, \delta_3\right)\cdot\tau = \left(\delta_1, \delta_2, \delta_3\right)\cdot\left( \begin{array}{ccc} \tau_{1, 1} & \tau_{1, 2} & \tau_{1, 3}\\ \tau_{2, 1} & \tau_{2, 2} & \tau_{2, 3} \\ \tau_{3, 1} & \tau_{3, 2} & \tau_{3, 3} \end{array}\right) \end{align} \]
Ziel dieses Abschnittes ist es, die $\tau_{i, j}$ eindeutig zu bestimmen. Motiviert durch Glg. (29.1) fordert man die Symmetrie von $\tau$, also
\[ \begin{align} \tau_{i, j} = \tau_{j, i}. \end{align} \]
Der Tensor $\tau$ wird in mehreren modifizierenden Interationen entwickelt. Man setzt zunächst an
\[ \begin{align} \tau = K\left( \begin{array}{ccc} 2\delta_1u_1 & \delta_2u_1 + \delta_1u_2 & \delta_3u_1 + \delta_1u_3\\ \delta_2u_1 + \delta_1u_2 & 2\delta_2u_2 & \delta_3u_2 + \delta_2u_3 \\ \delta_3u_1 + \delta_1u_3 & \delta_3u_2 + \delta_2u_3 & 2\delta_3u_3 \end{array}\right). \end{align} \]
Hieraus folgt
\[ \begin{align} \frac{\partial u_{1, \text{diff}}}{\partial t} &= K\left[\delta_1\left(2\delta_1u_1\right) + \delta_2\left(\delta_2u_1 + \delta_1u_2\right) + \delta_3\left(\delta_3u_1 + \delta_1u_3\right)\right]\nonumber\\ &= K\left(\delta_{1, 1}u_1 + \delta_{1, 1}u_1 + \delta_{2, 2}u_1 + \delta_{3, 3}u_1 + \delta_{2, 1}u_2 + \delta_{3, 1}u_3\right)\nonumber\\ &= K\left[\delta_{1, 1}u_1 + \delta_{2, 2}u_1 + \delta_{3, 3}u_1 + \delta_1\left(\delta_1u_1 + \delta_2u_2 + \delta_3u_3\right)\right] = K\left(\frac{3}{2}\Delta u_1 + \frac{3}{2}\delta_1D\right). \end{align} \]
Somit folgt aus
\[ \begin{align} \tau = K\left( \begin{array}{ccc} \frac{4}{3}\delta_1u_1 & \frac{2}{3}\left(\delta_2u_1 + \delta_1u_2\right) & \frac{2}{3}\left(\delta_3u_1 + \delta_1u_3\right)\\ \frac{2}{3}\left(\delta_2u_1 + \delta_1u_2\right) & \frac{4}{3}\delta_2u_2 & \frac{2}{3}\left(\delta_3u_2 + \delta_2u_3\right) \\ \frac{2}{3}\left(\delta_3u_1 + \delta_1u_3\right) & \frac{2}{3}\left(\delta_3u_2 + \delta_2u_3\right) & \frac{4}{3}\delta_3u_3 \end{array}\right) \end{align} \]
für die Reibungsbeschleunigung
\[ \begin{align} \frac{\partial u_{1, \text{diff}}}{\partial t} &= K\left(\Delta u_1 + \delta_1D\right). \end{align} \]
Der Term proportional zu $\delta_1D$ tritt im Rahmen der shallow water equations, da diese ein inkompressibles Medium beschreiben, nicht auf und ist daher möglichst zu eliminieren oder wenigstens zu reduzieren. Die Terme welche nicht auf der Diagonalen stehen dürfen nicht modifiziert werden, um die Symmetrie des Tensors nicht zu zerstören. Daher setzt man mit einem Parameter $0 \leq \alpha \leq 1$ an
\[ \begin{align} \tau &= K\left( \begin{array}{ccc} \frac{4}{3}\delta_1u_1 - \alpha D & \frac{2}{3}\left(\delta_2u_1 + \delta_1u_2\right) & \frac{2}{3}\left(\delta_3u_1 + \delta_1u_3\right)\\ \frac{2}{3}\left(\delta_2u_1 + \delta_1u_2\right) & \frac{4}{3}\delta_2u_2 - \alpha D & \frac{2}{3}\left(\delta_3u_2 + \delta_2u_3\right) \\ \frac{2}{3}\left(\delta_3u_1 + \delta_1u_3\right) & \frac{2}{3}\left(\delta_3u_2 + \delta_2u_3\right) & \frac{4}{3}\delta_3u_3 - \alpha D \end{array}\right).\nonumber\\ &= K\left( \begin{array}{ccc} \frac{2}{3}\left(2 - \alpha\right)\delta_1u_1 - \alpha\frac{2}{3}\left(\delta_2u_2 + \delta_3u_3\right) & \frac{2}{3}\left(\delta_2u_1 + \delta_1u_2\right) & \frac{2}{3}\left(\delta_3u_1 + \delta_1u_3\right)\\ \frac{2}{3}\left(\delta_2u_1 + \delta_1u_2\right) & \frac{2}{3}\left(2 - \alpha\right)\delta_2u_2 - \alpha\frac{2}{3}\left(\delta_1u_1 + \delta_3u_3\right) & \frac{2}{3}\left(\delta_3u_2 + \delta_2u_3\right) \\ \frac{2}{3}\left(\delta_3u_1 + \delta_1u_3\right) & \frac{2}{3}\left(\delta_3u_2 + \delta_2u_3\right) & \frac{2}{3}\left(2 - \alpha\right)\delta_3u_3 - \alpha\frac{2}{3}\left(\delta_1u_1 + \delta_2u_2\right) \end{array}\right). \end{align} \]
Hieraus folgt nun
\[ \begin{align} \frac{\partial u_{1, \text{diff}}}{\partial t} &= K\left[\Delta u_1 + \left(1 - \alpha\right)\delta_1D\right]. \end{align} \]
Der Faktor $1 - \alpha$ sollte möglichst nahe an Null sein, also sollte $\alpha$ möglichst nahe bei Eins sein. Nach Glg. (8.66) gilt für die Dissipationsrate
\[ \begin{align} \frac{3}{2}h\epsilon &= \frac{2}{3}\left(2 - \alpha\right)\left(\delta_1u_1\right)^2 - \alpha\frac{2}{3}\left(\delta_1u_1\delta_2u_2 + \delta_1u_1\delta_3u_3\right) + \frac{2}{3}\left(\delta_1u_2\delta_2u_1 + \delta_1u_2\delta_1u_2\right) + \frac{2}{3}\left(\delta_1u_3\delta_3u_1 + \delta_1u_3\delta_1u_3\right)\nonumber\\ &+ \frac{2}{3}\left(\delta_2u_1\delta_2u_1 + \delta_2u_1\delta_1u_2\right) + \frac{2}{3}\left(2 - \alpha\right)\left(\delta_2u_2\right)^2 - \alpha\frac{2}{3}\left(\delta_2u_2\delta_1u_1 + \delta_2u_2\delta_3u_3\right) + \frac{2}{3}\left(\delta_2u_3\delta_3u_2 + \delta_2u_3\delta_2u_3\right)\nonumber\\ &+ \frac{2}{3}\left(\delta_3u_1\delta_3u_1 + \delta_3u_1\delta_1u_3\right) + \frac{2}{3}\left(\delta_3u_2\delta_3u_2 + \delta_3u_2\delta_2u_3\right) + \frac{2}{3}\left(2 - \alpha\right)\left(\delta_3u_3\right)^2 - \alpha\frac{2}{3}\left(\delta_3u_3\delta_1u_1 + \delta_3u_3\delta_2u_2\right)\nonumber\\ &= \frac{2}{3}\left(2 - \alpha\right)\left(\delta_1u_1\right)^2 - \alpha\frac{2}{3}\left(\delta_1u_1\delta_2u_2 + \delta_1u_1\delta_3u_3\right) + \textcolor{red}{\frac{2}{3}\left(\delta_2u_1 + \delta_1u_2\right)^2} + \textcolor{red}{\frac{2}{3}\left(\delta_3u_1 + \delta_1u_3\right)^2}\nonumber\\ &+ \frac{2}{3}\left(2 - \alpha\right)\left(\delta_2u_2\right)^2 - \alpha\frac{2}{3}\left(\delta_2u_2\delta_1u_1 + \delta_2u_2\delta_3u_3\right) + \textcolor{red}{\frac{2}{3}\left(\delta_3u_2 + \delta_2u_3\right)^2}\nonumber\\ &+ \frac{2}{3}\left(2 - \alpha\right)\left(\delta_3u_3\right)^2 - \alpha\frac{2}{3}\left(\delta_3u_3\delta_1u_1 + \delta_3u_3\delta_2u_2\right). \end{align} \]
Die rot markierten Terme entstehen aus den Nicht-Diagonalelementen und sind alle quadratisch und somit nichtnegativ. Für die verbleibenden Terme gilt
\[ \begin{align} &\frac{2}{3}\left(2 - \alpha\right)\left(\delta_1u_1\right)^2 - \alpha\frac{2}{3}\left(\delta_1u_1\delta_2u_2 + \delta_1u_1\delta_3u_3\right) + \frac{2}{3}\left(2 - \alpha\right)\left(\delta_2u_2\right)^2 - \alpha\frac{2}{3}\left(\delta_2u_2\delta_1u_1 + \delta_2u_2\delta_3u_3\right)\nonumber\\ &+ \frac{2}{3}\left(2 - \alpha\right)\left(\delta_3u_3\right)^2 - \alpha\frac{2}{3}\left(\delta_3u_3\delta_1u_1 + \delta_3u_3\delta_2u_2\right)\nonumber\\ &= \frac{2}{3}\left(2 - \alpha\right)\left(\delta_1u_1\right)^2 + \frac{2}{3}\left(2 - \alpha\right)\left(\delta_2u_2\right)^2 + \frac{2}{3}\left(2 - \alpha\right)\left(\delta_3u_3\right)^2 - \alpha\frac{4}{3}\left(\delta_1u_1\delta_2u_2 + \delta_1u_1\delta_3u_3 + \delta_2u_2\delta_3u_3\right)\nonumber\\ &\hastobe \alpha\frac{2}{3}\left(\delta_1u_1 - \delta_2u_2\right)^2 + \alpha\frac{2}{3}\left(\delta_1u_1 - \delta_3u_3\right)^2 + \alpha\frac{2}{3}\left(\delta_2u_2 - \delta_3u_3\right)^2. \end{align} \]
Die Forderung im letzten Schritt begründet sich darauf, eine nichtnegative Dissipationsrate erhalten zu wollen. Es ergibt sich
\[ \begin{align} \frac{2}{3}\left(2 - \alpha\right) = \frac{4}{3}\alpha \Leftrightarrow \frac{4}{3} - \alpha\frac{2}{3} = \frac{4}{3}\alpha \Rightarrow \frac{4}{3} = \frac{6}{3}\alpha \Leftrightarrow\alpha = \frac{2}{3}. \end{align} \]
Mit dieser Wahl von $\alpha$ tragen auch die Terme auf der Hauptdiagonalen nichtnegativ zur Dissipationsrate bei und es gilt
\[ \begin{align} \epsilon \geq 0. \end{align} \]
Das Ergebnis für den Reibungstensor lautet somit
\[ \begin{align} \tau &= K\left( \begin{array}{ccc} \frac{4}{9}\left[2\delta_1u_1 - \left(\delta_2u_2 + \delta_3u_3\right)\right] & \frac{2}{3}\left(\delta_2u_1 + \delta_1u_2\right) & \frac{2}{3}\left(\delta_3u_1 + \delta_1u_3\right)\\ \frac{2}{3}\left(\delta_2u_1 + \delta_1u_2\right) & \frac{4}{9}\left[2\delta_2u_2 - \left(\delta_1u_1 + \delta_3u_3\right)\right] & \frac{2}{3}\left(\delta_3u_2 + \delta_2u_3\right)\\ \frac{2}{3}\left(\delta_3u_1 + \delta_1u_3\right) & \frac{2}{3}\left(\delta_3u_2 + \delta_2u_3\right) & \frac{4}{9}\left[2\delta_3u_3 - \left(\delta_1u_1 + \delta_2u_2\right)\right] \end{array}\right). \end{align} \]
Die Verallgemeinerung
\[ \begin{align} \tau &= \left( \begin{array}{ccc} K_c\frac{4}{9}\left[2\delta_1u_1 - \left(\delta_2u_2 + \delta_3u_3\right)\right] & K_1\frac{2}{3}\left(\delta_2u_1 + \delta_1u_2\right) & K_2\frac{2}{3}\left(\delta_3u_1 + \delta_1u_3\right)\\ K_1\frac{2}{3}\left(\delta_2u_1 + \delta_1u_2\right) & K_c\frac{4}{9}\left[2\delta_2u_2 - \left(\delta_1u_1 + \delta_3u_3\right)\right] & K_3\frac{2}{3}\left(\delta_3u_2 + \delta_2u_3\right)\\ K_2\frac{2}{3}\left(\delta_3u_1 + \delta_1u_3\right) & K_3\frac{2}{3}\left(\delta_3u_2 + \delta_2u_3\right) & K_c\frac{4}{9}\left[2\delta_3u_3 - \left(\delta_1u_1 + \delta_2u_2\right)\right] \end{array}\right)\tag{29.15}\label{eq:stress_tensor_hex_regular} \end{align} \]
zerstört keine der überprüften Aussagen.
Die Divergenz des Tensors Glg. (29.15) erfüllt noch nicht die Thuburn-Bedingung. Motiviert durch Abschn. 27.6.1 macht man für die Reibungsbeschleunigung an der Kante $e$ daher den Ansatz
\[ \begin{align} h_e\frac{\partial u_{e, \text{diff}}}{\partial t} &= \left[\nabla\left(K_cE_{c(e)}\right)\right]_e + \frac{1}{l_e}\left(\frac{1}{3}\sum_{r\in t_1}K_rF_{r(e)} - \frac{1}{3}\sum_{r\in t_2}K_rF_{r(e)}\right)\nonumber\\ &= \left[\nabla\left(K_cE_{c(e)}\right)\right]_e + \frac{1}{3l_e}\left(\sum_{r\in t_1}K_rF_{r(e)} - \sum_{r\in t_2}K_rF_{r(e)}\right). \end{align} \]
Hierbei sind $E_c$ und $F_r$ Deformationen, welche in den Zellzentren bzw. an den Kanten lokalisiert sind. Die $E_c$ sind auf der Diagonalen des Stress-Tensors lokalisiert, sie werden als Schubdeformationen (engl. strain deformation) bezeichnet. Die $F_r$ sind nicht auf der Diagonalen des Stress-Tensors lokalisiert, sie werden als Scherdeformationen (engl. shear deformation) bezeichnet. Man lässt hierbei zu, dass beide außer von ihrer Position auch von der Kante abhängen können, an der man die Beschleunigung berechnen möchte, also $E_{c(e)}$ bzw. $F_{r(e)}$. Das Dreieck $t_1$ ist in Bezug auf die Kante $e$ in der Richtung $\mathbf{k}\times\mathbf{n}_e$ lokalisiert, die Kante $t_2$ in Richtung $-\mathbf{k}\times\mathbf{n}_e$. $K_c$ und $K_r$ sind nichtnegative Diffusionskoeffizienten, welche in diesem Kapitel offengelassen werden. O. B. d. A. nimmt man an, dass $\mathbf{n}_e$ parallel zu $\mathbf{i}_1$ ausgerichtet ist.
Auf dem regulären Gitter erhält man für die Schubdeformation
\[ \begin{align} E_{c(e)} &= \frac{4}{9}\left[2\delta_1u_1 - \left(\delta_2u_2 + \delta_3u_3\right)\right] = \frac{4}{9}\left[3\delta_1u_1 - \left(\delta_1u_1 + \delta_2u_2 + \delta_3u_3\right)\right] = \frac{4}{3}\delta_1u_1 - \frac{4}{9}\left(\delta_1u_1 + \delta_2u_2 + \delta_3u_3\right)\nonumber\\ &= -\frac{4}{9}\left(\delta_1u_1 + \delta_2u_2 + \delta_3u_3\right) + \frac{4}{3}\delta_1u_1 = -\frac{2}{3}\frac{2}{3}\left(\delta_1u_1 + \delta_2u_2 + \delta_3u_3\right) + \frac{4}{3}\delta_1u_1\nonumber\\ &\stackrel{\href{ch-26-probleme-eines-dreielementigen-erzeugend.html#eq:div_3elements2d}{\text{Glg. (27.33)}}}{=} -\frac{2}{3}D + \frac{4}{3}\delta_1u_1 \end{align} \]
Für die Scherdeformationen muss laut Glg. (29.15)
\[ \begin{align} &\delta_2\left[\frac{2}{3}\left(\delta_2u_1 + \delta_1u_2\right)\right] + \delta_3\left[\frac{2}{3}\left(\delta_3u_1 + \delta_1u_3\right)\right] \hastobe \frac{1}{3}\delta_y\left(F_{1(e)} + F_{2(e)} + F_{3(e)}\right)\nonumber\\ &\stackrel{\href{ch-26-probleme-eines-dreielementigen-erzeugend.html#eq:ddy_hex}{\text{Glg. (27.31)}}}{=} \frac{1}{3\sqrt{3}}\left(\delta_2 - \delta_3\right)\left(F_{1(e)} + F_{2(e)} + F_{3(e)}\right)\tag{29.18}\label{eq:deriv_hex_diff_deformed_0} \end{align} \]
gelten. Die Wahl von $F_{1(e)}$ ist hierbei willkürlich, da dieser Term zunächst addiert und dann wieder abgezogen wird. Man wählt daher naheliegenderweise
\[ \begin{align} F_{1(e)} &= \zeta_1. \end{align} \]
Zunächst behandelt man Glg. (29.18) in der Form
\[ \begin{align} &\delta_2\left[\frac{2}{3}\left(\delta_2u_1 + \delta_1u_2\right)\right] + \delta_3\left[\frac{2}{3}\left(\delta_3u_1 + \delta_1u_3\right)\right] = \frac{1}{\sqrt{3}}\left(\delta_2F_{3(e)} - \delta_3F_{2(e)}\right). \end{align} \]
Hiermit erhält man
\[ \begin{align} F_{2(e)} &= -\frac{2}{\sqrt{3}}\left(\delta_3u_1 + \delta_1u_3\right) = \frac{2}{\sqrt{3}}\left(-\delta_3u_1 - \delta_1u_3\right) = \frac{2}{\sqrt{3}}\left[\left(\delta_3u_1 - \delta_1u_3\right) - 2\delta_3u_1\right] \stackrel{\href{ch-26-probleme-eines-dreielementigen-erzeugend.html#eq:hex_curl_rhombus_1}{\text{Glg. (27.171)}}}{=} \zeta_2 - \frac{4}{\sqrt{3}}\delta_3u_1,\\ F_{3(e)} &= \frac{2}{\sqrt{3}}\left(\delta_2u_1 + \delta_1u_2\right) = \frac{2}{\sqrt{3}}\left[\left(\delta_1u_2 - \delta_2u_1\right) + 2\delta_2u_1\right] \stackrel{\href{ch-26-probleme-eines-dreielementigen-erzeugend.html#eq:hex_curl_rhombus_2}{\text{Glg. (27.172)}}}{=} \zeta_3 + \frac{4}{\sqrt{3}}\delta_2u_1. \end{align} \]
Analog zu Abschn. 27.6.1 verwendet man nun die gemittelten Scherdeformationen, was auf Glg. (29.18) führt. Das Verfahren zur Berechnung des Gradienten der Komponente $u_1$ wird in Abschn. 29.2.0.1 im Rahmen des deformierten Gitters vorgestellt.
In diesem Abschnitt werden die Resultate des vorherigen Abschnittes auf das deformierte hexagonale C-Gitter verallgemeinert. Auch hier bleibt man im Rahmen der shallow water equations. Im Wesentlichen kann man die gleichen Formeln wie auf dem regulären Gitter verwenden, jedoch können die Geschwindigkeitsgradienten etwas anders notiert werden. Wegen
\[ \begin{align} F_{2(e)} + F_{3(e)} &= \zeta_2 - \frac{4}{\sqrt{3}}\delta_3u_1 + \zeta_3 + \frac{4}{\sqrt{3}}\delta_2u_1 = \zeta_2 + \zeta_3 + 4\frac{1}{\sqrt{3}}\left(\delta_2 - \delta_3\right)u_1 = \zeta_2 + \zeta_3 + 4\delta_yu_1\nonumber\\ &= \zeta_2 + 2\delta_yu_1 + \zeta_3 + 2\delta_yu_1 \end{align} \]
kann man die Deformationen auch in der Form
\[ \begin{align} E_{c(e)} &= -\frac{2}{3}D + \frac{4}{3}\delta_xu_1\vline_c,\\ F_{1(e)} &= \zeta_1,\\ F_{2(e)} &= \zeta_2 + 2\delta_yu_1\vline_2,\\ F_{3(e)} &= \zeta_3 + 2\delta_yu_1\vline_3 \end{align} \]
aufschreiben. Hierbei werden die Geschwindigkeitsgradienten über Zellen bzw. Rhomben berechnet.
Für die Dissipationsrate $\epsilon_\text{hex}$ kann einfach
\[ \begin{align} \epsilon_\text{hex} = -\rho\frac{\partial\mathbf{v}_\text{diff}}{\partial t}\cdot\mathbf{v} \end{align} \]
angesetzt werden. Dies entspricht zwar nicht exakt der molekularen Dissipationsrate Glg. (8.66), jedoch sind die globalen Integrale über beide Ansätze nach Glg. (8.96) identisch. Somit ist energetische Selbstkonsistenz gegeben.