Forschung: Aktueller Forschungsschwerpunkt 2000

(Abteilung: Klimaprozesse)

 

Grobstruktursimulation - eine Methode zur Berechnung turbulenter atmosphärischer Strömungen

 

Die meisten in der Natur und bei technischen Anwendungen auftretenden Strömungen sind turbulent. Insbesondere strömt die untere Atmosphäre turbulent, wenn als turbulente Strömungen in der Atmosphäre diejenigen aperiodischen, dreidimensionalen und stochastischen Bewegungen bezeichnet werden, die die Luft und die mit ihr verbundenen Eigenschaften im Übergang zwischen größerskaliger Advektion und molekularer Diffusion vermischen. Von besonderer Bedeutung in diesem Zusammenhang sind der Austausch von Impuls, Wärme, Wasser und anderer Stoffbei-mengungen zwischen Atmosphäre und dem unteren Rand. Diese werden von der Turbulenz kontrolliert. Die Kenntnis der turbulenten Transportvorgänge in der atmosphärischen Grenzschicht ist daher von großer Bedeutung für das Verständnis des Energie- und Stoffhaushalts der gesamten Atmosphäre.

 

Fortschritte beim Verständnis und bei der Vorhersage turbulenter atmosphärischer Strömungen können nur im Verbund von Experiment, Theorie und numerischer Modellierung erzielt werden. Hierbei kommt der numerischen Simulation turbulenter Strömungen eine immer größere Bedeutung zu, besonders bei Problemen, die komplexe Wechselwirkungen zwischen Strömungsfeld, Wolken- und Strahlungs-prozessen aufweisen. Der größte Vorteil der numerischen Simulation von Strömungen besteht darin, daß unter kontrollierbaren Bedingungen relativ schnell und preiswert Ergebnisse erhalten werden können.

 

Es gibt zwei Methoden die Strukturen und mikroskaligen Prozesse turbulenter hydrodynamischer Vorgänge durch eine numerische Lösung der grundlegenden Gleichungen zu untersuchen: die Direkte Numerische Simulation (DNS) und die Grobstruktursimulation (engl: “Large Eddy Simulation”, LES).

 

Den genauesten Zugang zur Turbulenzsimulation liefert die DNS, weil dabei die ungemittelten Navier-Stokes Gleichungen mit Ausnahme der unvermeidlichen numerischen Diskretisierung ohne zusätzliche Annahmen mit einer vorgegebenen zeitlichen und räumlichen Auflösung der Turbulenzelemente direkt gelöst werden. Sie erfaßt alle Bewegungen, indem alle Skalen der Strömung bis hin zur Dissipationslänge h aufgelöst werden. Bei einer DNS sind zur vollständigen Berechnung einer Strömung sowohl die Mindestgröße des Modellgebiets LD als auch die Maschenweite D des numerischen Gitters festgelegt. Die größten turbulenten Strukturen des simulierten Phänomens müssen kleiner als LD und die kleinsten Wirbel, gegeben durch die Kolmogorov-Länge h, größer als D sein. Die Zahl der notwendigen Gitterpunkte pro Raumrichtung beträgt also LD/D. Dieses Größenverhältnis kann auch in Beziehung zur Reynolds-Zahl (Re-Zahl) gesetzt werden, die das Verhältnis von Trägheits- zu Zähigkeitskräften ausdrückt. Dieser Zusammenhang kann in der Form LD/D ~ Ren angegeben werden, wobei der Exponent n je nach Strömungsart und Turbulenzintensität Werte im Intervall 3/4 £ n £ 1 annimmt. Für eine direkte dreidimensionale Simulation ist die erforderliche Gitterpunktszahl also proportional zu Re3n, so daß eine Verdoppelung dieses Parameters mindestens die fünffache Anzahl an Gitterpunkten erfordert. Hieraus wird deutlich, daß selbst mit den leistungsfähigsten Computern mit massiv-parallelen Rechnerarchitekturen voll turbulente Strömungen hoher Re-Zahl nicht direkt berechnet werden können. Als obere Grenze für die Zahl der Gitterpunkte wird heute der Wert 5123 (~ 1.35×108) angesehen.

 

Atmosphärische Strömungen weisen sehr hohe Re-Zahlen auf (Re ~ 108), was bei einer DNS eine große Anzahl von Gitterpunkten zur Auflösung des Modellgebiets erfordert. In der konvektiven Grenzschicht beispielsweise sind die größten turbulenten Wirbel durch die Mächtigkeit der Grenzschicht beschränkt und weisen daher Abmessungen von etwa 103 Metern auf. Die kleinsten Wirbel hingegen im Dissipationsbereich des Spektrums haben eine Ausdehnung von nur etwa 10-3 Metern. Das Spektrum der turbulenten Wirbel in der konvektiven Grenzschicht umfaßt demnach einen Bereich von einem Millimeter bis zu einigen Kilometern. Eine numerische Integration der Navier-Stokes Gleichungen würde daher 1018 Gitterpunkte umfassen. Dies ist heute (und in absehbarer Zukunft) nicht durchführbar. Selbst ein Teraflop-Computer (1012 Gleitkommaoperationen pro Sekunde) hätte - heutige Lösungsalgorithmen vorausgesetzt - für eine Stunde Simulationszeit einige tausend Jahre zu rechnen.

 

Glücklicherweise weist die Turbulenz in der Grenzschicht folgende Eigenschaft auf: Unter einem energetischen Aspekt betrachtet sind nicht alle Skalen des Turbulenzspektrums gleich wichtig. Vielmehr steckt die meiste kinetische Energie in den größten Skalen der Bewegung, und in der Regel nimmt die Energie als Funktion der Wellenzahl nach einem Potenzgesetz hin zu größeren Wellenzahlen ab. Darüber hinaus sind auch für die Impuls-, Wärme- und Stofftransporte weitgehend die großen Strukturen verantwortlich, während die kleineren Wirbel einen eher dissipativen Charakter aufweisen, indem sie einen Energietransfer von niedrigen zu hohen Wellenzahlen bewirken. Diese Tatsachen bilden die Basis des Konzepts der Grobstruktursimulation. Die Grobstruktursimulationstechnik erfaßt bzw. simuliert nur die großen, energiereichen Wirbel und berücksichtigt die Effekte nicht aufgelöster, subskaliger Dynamik durch eine einfache Parameterisierung. Demnach wird bei der LES-Technik statt der ursprünglichen Navier-Stokes Gleichungen ein Gleichungssystem gelöst, welches über den hochfrequenten Teil der turbulenten Schwankungen gemittelt ist. Man versucht also gar nicht erst die Bewegung der kleinsten Wirbel zu berechnen. Vielmehr simuliert man explizit nur die großen Turbulenzballen und versucht, den Einfluß der kleinen Wirbel angemessen zu berücksichtigen. Die gemittelten Felder der Geschwindigkeit und der anderen Zustandsgrößen variieren viel langsamer als die ursprünglichen, weshalb bei der LES-Technik ein weitaus gröberes Rechengitter als bei der DNS-Technik zur Simulation von Strömungen hoher Turbulenzgrade erforderlich ist. Um die sogenannte Grobstruktur, die die vom Gitter aufgelösten Bewegungsvorgänge enthält, von der Feinstruktur, den sogenannten subskaligen Beiträgen (engl.: “subgrid-scale” (SGS) - Anteil) zu unterscheiden, müssen die in Raum und Zeit kontinuierlichen Navier-Stokes Gleichungen räumlich gefiltert werden. Die Bezeichnung LES ist üblicherweise den Anwendungen vorbehalten, bei denen die Skala der Filteroperation (Filterbreite) innerhalb oder zumindest sehr nahe am Trägheitsunterbereich der dreidimensionalen Turbulenz liegt. In diesem Fall kann für den überwiegenden Teil der Grenzschicht der dissipative Effekt der kleinen Wirbel mit einfachen Turbulenzmodellen erfaßt werden, die auf der Annahme der Gesetzmäßigkeiten des Trägheitsunterbereichs im Turbulenzspektrum basieren. Nur in unmittelbarer Wandnähe und innerhalb der Entrainmentzone am Oberrand der Grenzschicht, wo die energieintensiven Wirbel kleinere Abmessungen als die Filterskala annehmen, wird die Annahme einer auf dem Trägheitsunterbereich basierten Theorie fragwürdig. Da in diesem Höhenbereich nahezu das gesamte Bewegungsspektrum subskalig ist, sind zur Parameterisierung dieser Effekte aufwendigere Feinstrukturmodelle erforderlich. Trotz dieser konzeptionellen Schwäche stellen LES-Modelle die geeignetsten Instrumente dar, um Grundlagenfragen zur Turbulenz zu beantworten, und sie können darüber hinaus dazu verwendet werden, komplexere Schließungsansätze für größerskalige Modelle weiterzuentwickeln, zu testen und zu entwerfen.

 

Am Max-Planck-Institut für Meteorologie in Hamburg wurde ein Grobstruktursimulationsmodell entwickelt, mit dem sich turbulente atmosphärische Strömungen beschreiben lassen. Das verwendete Modell zeichnet sich dadurch aus, daß Phasenübergänge des Wassers bei der Wolkenbildung und -dissipation, die langwellige Strahlungsabkühlung in Wolken und der Einfluß einer großräumigen Vertikalbewegung berücksichtigt werden können. Die Eignung der Grobstruktursimulationsmethode zur realitätsnahen Beschreibung der turbulenten Strukturen in der Atmosphäre wird im folgenden anhand von drei Beispiel demonstriert. Die Anwendungsfälle umfassen unterschiedliche Strömungsregime, bei denen (a) die Heizung von der unteren Berandung (konvektive Grenzschicht während eines Kaltluftausbruches), (b) die Strahlungsabkühlung an der Obergrenze einer Wolke (stratocumulusbedeckte Grenzschicht) und (c) die latente Wärmefreisetzung den wesentlichen Antrieb für die Umschichtung der Wolken liefern.

 

Das erste Beispiel befaßt sich mit der Simulation der konvektiven Grenzschicht während eines Kaltluftausbruches. In Anlehnung an einen während der Meßkampagne ARKTIS 1988 in der Grönlandsee beobachteten Fall wird versucht, die sich unter Bedingungen eines Kaltluftausbruches vollziehende Entwicklung der Grenzschicht nachzuvollziehen. Zur Veranschaulichung der berechneten Wirbelstrukturen sind in Abbildung 1 typische Realisationen der Sekundärzirkulation in Horizontal- und Vertikalschnitten zu verschiedenen Zeitpunkten der Entwicklung dargestellt. Gezeigt werden Isolinien des normierten Flüssigwassergehalts, wobei zur Skalierung das in der betreffenden Schnittebene auftretende Betragsmaximum verwendet wird. Die Ergebnisse der dreidimensionalen Modellrechnungen belegen, daß bei Vorgabe des beobachteten Temperaturunterschiedes von etwa -5 K zwischen Luft und Wasser sich im Modell als Reaktion auf eine kurzzeitige Störung des anfänglich horizontal homogenen Temperaturfeldes bereits nach ca. 15 Minuten Grenzschichtrollen und Wolkenstraßen ausbilden, die bezüglich der Richtung des geostrophischen Windes einen Orientierungswinkel von etwa 500 aufweisen. Nach einer Simulationszeit von 30 Minuten stimmen sowohl die Wolkenhöhe von 550 m als auch die Wellenlänge der entstehenden Wolkenbänder von etwa 1300 m gut mit den Beobachtungen überein. Die Mächtigkeit der Wolkenschicht und der Bedeckungsgrad nehmen mit zunehmendem Abstand von der Eiskante zu. Nach zwei Stunden Integrationszeit hat sich eine 450 m mächtige geschlossene Wolkendecke gebildet, in der die rollenförmige Strömungsstruktur, deren Wirbelachsen fast parallel zur Richtung des geostrophischen Windes ausgerichtet sind, jedoch noch klar erkennbar bleibt. Eine mit dem Modell durchgeführte Parameterstudie zeigt, daß sowohl die Intensität der Sekundärzirkulation als auch der Bedeckungsgrad sehr empfindlich auf die Struktur des großräumigen Feuchtefeldes reagiert. Abbildung 2 zeigt eine Simulation der Grenzschichtentwicklung während eines Kaltluftausbruches bei dem das synoptische Feld im Gegensatz zum Kontrollauf eine lineare Abnahme der spezifischen Feuchte mit der Höhe aufweist. Dargestellt sind die Isolinien der Vertikalgeschwindigkeit (oben) und des Flüssigwassergehalts (unten) im Vertikal- (links) und Horizontalschnitt (rechts) durch das Modellgebiet nach zwei Stunden Integrationszeit. Man erkennt, daß sich in einer trockeneren Atmosphäre, die oberhalb der Inversion eine relative Feuchte von etwa 50% aufweist, eine flachere, durchbrochene Wolkenschicht ausbildet, da die durch das Einmischen trockener Luft bedingte Verdunstung von Wolkenwasser die konvektive Entwicklung abschwächt und den mit Wolken bedeckten Flächenanteil erheblich reduziert.

 

Als weiteres Beispiel wird die marine stratocumulusbedeckte Grenzschicht betrachtet. Die Untersuchung soll zeigen, ob das Grobstruktursimulationsmodell in der Lage ist, die beobachtete Turbulenzstruktur einer wolkenbedeckten Grenzschicht richtig wiederzugeben. Darüber hinaus werden die Ergebnisse der LES-Modellierung einer umfangreichen Fehleranalyse unterzogen, um einige Aspekte der parametrischen und strukturellen Unsicherheit zu beleuchten.

 

Zur Bestimmung der statistischen Signifikanz der vom LES-Modell erzeugten Datensätze werden Ensemble-Läufe des Modells durchgeführt, wobei das Ensemble 21 Realisationen des Prozesses umfaßt. Bei der Festlegung der externen Parameter, die den großräumigen atmosphärischen Zustand charakterisieren, sowie bei der Spezifikation der mittleren Anfangsbedingungen werden Daten verwendet, die während einer Intensivmeßphase (Meßflug RF06 am 13.6.1992) des Atlantic Stratocumulus Transition Experiments (ASTEX) 1992 gewonnen wurden. Die Modellintegrationen beginnen mit der eben beschriebenen Gleichgewichtskonfiguration, der zusätzlich eine zufällige Temperaturstörung kleiner Amplitude überlagert wird. Die einzelnen Realisationen des Ensembles unterscheiden sich nur dadurch, daß verschiedene Sätze von räumlich unkorrelierten Zufallsstörungen zur Initialisierung des Temperaturanomaliefeldes benutzt werden. Die im folgenden gezeigten Profile stellen Momente zweiter Ordnung dar und werden durch Mittelwertbildung über Horizontalebenen gewonnen. Als Vertikalkoordinate wird die mit der Inversionshöhe zi skalierte Höhe verwendet. Abbildung 3 zeigt mittlere Vertikalprofile der (a) totalen vertikalen Geschwindigkeitsvarianz, (b) des totalen Auftriebsflusses, (c) des totalen Gesamtwasserflusses und (d) des Niederschlagsflusses für die verschiedenen Ensembleläufe zum Zeitpunkt t = 9000 s. Das Profil des jeweiligen Ensenblemittelwertes ist durch eine dicke Linie markiert. Aus Abbildung 3 wird deutlich, daß die Profile aller Momente eine große Streuung aufweisen. Bei z /zi = 0.7 beträgt der relative Fehler für alle Momente etwa 40%. Dies deutet darauf hin, daß der instantane horizontale Mittelwert bei der hier verwendeten Modellgebietsgröße nicht geeignet ist, um statistisch hoch signifikante Schätzwerte für die Momente zweiter Ordnung zu liefern. Der Einfluß der zeitlichen Mittelung auf die Güte der abgeleiteten statistischen Kenngrößen wird in Abbildung 4 dargelegt. Dargestellt sind Profile derselben Größen wie in Abbildung 3, jedoch wird - zentriert um t = 9000 s - zusätzlich eine einstündige Mittlelung vorgenommen. Augenscheinlich wird der statistische Stichprobenfehler bei der Bestimmung der Momente durch die zusätzliche Zeitmittelung erheblich reduziert. Dennoch weisen die individuellen Realisationen des Ensembles immer noch eine beträchtliche Abweichung vom Ensemblemittelwert auf; bei z /zi = 0.7 beträgt der relative Fehler für alle Momente etwa 10%. Zur weiteren Reduktion des Stichprobenfehlers müßte demnach eine längeres zeitliches Mittelungsintervall - sofern dies sinnvoll erscheint - oder aber ein größeres Modellgebiet bei den Simulationsrechnungen verwendet werden.

 

Die parametrische Unsicherheit der Grobstruktursimulationsresultate, die sich aus der unvollständigen Kenntnis der Eingabe- und Eingangsgrößen des Modells, d. h. der Anfangs-, Rand- und großräumigen Umgebungsbedingungen, ergibt, wird im folgenden untersucht und quantifiziert. Bei dieser Analyse wird herausgearbeitet, wie die Modellresultate auf kleine Variationen eines ungenau bekannten Faktors bei ansonsten festgehaltenen Eingangsgrößen reagieren. Mit Hilfe einfacher Fehlerfortpflanzungsprozeduren werden die Ergebnisse dieser Parameterstudie dazu verwendet, die Konfidenzintervalle für die berechneten statistischen Kenngrößen als Funktion der relevanten externen Einflußgrößen zu ermitteln. Abbildung 5 zeigt stündliche Mittelwerte der Vertikalprofile des (a) Flüssigwassergehalts, (b) des totalen Auftriebsflusses, (c) des totalen Gesamtwasserflusses, (d) des Niederschlagsflusses, (e) der totalen Geschwindigkeitsvarianz und (f) der totalen turbulenten kinetischen Energie für einen Modellauf, mit dem die während ein Meßtages des ASTEX-Experiments beobachtete Turbulenzstruktur nachgebildet werden soll. Zum Vergleich sind daher mit Rauten markierte Meßwerte dieser Größen in die Abbildung eingetragen. Zusätzlich sind in Abbildung 5 in ausgewählten Höhen die 90%-Konfidenzintervalle, die mit Hilfe der parametrischen Sensitivitätsstudie ermittelt wurden, für die verschiedenen dargestellten Größen eingetragen. Man erkennt, daß mit Ausnahme des Niederschlagsflusses die modellierten thermodynamischen, dynamischen und wolkenphysikalischen Zustandsgrößen der stratocumulusbedeckten Grenzschicht eine recht gute Übereinstimmung mit den aus Meßdaten abgeleiteten Werten aufweisen. Die Unterschiede zwischen den Modellvorhersagen und den Meßdaten sind (mit Ausnahme des Niederschlagsflusses) stets kleiner als die durch die parametrische Modellunsicherheit bewirkte Streuung in den Zielvariablen. Dabei tragen die Unsicherheiten in der Kenntnis der Inversionsstärke des Gesamtwassergehalts und in der langwelligen Strahlungsabkühlung die weitaus größten Anteile zur modellierten Varianz bei. Eine mögliche Erklärung für das Unvermögen des Modells, den Niederschlagsfluß vollständig richtig wiederzugeben, besteht darin, daß der Niederschlagsprozeß einen hochgradig in Raum und Zeit variierenden und intermittierenden Vorgang darstellt, dessen Simulation aufwendigere Initialisierungsprozeduren und wirklichkeitsgetreuere Antriebs-felder für das Modell erfordert, als die, die in dieser Studie verwendet werden.

 

Im dritten Beispiel werden schließlich die dynamischen, wolkenphysikalischen und strahlungsbedingten Vorgänge im Nachlaufbereich eines Kondensstreifens untersucht. Der Modelleinsatz erfolgt in der sogenannten Dispersionsphase im Lebenszyklus eines Kondens-streifens, d. h. in demjenigen Zeitabschnitt, in dem die vom Flugzeug induzierte Turbulenz weitgehend abgeklungen ist und die Emissionen sich mit der Umgebungsluft vermischt und über ein gewisses Volumen verteilt haben. Bei der Festlegung der externen Parameter, die den großräumigen atmosphärischen Grundzustand charakterisieren, werden Bedingungen gewählt, die typisch für das Auftreten langlebiger Kondensstreifen in Reiseflughöhe sind. Dies impliziert einen stabil geschichteten atmosphärischen Grundzustand, der bezüglich der Wasserphase untersättigt jedoch bezüglich der Eisphase leicht übersättigt ist; für die Umgebungstemperatur im Flugniveau wird eine Temperatur von -55o C angenommen.

 

Zur Veranschaulichung der berechneten Wirbelstruktur werden in den Abbildungen 6 typische Realisationen des Strömungsfeldes in vertikalen Schnittebenen senkrecht zur Flugrichtung gezeigt, die die zeitliche Entwicklung des Kondensstreifens über eine Zeitspanne von 30 Minuten verfolgen. Dargestellt sind Isolinien des Eisgehalts sowie die Vektoren des Geschwindigkeitsfeldes (Maximum angegeben). Man erkennt, daß nach zwei Minuten Simulationszeit das Strömungsfeld aus zwei gegenläufig rotierenden Wirbeln besteht, deren Zentren sich am Rande des Kondensstreifens befinden. Die Wirbel mischen feuchte, übersättigte Umgebungsluft in den Kondensstreifen ein, deren überschüssiger Wasserdampf nachfolgend infolge von Phasenflüssen an den bereits existierenden Partikeln zu Eis sublimiert wird. Nach 30 Minuten Simulationszeit nimmt der Kondensstreifen einen Großteil des horizontalen Modellgebiets ein. Die Maxima des Eisgehalts sind im Zentrum des Kondensstreifens zu finden, wobei Werte von 20 mg kg-1 anzutreffen sind. Während der Entwicklung des Kondensstreifens weist das Strömungsfeld eine zunehmend komplexere Struktur auf, die sich aus mehreren Wirbeln zusammensetzt. Die maximale Amplitude des Vertikalgeschwindigkeitsfeldes beträgt etwa 0.1 m s-1 und ist am oberen Rand des Kondensstreifens zu finden. Eine eingehende Analyse ergibt, daß die sich im Kondensstreifen entwickelnden Wirbel im wesentlichen durch die durch Auftrieb produzierte kinetische Energie gespeist werden. Dabei spielt die durch Phasenumwandlung freigesetzte latente Wärme die entscheidende Rolle; hingegen übt die durch differentielle Strahlungserwärmung hervorgerufene kinetische Energieproduktion wegen der geringen optischen Dicke der Kondensstreifen nur einen geringen Einfluß aus. Eine Parameterstudie zeigt, daß langlebige, d. h. sich infolge der Eigendynamik ausbreitende Kondensstreifen, sich nur in einer in bezug auf die Eisphase übersättigten Atmosphäre ausbilden können. Weiterhin kann demonstriert werden, daß neben der relativen Feuchte und der Temperatur der umgebenden Atmosphäre die Dichteschichtung des Mediums die entscheidende Einflußgröße darstellt, die das Ausbreitungsverhalten des Kondensstreifens maßgeblich bestimmt. Die Baroklinität der Atmosphäre und damit die Scherung des Grundstromprofils ist nur von untergeordneter Bedeutung, da sie die Geschwindigkeit der geometrischen Querschnittsänderung des Kondensstreifens nur geringfügig modifiziert.

 

Die hier vorgestellten Simulationsrechnungen belegen, daß mit Grobstruktursimulationsmodellen ein leistungsfähiges Werkzeug zur Erforschung atmosphärischer Turbulenz und ihrer Eigenschaften zur Verfügung steht. Mit der Weiterentwicklung der Leistungsfähigkeit moderner Computer ist zu erwarten, daß in der Zukunft mit der LES-Technik in Bereiche der Turbulenzdynamik vorgedrungen werden kann, die bisher nur ansatzweise bearbeitet werden konnten. So werden zukünftig beispielsweise Langzeitsimulationen realisierbar sein, wo bei sich ändernden Randbedingungen die tageszeitliche Entwicklung einer wolkenbedeckten Grenzschicht nachvollzogen wird.


Abbildungen: