Der Baum liefert Zellen und Nachbarschaften. Der Hydrolöser muss zusätzlich festlegen, was erhalten bleibt, welche Wellen über eine Kante laufen und wie aus Zellmittelwerten brauchbare Randzustände entstehen.
1. Welche Größen werden gespeichert?
Für ein ideales kompressibles Gas ohne Viskosität und Wärmeleitung speichern wir Dichte, beide Impulskomponenten und Gesamtenergiedichte. Geschwindigkeit und Druck werden daraus abgeleitet:
q ist die innere Energiedichte; γ>1 der konstante Adiabatenexponent. Zulässige Zustände benötigen ρ>0 und q>0. Negative Drücke sind in diesem Modell kein physikalischer Zustand. Wenn kinetische und gesamte Energie fast gleich groß sind, ist ihre Differenz numerisch empfindlich.
In einem ebenen Test kann ρ eine Volumendichte einer in der dritten Richtung invarianten Strömung bedeuten. In unseren selbstgravitierenden Massenschichten bezeichnet sie eine Flächendichte; Druck und Energie sind dann ebenfalls vertikal integrierte Modellgrößen. Die identische algebraische Form ersetzt keine Entscheidung über die Einheiten.
2. Von der Erhaltungsgleichung zum Kantenfluss
In Richtung einer äußeren Einheitsnormalen n und mit u=v·n lautet der physikalische Fluss:
Die ersten drei Komponenten transportieren Masse und Impuls; der Druck übt zusätzlich eine Kraft auf die Kante aus. Der Energiefluss enthält neben dem Energietransport die Druckarbeit. Integration über eine Zelle K und Anwendung des Divergenzsatzes führen zu:
Gespeichert wird der Zellmittelwert Ū, kein beliebiger Punktwert. Aₖ ist die Fläche, |f| die Kantenlänge. Der numerische Fluss erhält zwei Zustände: einen aus der Zelle und einen vom Nachbarn oder der Randbedingung. Im Projekt wird er am Kantenmittelpunkt ausgewertet.
Konservativität entsteht durch gemeinsame Buchführung: Jede innere Kante besitzt nur einen orientierten Fluss. Er geht mit entgegengesetztem Vorzeichen in die beiden angrenzenden Zellen ein. Summiert man alle Zellbilanzen, heben sich die inneren Beiträge auf. Übrig bleiben Randflüsse und Quellen.
Die erste Identität ist Konsistenz, die zweite die passende Umkehrung der Orientierung. Ein kleiner globaler Bilanzrest beweist allerdings keine scharfe Schockauflösung: Auch eine stark verschmierte Lösung kann konservativ sein.
3. Warum ein Riemann-Löser?
Links und rechts einer Kante liegen meist verschiedene Zustände. Das lokale Riemannproblem fragt nach der Entwicklung dieses Sprungs. In Normalenrichtung besitzt Euler die Wellengeschwindigkeiten u−c, u, u und u+c, mit c=√(γp/ρ). Die äußeren Familien sind akustisch; die mittleren tragen Kontakt- und Scherinformation.
Ein Schock muss die Rankine–Hugoniot-Bedingung erfüllen:
Zusätzlich benötigt eine physikalische schwache Lösung die passende Entropiebedingung. Positiver Druck allein beweist diese nicht. Ein approximativer Riemann-Löser ist also nicht bloß eine naive zentrale Diskretisierung der Differentialgleichung.
Rusanov, Roe und HLLC
Rusanov verwendet eine einzige schnelle Dissipationsskala a=max(|uₗ|+cₗ, |uᵣ|+cᵣ):
Das ist einfach, kann aber langsamere Kontaktwellen stark verschmieren. Roe zerlegt den Sprung in charakteristische Wellen eines geeigneten gemittelten Flussjacobis und dämpft sie unterschiedlich. Roe wurde im historischen Pythoncode als Vergleich untersucht, ist jedoch nicht der heutige native Standard.
Der native Hydrooperator verwendet HLLC. Zwischen geschätzten äußeren Geschwindigkeiten Sₗ und Sᵣ liegt eine Kontaktwelle Sₘ. Mit aₗ=ρₗ(Sₗ−uₗ) und aᵣ=ρᵣ(Sᵣ−uᵣ) ergibt sich aus den Sprungbedingungen:
Diese Größen bestimmen die beiden Sternzustände. Der Zustand, dessen Wellenbereich die ruhende Kante enthält, liefert den Fluss. Unzulässige oder degenerierte Zwischenzustände erfordern einen kontrollierten Ersatzfluss. Die konkret verwendeten Geschwindigkeitsschätzungen und Sicherheitsprüfungen sind Teil der Implementierung, kein allgemeiner Positivitätsbeweis für alle HLLC-Varianten.
4. MUSCL: mehr als ein besserer Fluss
Ein Riemann-Löser kann nur die Zustände verarbeiten, die er erhält. P0 bietet überall in einer Zelle denselben Mittelwert an. Die lineare Rekonstruktion bietet stattdessen für jede konservative Komponente z ein affines Feld an:
xₖ ist der Flächenschwerpunkt. Deshalb integriert sich der lineare Zusatz über das Dreieck zu null: Der gespeicherte Mittelwert bleibt erhalten. Dies ist die räumliche Grundidee von MUSCL. Die vollständige Methode benötigt außerdem Begrenzer, Fluss und Zeitintegration.
Aus den Schwerpunktabständen dⱼ und Mittelwertdifferenzen bⱼ wird eine ungewichtete Least-Squares-Steigung bestimmt:
Bei Rang zwei und affinen Daten reproduziert der unbegrenzte Ansatz die richtige Steigung bis auf Rundung. C++ verwendet eine skalierte Zweispalten-QR; die geometrischen Gewichte werden bis zum Gitterwechsel gespeichert. Bei Rangmangel wird die Steigung auf null gesetzt.
In der Nähe eines Sprungs kann eine Gerade überschießen. Ein komponentenweiser Begrenzer verkleinert deshalb die Steigung anhand der Werte in den Ecksternen. Eine zusätzliche gemeinsame Skalierung aller Komponenten schützt die Dichte und innere Energie der rekonstruierten Zustände. Der Zellmittelwert wird dabei nicht durch einen Floor verändert.
5. Zwei Zeitstufen und eine CFL-Grenze
Die räumliche Diskretisierung definiert einen Operator L. Auf festem Gitter verwendet der Code SSP-RK2:
Gradienten, Begrenzer und Flüsse werden für die zweite Stufe neu berechnet. Für glatte Lösungen ergibt sich Zeitordnung zwei. Die SSP-Darstellung bewahrt geeignete Stabilitätseigenschaften eines Vorwärtsschritts nur dann, wenn dieser sie tatsächlich besitzt.
Kleine Zellen und hohe Wellengeschwindigkeiten verkürzen den globalen Schritt. Die zweite Stufe wird erneut geprüft. Bei unzulässigen Zuständen oder verletzter Grenze wird der Versuch verworfen und die Schrittweite halbiert, bis eine begrenzte Anzahl Versuche erreicht ist. Das ist eine Sicherung, kein Beweis, dass jede Problemsituation durch kleinere Schritte lösbar ist.
6. Welche Tests sagen was aus?
- Sod: Vergleich mit einer planaren exakten Riemannlösung. Frontlage, Kontakt und Verdünnungsfächer müssen getrennt betrachtet werden. Ein gedrehter Sod-Test prüft Orientierungseinflüsse, aber noch keine vollständig zweidimensionale Strömung.
- Glatter Transport: Beobachtete Ordnung ohne Schock. Eine scharfe Schockdarstellung allein weist keine zweite Ordnung nach.
- Kreisstoß und Sedov: Radiale Expansion in 2D. Ein radialer Mittelwert kann Winkelartefakte verdecken; Frontstreuung ist zusätzlich relevant. Die numerische endliche Energiequelle unterscheidet sich von der idealen Sedov-Punktquelle.
- AMR: Nach Gitterwechseln müssen die übertragenen Zellintegrale erhalten bleiben. Restriktion ist flächengewichtet; Prolongation benötigt eine begrenzte Rekonstruktion.
Quellen und Einordnung
- Van Leer (1979): A second-order sequel to Godunov’s method.
- Roe (1981): Approximate Riemann solvers, parameter vectors, and difference schemes.
- Toro, Spruce & Speares (1994): Restoration of the contact surface in the HLL-Riemann solver.
- Barth & Jespersen (1989): Upwind schemes on unstructured meshes.
- Projektvertrag, Herleitungen, weitere Literatur und Zahlen: Arbeit 0.3, Kapitel 6. Die konkrete Kombination wird hier nicht als allgemeiner Entropie- oder Positivitätssatz ausgegeben.