21 Ozeanzirkulation

Obwohl dies ein Buch über Meteorologie ist, gibt es aus drei Gründen ein Kapitel über den Ozean:

• Zwei Drittel der Erdoberfläche sind vom Ozean bedeckt. Der Ozean bildet somit einen großen Teil der unteren Randbedingung der Atmosphäre und ist naturgemäß entscheidend für den Wasserhaushalt.

• Ozean und Atmosphäre gehören zusammen zum Klimasystem und wechselwirken miteinander. Die Ozeanzirkulation beeinflusst das Klima der Atmosphäre, indem sie Wärme transportiert und moduliert somit die Position der Fronten, die Zyklogenese und meridionale klimatologische Unterschiede.

• Die theoretische Vorarbeit, die bereits für die Atmosphäre gemacht wurde, lässt sich recht leicht auf den Ozean anwenden, weshalb es schade wäre, dies nicht zu tun.

21.1 Zustandsgrößen des Ozeans

Die relevanten Zustandsvariablen im Ozean sind die drei Geschwindigkeitskomponenten $u$, $v$, $w$, die Oberflächenauslenkung $\eta$, die thermodynamischen Variablen $p$ und $T$ sowie der Salzgehalt (auch: Salinität) $S$. Alternativ kann man auch hieraus bijektiv ableitbare Größen verwenden, wie z. B. generalisierte Geschwindigkeiten, die potentielle Dichte oder die potentielle Temperatur. Dies ist analog zu den Zustandsvariablen der feuchten Atmosphäre ohne Kondensate. Die Rolle der Feuchte wird dabei von der Salinität übernommen. Ein Unterschied zur Atmosphäre besteht darin, dass die diabatischen Terme insbesondere in tieferen Schichten aufgrund der schwachen Strahlung, geringen Strömungsgeschwindigkeit und somit geringen Reibung sowie der Abwesenheit von Phasenübergängen sehr schwach sind. Dies führt dazu, dass sich Teilchen im Wesentlichen entlang von Isopyknen bewegen, dies sind Flächen gleicher potentieller Dichte.

21.2 Thermische Zustandsgleichung

Für Flüssigkeiten existieren im Gegensatz zu idealen Gasen keine einfachen, analytisch herleitbaren und hinreichend genauen thermischen Zustandsgleichungen. Bei Seewasser kommt die zusätzliche Komplikation des Salzes hinzu. Daher existiert eine Reihe empirischer Zustandsgleichungen für Seewasser. Problematisch bei empirischen, also an Messdaten approximierten Ausdrücken ist jedoch ihre fehlende thermodynamische Konsistenz: angenommen, man hat einen empirischen Ausdruck für die Temperatur als Funktion von Volumen, Entropie und Teilchenzahl und einen anderen empirischen Ausdruck für den Druck als Funktion von Volumen, Entropie und Teilchenzahl, ist nicht gewährleistet, dass die Maxwell-Relation Glg. (5.97) gilt. Um dieses Problem zu lösen, ist man für Seewasser dazu übergegangen, eine Formel für ein thermodynamisches Potential zu verwenden, und andere thermodynamische Größen durch partielles Differenzieren dieser Formel zu berechnen. Dies führt automatisch zu thermodynamischer Konsistenz. Die moderne standardisierte Gleichung, die sich dieses Prinzips bedient, ist TEOS-10 (Thermodynamic Equation of Seawater 2010) [32].

Diese Gleichung verwendet das Gibbs-Potential $G = G\left(T,p,m\right)$ bzw. das spezifische Gibbs-Potential $g\left(T,p\right) = G\left(T,p,m\right)/m$. Diese Transformation ist möglich, da das Gibbs-Potential wegen Glg. (5.85) eine extensive Zustandsgröße ist. Hieraus lässt sich mittels Glg. (5.95) die Massendichte berechnen:

\[ \begin{align} \rho = \frac{m}{V} = \frac{1}{\frac{V}{m}} = \frac{1}{\frac{\partial g}{\partial p}} \end{align} \]

Hierbei wurde der Index $T,N$ bzw. $T,m$ weggelassen, da diese Größen wegen $g = g\left(T,p\right)$ ohnehin bei dieser partiellen Ableitung konstant sind.

21.3 Windgetriebene Zirkulation

Die in Abschn. 17.3 durchgeführte Herleitung der Ekman-Spirale kann man ohne Modifikationen auf einen ebenen Ozeanboden übertragen. An der Ozeanoberfläche hingegen fehlt die Adhäsionsbedinung $\lim_{\mathbf{r} \to \partial V}\mathbf{v} = \mathbf{0}$, vielmehr wäre die Annahme sinnvoll, dass an der Ozeanoberfläche Windgeschwindigkeit und Strömungsgeschwindigkeit gegen einen gemeinsamen Grenzwert konvergieren. Man müsste dann die Ekman-Spirale simultan in Ozean und Atmosphäre lösen.

21.3.1 Sverdrup-Balance

Es gibt jedoch eine einfachere Möglichkeit, den Einfluss des Windes auf die Zirkulation im Ozean grundlegend zu untersuchen. Hierzu setzt man zunächst das Kräftegleichgewicht der Glg.en (17.35) - (17.36) an:

\[ \begin{align} f\mathbf{k}\times\mathbf{v}_{h} = -\nabla\phi + \frac{1}{\rho_0}\frac{\partial\mathbf{\tau}}{\partial z} \end{align} \]

Man geht nun davon aus, dass die vertikalen Dichtegradienten so klein sind, dass man den Term $\frac{1}{\rho_0}$ in die partielle Ableitung hereinziehen kann, also

\[ \begin{align} f\mathbf{k}\times\mathbf{v}_{h} = -\nabla\phi + \frac{\partial}{\partial z}\left(\frac{\mathbf{\tau}}{\rho_0}\right). \end{align} \]

Hierauf wendet man nun die Rotation an, also

\[ \begin{align} \nabla\times f\mathbf{k}\times\mathbf{v}_{h} = \nabla\times\frac{\partial}{\partial z}\left(\frac{\mathbf{\tau}}{\rho_0}\right).\tag{21.4}\label{eq:wind-driven_circ_deriv_1} \end{align} \]

Mit Glg. (B.53) folgt

\[ \begin{align} \nabla\times f\mathbf{k}\times\mathbf{v}_{h} &= \left(\mathbf{v}_{h}\cdot\nabla\right)f\mathbf{k} - \mathbf{v}_{h}\left(\nabla\cdot f\mathbf{k}\right) + f\mathbf{k}\nabla\cdot\mathbf{v}_{h} - \left(f\mathbf{k}\cdot\nabla\right)\mathbf{v}_{h}. \end{align} \]

Projektion auf die vertikale Richtung ergibt

\[ \begin{align} \mathbf{k}\cdot\nabla\times f\mathbf{k}\times\mathbf{v}_{h} &= \mathbf{k}\cdot\left(\mathbf{v}_{h}\cdot\nabla\right)f\mathbf{k} + f\nabla\cdot\mathbf{v}_{h}\nonumber\\ \Rightarrow \mathbf{k}\cdot\nabla\times f\mathbf{k}\times\mathbf{v}_{h} &= v\beta - f\frac{\partial w}{\partial z}. \end{align} \]

Integration über das Höhenintervall $\left[z_B, z_T\right]$ führt somit auf

\[ \begin{align} \beta\int_{z_B}^{z_T}vdz - f\left(w_T - w_B\right) = \mathbf{k}\cdot\frac{1}{\rho_0}\nabla\times\left(\mathbf{\tau}_T - \mathbf{\tau}_B\right). \end{align} \]

Mit den Definitionen

\[ \begin{align} D \coloneqq z_T - z_B, & {} & \newoverline{v} \coloneqq \frac{1}{D}\int_{z_B}^{z_T}vdz \end{align} \]

kann man dies kürzer als

\[ \begin{align} \beta D\newoverline{v} - f\left(w_T - w_B\right) = \mathbf{k}\cdot\frac{1}{\rho_0}\nabla\times\left(\mathbf{\tau}_T - \mathbf{\tau}_B\right) \end{align} \]

notieren. Da es hier um die ganze Wassersäule geht, kann man die einfache Randbedingung

\[ \begin{align} w_T = w_B = 0 \end{align} \]

anwenden. Damit erhält man

\[ \begin{align} \beta D\newoverline{v} = \mathbf{k}\cdot\frac{1}{\rho_0}\nabla\times\left(\mathbf{\tau}_T - \mathbf{\tau}_B\right). \end{align} \]

Vernachlässigt man hier auch noch die Bodenreibung, die aufgrund der geringen Strömungsgeschwindigkeit im tiefen Ozean deutlich kleiner ist als der Wind stress, erhält man die sogenannte Sverdrup-Balance, oder auch Sverdrup-Relation:

Sprachlich lautet die Sverdrup-Balance: „Die vertikal integrierte meridionale Geschwindigkeit ist gleich der Rotation der Windschubspannung dividiert durch den Rossby-Parameter.“

21.3.1.1 Stommel-Modell

Um den Reibungsterm am Grund zu parametrisieren, verwendet man das Stommel-Modell, welches lautet

\[ \begin{align} \mathbf{\tau}_B = r\newoverline{\mathbf{v}_{h}} \end{align} \]

mit einer Konstante $r > 0$. Mit $w = 0$ gilt auch $\nabla\cdot\mathbf{v}_{h} = 0$, weshalb man den gemittelten Wind durch eine Stromfunktion $\psi$ ausdrücken kann:

\[ \begin{align} \newoverline{u} = -\frac{\partial\psi}{\partial y}, & {} & \newoverline{v} = \frac{\partial\psi}{\partial x}\tag{21.14}\label{eq:stommel_stream} \end{align} \]

Somit kann man notieren

\[ \begin{align} \beta D\frac{\partial\psi}{\partial x} + \mathbf{k}\cdot\frac{1}{\rho_0}\nabla\times\mathbf{\tau}_B &= \mathbf{k}\cdot\frac{1}{\rho_0}\nabla\times\mathbf{\tau}_T\nonumber\\ \Leftrightarrow \beta D\frac{\partial\psi}{\partial x} + \frac{r}{\rho_0}\Delta\psi &= \mathbf{k}\cdot\frac{1}{\rho_0}\nabla\times\mathbf{\tau}_T\nonumber \end{align} \]

mit

\[ \begin{align} \epsilon \coloneqq \frac{\rho_0\beta D}{r}. \end{align} \]

Glg. (21.15) ist die Bewegungsgleichung des Stommel-Modells. Hieraus lässt sich die Ozeanströmung bei gegebener Windeinwirkung diagnostisch ableiten. Hierfür wird Glg. (21.15) auf der Menge $\left[0, L_x\right] \times \left[0, L_y\right]$ mit $L_x,L_y > 0$ gelöst, also für ein rechteckiges Becken. Hierzu geht man von einem wind stress der Form

\[ \begin{align} \tau_T^{(x)} &= -C\rho_aU^2\cos\left(\pi\frac{y}{L_y}\right),\\ \tau_T^{(y)} &= 0 \end{align} \]

aus. Hierbei ist $\rho_a$ die Dichte der Atmosphäre. Hieraus folgt

\[ \begin{align} \frac{1}{r}\mathbf{k}\cdot\nabla\times\mathbf{\tau}_T = -C\frac{\pi}{rL_y}\rho_aU^2\sin\left(\pi\frac{y}{L_y}\right). \end{align} \]

Setzt man dies in Glg. (21.15) ein, erhält man

\[ \begin{align} \Delta\psi + \epsilon\frac{\partial\psi}{\partial x} &= -C\frac{\pi}{rL_y}\rho_aU^2\sin\left(\pi\frac{y}{L_y}\right).\tag{21.20}\label{eq:stommel_dynamics_mod} \end{align} \]

Die Randbedingungen sind

\[ \begin{align} u\left(0,y\right) = u\left(L_x,y\right) &= 0,\nonumber\\ v\left(x,0\right) = v\left(x,L_y\right) &= 0, \end{align} \]

da kein Wasser das Becken verlassen soll. Hieraus folgt

\[ \begin{align} \psi\left(0,y\right) = \psi\left(L_x,y\right) = \psi\left(x,0\right) = \psi\left(x,L_y\right) = 0. \end{align} \]

Nun geht man weiterhin davon aus, dass die Tiefe homogen ist, sowie dass $f$ linear von der Breite abhängt, womit $\beta$ homogen ist ($\beta-$Ebene). Glg. (21.20) ist linear, d. h. Summen von Lösungen sind wieder Lösungen. Weiterhin ist diese Gleichung inhomogen. Zunächst löst man den homogenen Anteil

\[ \begin{align} \Delta\psi + \epsilon\frac{\partial\psi}{\partial x} &= 0. \end{align} \]

Als Ansatz macht man

\[ \begin{align} \psi\left(x,y\right) = \left[A\exp\left(\lambda_1x\right) + B\exp\left(\lambda_2x\right)\right]\sin\left(\frac{m\pi y}{L_y}\right). \end{align} \]

Einsetzen ergibt

\[ \begin{align} &\left[A\lambda_1^2\exp\left(\lambda_1x\right) + B\lambda_2^2\exp\left(\lambda_2x\right)\right]\sin\left(\frac{m\pi y}{L_y}\right) - \frac{m^2\pi^2}{L_y^2}\left[A\exp\left(\lambda_1x\right) + B\exp\left(\lambda_2x\right)\right]\sin\left(\frac{m\pi y}{L_y}\right)\nonumber\\ &+ \epsilon\left[A\lambda_1\exp\left(\lambda_1x\right) + B\lambda_2\exp\left(\lambda_2x\right)\right]\sin\left(\frac{m\pi y}{L_y}\right) = 0\nonumber\\ &\Leftrightarrow\left[A\lambda_1^2\exp\left(\lambda_1x\right) + B\lambda_2^2\exp\left(\lambda_2x\right)\right] - \frac{m^2\pi^2}{L_y^2}\left[A\exp\left(\lambda_1x\right) + B\exp\left(\lambda_2x\right)\right] + \epsilon\left[A\lambda_1\exp\left(\lambda_1x\right) + B\lambda_2\exp\left(\lambda_2x\right)\right] = 0 \end{align} \]

Da $\exp\left(\lambda_1x\right)$ und $\exp\left(\lambda_2x\right)$ linear unabhängig voneinander sind, müssen die Summen ihrer Vorfaktoren jeweils gleich Null sein, also

\[ \begin{align} \lambda_1^2 + \epsilon\lambda_1 - \frac{m^2\pi^2}{L_y^2} &= 0,\\ \lambda_2^2 + \epsilon\lambda_2 - \frac{m^2\pi^2}{L_y^2} &= 0. \end{align} \]

Hieraus folgen mit der pq-Formel

\[ \begin{align} \lambda_1 &= -\frac{\epsilon}{2} + \sqrt{\frac{\epsilon^2}{4} + \frac{m^2\pi^2}{L_y^2}},\\ \lambda_2 &= -\frac{\epsilon}{2} - \sqrt{\frac{\epsilon^2}{4} + \frac{m^2\pi^2}{L_y^2}}. \end{align} \]

Nun benötigt man noch eine spezielle Lösung von Glg. (21.20) inklusive des inhomogenen Anteils. Hierfür wählt man

\[ \begin{align} \psi\left(x,y\right) = C\frac{L_y}{r\pi}\rho_aU^2\sin\left(\pi\frac{y}{L_y}\right). \end{align} \]

Die allgemeine Lösung lautet also

\[ \begin{align} \psi\left(x,y\right) = \left[A\exp\left(\lambda_1x\right) + B\exp\left(\lambda_2x\right)\right]\sin\left(\frac{m\pi y}{L_y}\right) + \frac{CL_y}{r\pi}\rho_aU^2\sin\left(\pi\frac{y}{L_y}\right). \end{align} \]

$A$, $B$ und $m$ können aus den Randbedingungen abgeleitet werden. Die Randbedingungen $\psi\left(x,0\right) = \psi\left(x,L_y\right) = 0$ sind erfüllt. Aus $\psi\left(0,y\right) = \psi\left(L_x,y\right) = 0$ für alle $0\leq y\leq L_y$ folgt zunächst $m=1$, da sich sonst die beiden sin-Terme nicht gegeneinander aufheben können. Es gilt also

\[ \begin{align} \psi\left(x,y\right) = \left[A\exp\left(\lambda_1x\right) + B\exp\left(\lambda_2x\right) + \frac{L_y}{r\pi}C\rho_aU^2\right]\sin\left(\pi\frac{y}{L_y}\right).\tag{21.32}\label{eq:stommel_deriv_2} \end{align} \]

Mit $\psi\left(0,y\right) = \psi\left(L_x,y\right) = 0$ erhält man nun das folgende lineare Gleichungssystem für $A$ und $B$:

\[ \begin{align} A + B &= -\frac{L_y}{r\pi}C\rho_aU^2,\tag{21.33}\label{eq:stommel_deriv_1}\\ A\exp\left(\lambda_1L_x\right) + B\exp\left(\lambda_2L_x\right) &= -\frac{L_y}{r\pi}C\rho_aU^2 \end{align} \]

Multipliziert man die erste Gleichung mit $\exp\left(\lambda_1L_x\right)$ und zieht das Resultat von der zweiten ab, erhält man

\[ \begin{align} B\left[\exp\left(\lambda_2L_x\right) - \exp\left(\lambda_1L_x\right)\right] &= \frac{L_y}{r\pi}C\rho_aU^2\left(\exp\left(\lambda_1L_x\right) - 1\right)\nonumber\\ \Leftrightarrow B&= \frac{L_y}{r\pi}C\rho_aU^2\frac{\exp\left(\lambda_1L_x\right) - 1}{\exp\left(\lambda_2L_x\right) - \exp\left(\lambda_1L_x\right)}. \end{align} \]

Setzt man dies in Glg. (21.33) ein, erhält man

\[ \begin{align} A &= -\frac{L_y}{r\pi}C\rho_aU^2 - B = -\frac{L_y}{r\pi}C\rho_aU^2 - \frac{L_y}{r\pi}C\rho_aU^2\frac{\exp\left(\lambda_1L_x\right) - 1}{\exp\left(\lambda_2L_x\right) - \exp\left(\lambda_1L_x\right)}\nonumber\\ &= -\frac{L_y}{r\pi}C\rho_aU^2\frac{\exp\left(\lambda_2L_x\right) - \exp\left(\lambda_1L_x\right)}{\exp\left(\lambda_2L_x\right) - \exp\left(\lambda_1L_x\right)} - \frac{L_y}{r\pi}C\rho_aU^2\frac{\exp\left(\lambda_1L_x\right) - 1}{\exp\left(\lambda_2L_x\right) - \exp\left(\lambda_1L_x\right)}\nonumber\\ &= -\frac{L_y}{r\pi}C\rho_aU^2\frac{\exp\left(\lambda_2L_x\right) - 1}{\exp\left(\lambda_2L_x\right) - \exp\left(\lambda_1L_x\right)}. \end{align} \]

Setzt man die Resultate für $A$ und $B$ in Glg. (21.32) ein und definiert weiter

\[ \begin{align} \gamma \coloneqq \exp\left(\lambda_2L_x\right) - \exp\left(\lambda_1L_x\right), \end{align} \]

folgt

\[ \begin{align} \psi\left(x,y\right) &= \frac{L_y}{r\pi}C\rho_aU^2\left[\frac{1 - \exp\left(\lambda_2L_x\right)}{\gamma}\exp\left(\lambda_1x\right) + \frac{\exp\left(\lambda_1L_x\right) - 1}{\gamma}\exp\left(\lambda_2x\right) + 1\right]\sin\left(\pi\frac{y}{L_y}\right).\tag{21.38}\label{eq:stommel_solution} \end{align} \]

Setzt man dies in Glg. (21.14) ein, erhält man

\[ \begin{align} u\left(x,y\right) &= -\frac{C\rho_aU^2}{r}\left[\frac{1 - \exp\left(\lambda_2L_x\right)}{\gamma}\exp\left(\lambda_1x\right) + \frac{\exp\left(\lambda_1L_x\right) - 1}{\gamma}\exp\left(\lambda_2x\right) + 1\right]\cos\left(\pi\frac{y}{L_y}\right),\tag{21.39}\label{eq:stommel_result_u}\\ v\left(x,y\right) &= \frac{L_y}{r\pi}C\rho_aU^2\left[\lambda_1\frac{1 - \exp\left(\lambda_2L_x\right)}{\gamma}\exp\left(\lambda_1x\right) + \lambda_2\frac{\exp\left(\lambda_1L_x\right) - 1}{\gamma}\exp\left(\lambda_2x\right)\right]\sin\left(\pi\frac{y}{L_y}\right).\tag{21.40}\label{eq:stommel_result_v} \end{align} \]

Die Lösungen dieser Gleichungen sind für

\[ \begin{align} \beta &= 7\cdot 10^{-12}\:\frac{1}{\text{ms}},\\ D &= 4.000\:\text{m},\\ U &= 10\:\frac{\text{m}}{\text{s}},\\ L_x &= 10.000\:\text{km},\\ L_y &= 6.000\:\text{km},\\ \rho_0 &= 1024\:\frac{\text{kg}}{\text{m}^3},\\ \rho_a &= 1,2\:\frac{\text{kg}}{\text{m}^3},\\ C &= 0,001,\\ r &= 5,0\:\frac{\text{kg}}{\text{m}^2\text{s}} \end{align} \]

in Abb. 21.1 dargestellt. Das Stommel-Modell erklärt also die Existenz westlicher Randströme (Golfstrom, Kuroshio).

Im Falle der $f-$Ebene gilt $\lambda_2 = -\lambda_1.$ Setzt man dies in Glg. (21.38) ein, erhält man

\[ \begin{align} \psi_f\left(x,y\right) &= \frac{L_y}{r\pi}C\rho_aU^2\left[\frac{1 - \exp\left(-\lambda_1L_x\right)}{\gamma}\exp\left(\lambda_1x\right) + \frac{\exp\left(\lambda_1L_x\right) - 1}{\gamma}\exp\left(-\lambda_1x\right) + 1\right]\sin\left(\pi\frac{y}{L_y}\right)\nonumber\\ \Rightarrow\psi_f\left(\frac{L_x}{2}+x,y\right) &= \frac{L_y}{r\pi}C\rho_aU^2\left[\frac{\exp\left(\lambda_1\frac{L_x}{2}\right) - \exp\left(-\lambda_1\frac{L_x}{2}\right)}{\gamma}\left(\exp\left(\lambda_1x\right) + \exp\left(-\lambda_1x\right)\right) + 1\right]\sin\left(\pi\frac{y}{L_y}\right). \end{align} \]

Die Lösung ist also in diesem Fall symmetrisch um $x = \frac{L_x}{2}$, es entsteht somit kein westlicher Randstrom. Die Breitenabhängigkeit des Coriolisparameters ist also konstitutiv für die Existenz der Randströme.

21.3.1.2 Munk-Modell

21.4 Thermohaline Zirkulation

Die thermohaline Zirkulation beschreibt die planetarische Zirkulation, die dadurch entsteht, dass sich Teilchen, solange sie keinen diabatischen Einwirkungen ausgesetzt sind, entlang von Isopyknen bewegen. Das Grundkonzept ist dabei wie folgt:

1. Durch diabatische Flüsse (insbesondere Verdunstung und daraus folgende Abkühlung und Salinitätserhöhung) verringert sich die potentielle Dichte einer Wassermasse nahe der Oberfläche.

2. Die Wassermasse sinkt ab (Tiefenwasserbildung) und die Teilchen bewegen sich entlang von Isopyknen durch den Ozean, bis die Isopykne wieder die Oberfläche schneidet.

Konkret findet ein großer Teil der Tiefenwasserbildung im Nordatlantik statt, nachdem sich das über den Golfstrom dorthin transportierte Wasser durch latente und sensible Wärmeflüsse abgekühlt hat. Man bezeichnet diese Zirkulation auch als ocean conveyor belt. Im Nordpazifik findet hingegen keine Tiefenwasserbildung statt, da

• durch die Passatwinde Süßwasser in Form von Niederschlag vom Atlantik zum Pazifik transportiert wird und der Atlantik dadurch salzreicher ist und

• das aus dem sehr salzigen Mittelmeer kommende Wasser die Tiefenwasserbildung im Atlantik verstärkt, was im Pazifik nicht der Fall ist.

Das so im Nordatlantik gebildete Wasser bezeichet man als North Atlantic deep water (NADW). Das schwerste und sich somit am weitesten unten befindende Wasser ist jedoch das Antarctic bottom water (AABW). Dieses Wasser ist nicht nur sehr kalt, sondern durch die Bildung von Eis und das dabei ausfrierende Salz auch sehr salzreich, weshalb es die größte potentielle Dichte hat.

21.5 Ozeanmodelle

Ozeanmodelle simulieren den Zustand des Ozeans, welcher durch die in Abschn. 21.1 beschriebenen Variablen festgelegt wird, analog zu Atmosphärenmodellen, die den Zustand der Atmosphäre simulieren. Gegenüber der Atmosphäre gibt es bei der Ozeanmodellen drei wesentliche Unterschiede:

• Die thermodynamische Zustandsgleichung ist wesentlich komplizierter. Dies ist numerisch allerdings kein großes Problem.

• Im Ozean gibt es keine Phasenübergänge. Dies ist eine wesentliche Vereinfachung, da damit das Problem der Berechnung der Übergangsraten und der latenten Wärme wegfällt undd weniger Tracer advehiert werden müssen.

• Während die Obergrenze der Atmosphäre mehr oder weniger willkürlich festgelegt wird, hat die Wassersäule eine klare vertikale Ausdehnung $[-H,\eta]$, hierbei sind $H$ die Tiefe und $\eta$ die Auslenkung der Wasseroberfläche aufgrund von Strömungen (Tiden und windgetriebene Zirkulation sowie Wellen). Die Schwierigkeiten hierbei sind, dass $H$ Null werden kann (Land), sowie dass $\eta$ zeitabhängig ist.

21.5.1 Prognostische Variablen

Die prognostischen Variablen von Ozeanmodellen sind im Wesentlichen die in Abschn. 21.1 beschriebenen, jedoch mit einigen Besonderheiten. Im Folgenden ist eine Liste typischerweise verwendeter prognostischer Variablen aufgeführt, wobei natürlich auch andere Variablen verwendet werden können, solange die in Abschn. 21.1 aufgelisteten daraus abgeleitet werden können:

• Die Horizontalgeschwindigkeiten und die Oberflächenauslenkung werden direkt als solche verwendet.

• Da in Ozeanmodellen fast immer von der hydrostatischen Approximation Gebrauch gemacht wird, wird die Vertikalgeschwindigkeit im $p-$System oder in Bezug auf eine andere generalisierte Vertikalkoordinate ausgedrückt.

• Anstatt der Temperatur wird die potentielle Temperatur verwendet, da diese ein passiver Tracer ist und eine einfache Kontinuitätsgleichung erfüllt.

• Die Dichte wird direkt als solche verwendet, da sie im Gegensatz zur potentiellen Dichte eine einfache Kontinuitätsgleichung erfüllt.

• Die Salinität wird entweder als spezifische Größe $q_s$ (Masse Salz pro Masse Meerwasser) oder als volumenbezogene Größe $\rho q_s$ (Masse Salz pro Volumen Meerwasser) verwendet.

Alle anderen thermodynamischen Größen können aus der potentiellen Temperatur, der Dichte und der Salinität abgeleitet werden.

21.5.2 Vertikalkoordinaten

Aufgrund der im vorherigen Abschnitt beschriebenen Herausforderung kommt der Vertikalkoordinate in der Ozeanmodellierung eine besondere Bedeutung zu. Allgemein ist eine Vertikalkoordinate $\kappa$ eine bijektive Abbildung

\[ \begin{align} z = z\left(\kappa\right). \end{align} \]

21.5.2.1 z-Koordinaten

Bei z-Koordinaten (oder auch geopotentielle Koordinaten genannt) sind die Koordinatenflächen genau horizontal. Hauptvorteil ist die Abwesenheit von metrischen Termen bei horizontalen Ableitungen. An den Küsten schneiden die Koordinatenflächen den Ozeanboden, sodass die Anzahl der Schichten von den Horizontalkoordinaten abhängt. Aufgrund der Oberflächenauslenkung $\eta$ ist außerdem die Dicke der ersten Schicht zeitabhängig. Die Bathymetrie wird dabei stufenförmig repräsentiert.

21.5.2.2 $\sigma_z$-Koordinaten

Bei $\sigma_z$-Koordinaten werden die Schichten an die Bathymetrie angepasst, sodass sie nicht den Ozeanboden schneiden. In flachen Gewässern werden die Schichten dabei sehr dünn und an Hängen verlaufen die Koordinatenflächen sehr steil. Die mathematische Formulierung ist analog zu den in Abschn. 12.3 beschriebenen geländefolgenden Koordinaten in der Atmosphäre.

21.5.2.3 $\rho$-Koordinaten

Innerhalb des tieferen Ozeans verlaufen die Bewegungen vorwiegend entlang von Isopyknen. Richtet man die Koordinatenflächen an diesen aus, ist dies numerisch vorteilhaft, da dann Vertikalbewegungen fast verschwinden. Dies ist jedoch nur möglich, wenn der Ozean stabil geschichtet ist (sonst ist die potentielle Dichte innerhalb einer Wassersäule nicht eindeutig auf die geometrische Höhe abbildbar), was nicht überall der Fall ist, weshalb diese Koordinaten nicht global anwendbar sind. Dies ist analog zu den in Abschn. 12.2 betrachteten isentropen Koordinaten in der Atmosphäre.

21.5.2.4 $z^\star$-Koordinaten

$z^\star$-Koordinaten („z-Stern-Koordinaten“) sind eine Verallgemeinerung von z-Koordinaten bei denen das Problem, dass die Dicke der obersten Schicht von der Oberflächenauslenkung abhängt, gelöst wird. Die Transformation auf geometrische Koordinaten lautet

\[ \begin{align} z = \eta + z^\star\frac{\eta + H}{H}. \end{align} \]

Die Koordinatenflächen sind hier also nicht mehr genau horizontal, sie schneiden aber weiterhin den Ozeanboden

21.5.2.5 Hybrid-Koordinaten

Alle bisher aufgelisteten Vertikalkoordinaten haben wesentliche Vor- und Nachteile. Hybrid-Koordinaten sind ein Versuch, die Vorteile der verschiedenen Koordinaten miteinenander zu kombinieren. Innerhalb der oberen paar Dekameter, der sogenannten Thermokline, ist der Ozean gut durchmischt und schwach geschichtet. Hier sind z- oder $z^\star$-Koordinaten sinnvoll. Weiter unten wird der Ozean von quasigeostrophischen Bewegungen enlang der Isopyknen (Flächen gleicher potentieller Dichte) dominiert. In Küstenregionen sind hingegen geländefolgende Koordinaten sinnvoll, um die Bathymetire abbilden zu können. Hybride Koordinaten sind dementsprechend eine Mischung von $\rho-$, $\sigma_z-$ und $z^\star-$Koordinaten. Sie geben dem Ozeanmodell HyCOM (Hybrid Coordinate Ocean Model) seinen Namen.

21.6 Seegangsvorhersage

21.6.1 Prognostische Variablen

Seegang ist ein anderes Wort für Wasseroberflächenwellen. Die prognostische Variable, auf die dabei abgezielt wird, ist die Auslenkung der Wasseroberfläche

\[ \begin{align} h = h\left(x, y, t\right) \end{align} \]

von der mittleren Position der Wasseroberfläche.Es kann nicht einfach $h$ als Abweichung vom Geoid definiert werden, da hier noch die dynamische Topographie überlagert wäre. Diese zählt nicht in die mit den Wellen verbundene Auslenkung hinein. Dementsprechend muss die Länge des Mittelungsintervalls gewählt werden. In vielen Situationen ist $h$ keine Funktion der horizontalen Koordinaten, wie zum Beispiel im Falle von Brandung oder bei aufgewühlter See, da in diesen Fällen die Position der Wasseroberfläche nicht mehr eindeutig festgelegt ist. Solche Effekte werden später in Form von Energie dissipierenden Quelltermen berücksichtigt.

Es könnte nun als Satz prognosticher Gleichungen die in Abschn. 13.8.1 hergeleiteten Flachwassergleichungen verwendet werden, eventuell mit einigen halb-empirischen Zusatztermen. Dies hat jedoch mehrere Nachteile:

• Um einen wesentlichen Anteil der Energie des Wellenfeldes zu erfassen, wären ein sehr kleiner Zeitschritt und eine sehr feine räumliche Auflösung erforderlich.

• Die Phasen der Wellen sind in den allermeisten Fällen nicht bekannt und auch nicht sehr relevant. Vielmehr sind Größen wie Richtung und Höhe der Wellen entscheidend.

Daher verwendet man zur Beschreibung von Wasseroberflächenwellen meistens eine Strahlungsübertragungsgleicung (radiative transfer equation (RTE)). Als prognostische Variable wird daher die spektrale Strahldichte

\[ \begin{align} N = N\left(\mathbf{k},\mathbf{r},t\right) \end{align} \]

verwendet. Hierbei sind $\mathbf{r}$ ein zweidimensionaler Ortsvektor und $\mathbf{k}$ ein zweidimensionaler Wellenvektor. Üblicherweise wird $\mathbf{k}$ nicht in kartesischen, sondern in Polarkoordinaten $\left(k,\theta\right)$ angegeben. Als Dispersionsrelation $\omega = \omega\left(k, \theta\right)$ wird diejenige der Wasseroberflächenwellen

\[ \begin{align} \omega^2 &= gk\tanh\left(kD\right) \Rightarrow \omega = \sqrt{gk\tanh\left(kD\right)}, \end{align} \]

verwendet, hierbei ist $D$ die mittlere Wassertiefe (Wassertiefe ohne Wellen). Hieraus folgen

\[ \begin{align} c_\text{ph} &= \frac{\omega}{k} = \sqrt{\frac{g\tanh\left(kD\right)}{k}},\\ c_\text{gr} &= \frac{1}{2\omega}\frac{\partial\omega^2}{\partial k} = \frac{1}{2\omega}\left[g\tanh\left(kD\right) + \frac{gkD}{\cosh^2\left(kD\right)}\right] = \frac{g\tanh\left(kD\right)}{2\omega}\left[1 + \frac{kD}{\sinh\left(kD\right)\cosh\left(kD\right)}\right]\nonumber\\ &= \frac{c_\text{ph}}{2}\left[1 + \frac{2kD}{2\sinh\left(kD\right)\cosh\left(kD\right)}\right] = \frac{c_\text{ph}}{2}\left[1 + \frac{2kD}{\sinh\left(2kD\right)}\right]. \end{align} \]

Die Gruppengeschwindigkeit ist also isotrop, aber nicht homogen,

\[ \begin{align} c_\text{gr} = c_\text{gr}\left(k,\mathbf{r},t\right). \end{align} \]

Die Orts- und Zeitabhängigkeit entsteht dabei über die Orts- und Zeitabhängigkeit von $D$.

21.6.2 Wave action equation

Hieraus kann man eine spektrale Strahlungsflussdichte

\[ \begin{align} N\left(k, \theta,\mathbf{r},t\right)\left(c_\text{gr}\left(k,\mathbf{r},t\right) + \mathbf{v}\left(\mathbf{r}, t\right)\right) \end{align} \]

herleiten, hierbei ist $\mathbf{v} = \mathbf{v}\left(\mathbf{r}, t\right)$ die Stromgeschwindigkeit. Die Strahlungsübertragungsgleichung ist eine Art Kontinuitätsgleichung für $N$, was konzeptionell mit der Energieerhaltung zusammenhängt:

Die Orts- und Zeitabhängigkeit wurde hierbei nicht mehr mitnotiert. Weiterhin wurden fünf zusätzliche Quellterme aufgenommen:

• $S_\text{nl}\left(k, \theta\right)$ ist die Energiequelle, die aus nichtlinearen Effekten resultiert. Dieser Term ist besonders relevant für die Wechselwirkung von Skalen. Ein Beispiel ist die Entstehung langer Wellen (sogenannte Dünung) aus der Windsee (unmittelbar aus Windeinwirkung entstehende Wellen), welche aus kürzeren Wellen besteht, über upscaling.

• $S_\text{ws}\left(k, \theta\right)$ ist die Energiequelle, die aus der Wechselwirkung mit dem Windfeld resultiert. Dies ist die Hauptenergiequelle des Wellenfeldes.

• $S_\text{wc}\left(k, \theta\right)$ ist die Energiequelle, die aus der Gischt (das sogenannte whitecapping) resultiert, dies ist ein dissipativer und somit üblicherweise negativer Term.

• $S_\text{diss}\left(k, \theta\right)$ ist die Energiequelle, die mit der inneren Reibung (Viskosität) zusammenhängt. Da die Viskosität quadratisch stärker auf die kleineren Skalen wirkt, führt dieser Effekt dazu, dass Windsee schneller als Dünung dissipiert wird.

• $S_\text{bd}\left(k, \theta\right)$ ist die Energiequelle, die aus der Wechselwirkung mit der Bathymetrie (bottom drag) zusammenhängt, auch dies ist ein dissipativer und somit üblicherweise negativer Term.

Abhängig von der konkreten Situation können weitere Quellterme aufgenommen werden. Glg. (21.60) bezeichnet man auch als wave action equation.

Die Gesamtenergie des Wellenspektrums $E$ an einem gewissen Ort zu einer gewissen Zeit ist das Integral über das gesamt Spektrum:

\[ \begin{align} E = \int_0^{2\pi}\int_0^\infty N\left(k,\theta\right)dkd\theta.\tag{21.61}\label{eq:wave_spectrum_total_energy} \end{align} \]

Das energiegewichtete spektrale Mittel einer Größe $\psi$ berechnet sich somit über

\[ \begin{align} \newoverline{\psi} \coloneqq \frac{1}{E}\int_0^{2\pi}\int_0^\infty\psi\left(k,\theta\right)N\left(k,\theta\right)dkd\theta.\tag{21.62}\label{eq:wave_spectral_average} \end{align} \]

21.6.3 Diagnostische Variablen

21.6.3.1 Signifikante Wellenhöhe

Die signifikante Wellenhöhe ist eine Art repräsentative Wellenhöhe, für diese gilt

\[ \begin{align} H_s = 4\sqrt{E}. \end{align} \]

21.6.3.2 Mittlere Wellenrichtung

Die mittlere Wellenrichtung $\theta_m$ berechnet sich gemäß Glg. (21.62) als

\[ \begin{align} \theta_m = \arctan2\left(b,a\right) \end{align} \]

mit

\[ \begin{align} a &\coloneqq \frac{1}{E}\int_0^{2\pi}\int_0^\infty\cos\left(\theta\right)N\left(k,\theta\right)dkd\theta,\nonumber\\ b &\coloneqq \frac{1}{E}\int_0^{2\pi}\int_0^\infty\sin\left(\theta\right)N\left(k,\theta\right)dkd\theta. \end{align} \]

21.6.3.3 Mittlere Wellenlänge

Für die mittlere Wellenlänge $L_m$ gilt

\[ \begin{align} L_m = \newoverline{\left(\frac{2\pi}{k}\right)} = 2\pi\newoverline{k^{-1}}. \end{align} \]

21.6.3.4 Mittlere Wellenperiode

Für die mittlere Wellenperiode gibt es verschiedene Möglichkeiten, diese zu berechnen:

\[ \begin{align} T_{m,1} &= \frac{2\pi}{\newoverline{\sigma}},\\ T_{m,2} &= \frac{2\pi}{\sqrt{\newoverline{\sigma^2}}},\\ T_{m,-1} &= \newoverline{\left(\frac{2\pi}{\sigma}\right)} = 2\pi\newoverline{\sigma^{-1}}. \end{align} \]

Hierbei ist $\sigma$ die relativ zum Meeresgrund gemessene Kreisfrequenz.

21.6.4 Seegangsvorhersagemodelle

Das gebräuchlichste Seegangsvorhersagemodell ist Wavewatch III. Modelle, die sich an Glg. (21.60) orientieren, sind sogenannte Seegangsvorhersagemodelle dritter Generation. Sie lösen diese Gleichung auf einem Gitter und sind daher Gitterpunktmodelle, auch wenn sie gelegentlich als Spektralmodelle bezeichnet werden, da ihre prognostische Variable spektrale Bedeutung hat. Sie lösen Glg. (21.60), wobei an jedem Gitterpunkt ein spektrales richtungsabhängiges Gitter im $\left(k,\theta\right)-$Raum aufgespannt wird. Meist wird zusätzlich eine recht große Vielfalt an halb-empirischen Quelltermen $Q_i$ aufgenommen.