Die schwierigsten Gleichungen in der CFD — CMPS-Strategien für drucklose Euler-Gleichungen (verdünnte Euler-Partikel)

1. Einleitung

Unter allen Gleichungssystemen der numerischen Strömungsmechanik zählen die drucklosen Euler-Gleichungen zu den am schwierigsten zu lösenden. Sie beschreiben die Bewegung einer verdünnten, dispersen Phase – typischerweise Feststoff- oder Flüssigkeitspartikel –, die als Kontinuum, jedoch ohne inneren Druckterm, behandelt wird. Durch das Fehlen von Druck entfällt der in der Gasdynamik vorhandene natürliche Wellenausbreitungsmechanismus, wodurch das System schwach hyperbolisch oder sogar nicht-hyperbolisch wird.

Das Fehlen von Hyperbolizität führt zu erheblichen numerischen Herausforderungen. Da die Gleichungen keine reelle Schallgeschwindigkeit zulassen, können herkömmliche Riemann-Solver keine korrekten Wellengeschwindigkeiten definieren. Selbst kleine Störungen können zu unphysikalischen Lösungen wie δ-Stoßwellen , Vakuumregionen oder undefinierten Kontaktdiskontinuitäten führen. Diese Verhaltensweisen machen das System extrem empfindlich gegenüber Diskretisierungsfehlern und häufig instabil unter Standard-CFD-Verfahren.

Im Verdünnungsgrenzfall, wenn der Partikelvolumenanteil und die Dichte gegen null streben, wird die Situation noch komplexer. Die Erhaltungsgleichungen müssen weiterhin gelten, doch die konvektiven Flüsse dominieren vollständig und können zu singulären Strukturen kollabieren. Daher gelten die drucklosen Euler-Gleichungen oft als die schwierigsten Gleichungen in der CFD .

CMPS begegnet diesen Problemen mit einem vollständig gekoppelten, vollständig impliziten Rahmenwerk , das eine kontrollierte künstliche Kompressibilität, regularisierte Flüsse und exakte Jacobi-Matrizen kombiniert und so ein nahezu schlecht gestelltes System in ein numerisch und physikalisch konsistentes Modell für den Transport verdünnter Partikel verwandelt.

2. Die mathematische Schwierigkeit der drucklosen Euler-Gleichungen

Das drucklose Euler-System im Kern in einer Raumdimension (ohne den Volumenanteilsterm) ist

$$ \frac{\partial \rho_p}{\partial t} + \frac{\partial (\rho_p u_p)}{\partial x} = 0, \qquad \frac{\partial (\rho_p u_p)}{\partial t} + \frac{\partial (\rho_p u_p^2)}{\partial x} = 0. $$

Die Fluss-Jacobi-Matrix dieses Systems ist

$$ A = \frac{\partial \mathbf{F}} {\partial \mathbf{U}} = \begin{bmatrix} 0 & 1 \\ -u_p^2 & 2u_p \end{bmatrix}, \qquad mit Eigenwerten \; \lambda_1 = \lambda_2 = u_p. $$

Beide Eigenwerte sind identisch, und die Matrix besitzt nur einen unabhängigen Eigenvektor, weshalb sie nicht diagonalisiert werden kann. Das System ist daher nicht streng hyperbolisch . Alle Charakteristika reduzieren sich auf die Partikelgeschwindigkeit \(u_p\), wodurch Wellentrennung aufgehoben und die Grundlage klassischer Upwind-CFD-Methoden untergraben wird. Die Folge sind bekannte Pathologien: δ-Stoßwellen , spontane Vakuumbildung und extreme Empfindlichkeit gegenüber Diskretisierung.

3. Wiederherstellung der Hyperbolizität in CMPS

CMPS stellt eine funktionsfähige Wellenstruktur wieder her, indem ein kleiner, kontrollierter künstlicher Druck für die Teilchenphase eingeführt wird.

$$ p_p = \varepsilon_p \; \rho_p^{\gamma_p}, $$

was zu einer endlichen Teilchen-Wellen-Geschwindigkeit führt

$$ c_d = \sqrt{\varepsilon_p \, \gamma_p \, \rho_p^{\gamma_p - 1}} . $$

Dadurch wird das entartete Spektrum in zwei verschiedene reelle Eigenwerte umgewandelt.

$$ \lambda_{1,2} = u_p \pm c_d, $$

Die numerische Hyperbolizität wird wiederhergestellt und stabile Aufwindflüsse ermöglicht. Die Kompressibilität ist lokalisiert und begrenzt : Eine Reibungsschwelle in der lokalen Konzentration stellt sicher, dass der künstliche Druck im wirklich verdünnten Bereich verschwindet, während ein Mach-basierter Limiter \(\varepsilon_p\) begrenzt, um unphysikalische Steifigkeit in langsamen Strömungen zu vermeiden.

3.1 Berechnung von \(\varepsilon_p\): Aktivierung und Grenzwert

Der künstliche Druck muss in sehr verdünnten Bereichen inaktiv und andernorts begrenzt sein. CMPS erreicht dies durch zwei Steuerungsmechanismen:

• Aktivierung durch Konzentration: Sei \(\alpha = \rho_p/\rho_{p,mat}\) der lokale Volumenanteil der Partikel (\(\rho_{p,mat}\) ist die Materialdichte). Definiere eine Reibungsschwelle \(\alpha_{fr}\) und eine maximale Packungsdichte \(\alpha_{max}\). Der Flächenmittelwert \(\alpha_{face} = (\alpha_L + \alpha_R)/2\) wird auf \([0,\alpha_{max}]\) begrenzt. Falls \(\alpha_{face} < \alpha_{fr}\), setze \(\varepsilon_p^{eff} = 0\), und die Partikelphase bleibt drucklos.
• Mach-basierter Begrenzer: Um die Pseudo-Schallgeschwindigkeit durch die lokale Teilchengeschwindigkeit zu begrenzen, wird Folgendes erzwungen: $$ c_d \le \frac{|\mathbf{u}_p|}{M_{min}} . $$ Unter Verwendung von \(c_d^2 = \varepsilon_p \, \gamma_p \, \rho_p^{\gamma_p-1}\) ergibt sich folgender effektiver Koeffizient: $$ \varepsilon_p^{eff} = \min\!\left( \varepsilon_p,\; \frac{|\mathbf{u}_p|^2}{\gamma_p\, \rho_p^{\gamma_p-1}\, M_{min}^2} \right). $$ Dies garantiert, dass die künstliche Wellengeschwindigkeit niemals die konvektive Skala dominiert.

Diese Regeln werden separat auf die linken/rechten Zustände an jeder Fläche angewendet, wodurch sich \(p_{p,L/R} = \varepsilon_{p,L/R}^{eff} \rho_{p,L/R}^{\gamma_p}\) und \(c_{d,L/R} = \sqrt{\varepsilon_{p,L/R}^{eff} \gamma_p \rho_{p,L/R}^{\gamma_p-1}} \) ergeben. In stark verdünnten Zonen fällt die Methode konstruktionsbedingt in den drucklosen Grenzfall zurück.

3.2 Bewertung des regularisierten Flusses

Mithilfe der künstlichen Wellengeschwindigkeiten wird der Flächenfluss mit einer robusten Formel vom Rusanov/HLL-Typ berechnet.

$$ \mathbf{F}_{\text{face}} = \tfrac{1}{2}(\mathbf{F}_L + \mathbf{F}_R) - \tfrac{1}{2}\, \lambda_{\max} (\mathbf{U}_R - \mathbf{U}_L), \quad \lambda_{\max} = \max\big( |u_L| + c_{d,L},\, |u_R| + c_{d,R} \big). $$

Dadurch wird die Dissipation an lokale Gradienten angepasst und geht stetig in den drucklosen Grenzfall über, in dem der künstliche Druck keine Wirkung mehr hat. In glatten Bereichen bleibt das Verfahren diffusionsarm; in der Nähe scharfer Gradienten liefert es die notwendige Dissipation, um δ-Stoßwellen und negative Dichten zu verhindern.

3.3 Kopplung und Erhaltung

Die Teilchengleichungen bleiben über den Impuls- und Wärmeaustausch zwischen den Phasen mit gleich großen und entgegengesetzten Quelltermen konsistent an die Gasphase gekoppelt, wodurch die strikte globale Erhaltung gewährleistet wird. Der künstliche Druck beeinflusst lediglich die Berechnung des Teilchenphasenflusses und stört nicht das Erhaltungsgleichgewicht zwischen den Phasen.

3.4 Implizite Stabilisierung und exakte Jacobische Matrizen

Alle Flüsse und Quellen werden linearisiert und in einem vollständig gekoppelten impliziten System unter Verwendung exakter, durch automatische Differentiation gewonnener Jacobimatrizen gelöst. Dies unterdrückt kurzwellige Instabilitäten und erhält die Newton-Konvergenz auch bei starker Phasenreibung und extremen Dichteverhältnissen.