Numerische Strömungsberechnung
Stefan Hickel - Angewandte Strömungssimulation 2
Parameter und Kennzahlen
Gleichungssystem
Turbulenzmodell
Randbedingungen
Diskrete Operatoren
Lösungsalgorithmen
Rechengitter
Programmierung
Berechnung
Visualisierung
Validierung
Physikalische Modellierung
Mathematische Modellierung
Numerische Modellierung
Lösung
Auswertung
CFX – Pre -> Parameterdatei
ICEM CFD -> Gitter
CFX – Solver -> Ergebnisdatei
CFX – Post -> Bilder und Erkenntnis
Finite Volumen Methode
Stefan Hickel - Angewandte Strömungssimulation 4
• Das Rechengebiet wird in nicht überlappende Bereiche := finite Volumina ( FV ) unterteilt.
• Jedes dieser FV stellt ein Kontrollvolumen ( KV ) dar, für das ein Mittelwert bilanziert wird.
• Jedem FV wird ein Knoten zugeordnet, an welchem die diskreten Werte gespeichert werden.
Zellknoten
Finites Volumen
Finite Volumen Methode
Stefan Hickel - Angewandte Strömungssimulation 5
• Die Erhaltungsgleichungen werden über die FV integriert
• Volumenintegrale können mit dem Satz von Gauss in Oberflächenintegrale umgewandelt werden.
∂∂tϕ =∇•Ψ →
∂∂t
ϕ dVV∫ = ∇•ΨdV
V∫
Erhaltungsgröße Fluss Ψ = Ψ( φ )
∇•ΨdVV∫ = n•Ψ
SV
∫ dS = ΨSrechts
∫ dS − ΨSlinks
∫ dS +−!
= Frechts −Flinks +−!
Finite Volumen Methode
Stefan Hickel - Angewandte Strömungssimulation 6
• Beispiel: Impulserhaltung in integraler Form
• Es werden Flüsse über die FV-Oberfläche S bilanziert!
∂∂t
ρu dVV∫ + ρuu
SV∫ ⋅ n dS = − p dS
SV∫ + τ ⋅ n
SV∫ dS + ρFV dVV∫
Zellknoten
Finites Volumen
Fluss
Finite Volumen Methode
Stefan Hickel - Angewandte Strömungssimulation 7
• Approximationsvorschriften werden zur numerischen Auswertung der Flüsse durch die Grenzflächen der Kontrollvolumina benötigt.
• Zur Approximation der Feldgrößen φ und Ψ an anderen Punkten als den Knotenpunkten (wo ja die Werte gespeichert sind) werden Interpolationsvorschriften verwendet.
Zellknoten
Finites Volumen
Fluss
Finite Volumen Methode
Stefan Hickel - Angewandte Strömungssimulation 8
• Approximationsvorschriften werden zur numerischen Auswertung der Flüsse durch die Grenzflächen der Kontrollvolumina benötigt.
• Zur Approximation der Feldgrößen φ und Ψ an anderen Punkten als den Knotenpunkten (wo ja die Werte gespeichert sind) werden Interpolationsvorschriften verwendet.
• Die diskretisierten Erhaltungsgleichungen werden für jedes FV ausgewertet.
• Zusammen mit geeigneten Randbedingungen entsteht ein algebraisches Gleichungssystem, welches numerisch gelöst wird.
Finite Volumen Methode
Stefan Hickel - Angewandte Strömungssimulation 9
Schwerpunkt dieser Vorlesung ist die notwendige zweifache Näherung: 1. Quadratur: Der Wert der Integrale der Flüsse Ψ über die FV-
Oberflächen wird angenähert durch die diskreten Werte der Variablen an einem oder mehreren Punkten auf den Zellwänden.
2. Interpolation: Die Werte von φ an den Zellwänden werden angenähert als Funktion der Werte von φ an den Zellknoten (FV Mitte).
Kompass Notation
Stefan Hickel - Angewandte Strömungssimulation 10
P EE E W
S SW SE
NE N NW
nw ne
se sw
e w
n
s
Oberflächenintegrale
Stefan Hickel - Angewandte Strömungssimulation 11
Mittelpunktsregel (2ter Ordnung):
• Das Integral wird angenähert durch das Produkt des Integranden am Mittelpunkt der Zellfläche mit der Größe der Zellfläche
• Der Fehler einer numerischen Vorschrift kann durch Taylorreihenentwicklung analysiert werden.
Fe = ΨdSSe
∫ =ΨeAe
Fe ≈ ΨeAe
ΨeAe
Ae = 1dSSe
∫
Ψ
x0-h/2 x0 x x0-h/2
:= Mittelwert Ψe =1Ae
ΨdSSe
∫
Tafelanschrieb:
Oberflächenintegrale
Stefan Hickel - Angewandte Strömungssimulation 12
Taylor-Analyse der Mittelpunktsregel
Ψ
x0-h/2 x0 x x0-h/2
Oberflächenintegrale
Stefan Hickel - Angewandte Strömungssimulation 13
• Die Ordnung bestimmt die Rate der Gitterkonvergenz bei hinreichend glatter Lösung
log Fehler
Anzahl Gitterpunkte
1. Ordnung
2. Ordnung
geringe Auflösung hohe Auflösung große Gitterweite kleine Gitterweite
4. Ordnung
Oberflächenintegrale
Stefan Hickel - Angewandte Strömungssimulation 14
Trapezregel (2ter Ordnung): Simpson Regel (4ter Ordnung):
Ψ
x0-h/2 x0 x x0-h/2
Fe = ΨdSSe
∫ ≈ AeΨne +Ψ se
2
Fe = ΨdSSe
∫ ≈ AeΨne + 4Ψe +Ψ se
6
Volumenintegrale
Stefan Hickel - Angewandte Strömungssimulation 15
• Üblich ist Integration 2ter Ordnung
• Diese Vorschrift ist exakt, wenn q konstant ist oder in V eine lineare Funktion vom Ort ist.
QP = qdVV∫ = qVP ≈ qPVP
Interpolation
Stefan Hickel - Angewandte Strömungssimulation 16
• Für die Auswertung der im Impulssatz auftretenden Integrale wird Ψ an den Zellwänden benötigt.
• Diese werden durch räumliche INTERPOLATION bestimmt Ψe = f (ΨW ,ΨP ,ΨE ,ΨEE ,...)
Interpolation
Stefan Hickel - Angewandte Strömungssimulation 17
Upwind Interpolation (Upwind Differencing Scheme - UDS):
Ψe =ΨP , falls (u ⋅ n)e > 0
ΨE , falls (u ⋅ n)e < 0
#$%
&%
P E W e
Interpolation
Stefan Hickel - Angewandte Strömungssimulation 18
Lineare Interpolation (Central Differencing Scheme - CDS):
Ψe =ΨEλe +ΨP (1−λe )
P E W e
λe =xe − xPxE − xP
Abbruchfehler – Numerische Diffusion
Stefan Hickel - Angewandte Strömungssimulation 19
• Der Fehler einer numerischen Vorschrift kann durch Taylorreihenentwicklung analysiert werden
• Abbruchfehler des UDS Verfahrens
• EUDS ist proportional zur Gitterweite (xe-xP) , d.h. es handelt sich um ein Verfahren 1. Ordnung.
Ψe =ΨP + (xe − xP )∂Ψ∂x
!
"#
$
%&P
+(xe − xP )
2
2∂2Ψ
∂x2!
"#
$
%&P
+ ...
EUDS ∝ (xe − xP )∂Ψ∂x
#
$%
&
'(P
(plus Terme höherer Ordnung)
Abbruchfehler – Numerische Diffusion
Stefan Hickel - Angewandte Strömungssimulation 20
• Beispielbetrachtung des Abbruchfehlers für UDS anhand der Advektionsgleichung (Transportgleichung):
• Wir nehmen oBdA an, dass U positiv ist.
• Welchen Einfluss hat der Diskretisierungsfehler auf die ursprünglich zu lösende Differentialgleichung?
∂ϕ∂t
+U ∂ϕ∂x
= 0
Abbruchfehler – Numerische Diffusion
Stefan Hickel - Angewandte Strömungssimulation 21
• Die Euler- / Upwind-Diskretisierung für die Ableitungen lautet
• Entwicklung von φ in eine zeitliche und eine räumliche Taylorreihe
ϕin+1 =ϕi
n +Δt ∂ϕ∂t
+Δt2
2∂2ϕ∂t2
+…
ϕi−1n =ϕi
n −Δx ∂ϕ∂x
+Δx2
2∂2ϕ∂x2
+…
n – Zeitpunkt
i – Gitterpunkt
∂ϕ∂t
≈ϕi
n+1 −ϕin
Δt∂ϕ∂x
≈ϕi
n −ϕi−1n
Δx
Abbruchfehler – Numerische Diffusion
Stefan Hickel - Angewandte Strömungssimulation 22
• Entwicklung von φ in eine zeitliche und eine räumliche Taylorreihe
• Mit der Upwind-Vorschrift für U>0 erhalten wir für die Ableitungen:
Fehler numerische Vorschrift was wir eigentlich berechnen wollten
n – Zeitpunkt
i – Gitterpunkt
∂ϕ∂t
=ϕi
n+1 −ϕin
Δt−Δt2∂2ϕ∂t2!
∂ϕ∂x
=ϕi
n −ϕi−1n
Δx+Δx2∂2ϕ∂x2!
Abbruchfehler – Numerische Diffusion
Stefan Hickel - Angewandte Strömungssimulation 23
• Einsetzen in die Advektionsgleichung liefert
• Die Vernachlässigung des Abbruchfehlers im Lösungsprozess liefert die diskretisierte Gleichung:
∂ϕ∂t
+U ∂ϕ∂x
= 0 ⇒
ϕin+1 −ϕi
n
Δt+U ϕi
n −ϕi−1n
Δx−Δt2∂2ϕ∂t2
+U Δx2∂2ϕ∂x2
Abbruchfehler! "### $###
= 0
ϕin+1 −ϕi
n
Δt+U ϕi
n −ϕi−1n
Δx= 0
Abbruchfehler – Numerische Diffusion
Stefan Hickel - Angewandte Strömungssimulation 24
• In Wirklichkeit berechnen wir die Lösung für die Modifizierte Differentialgleichung
• Weitere Umformung ergibt:
• Der Abbruchfehler der hier verwendeten numerischen Vorschrift hat den selben Charakter wie ein Diffusionsterm und wird daher auch “numerische Diffusion” genannt.
∂ϕ∂t
+U ∂ϕ∂x
+Δt2∂2ϕ∂t2
−U Δx2∂2ϕ∂x2
= 0
∂ϕ∂t
+U ∂ϕ∂x
+U 2Δt2
−UΔx2
$
%&
'
()∂2ϕ∂x2
= 0
Abbruchfehler – Numerische Diffusion
Stefan Hickel - Angewandte Strömungssimulation 25
• Der Abbruchfehler kommt allein durch das numerische Diskretisierungsschema zustande und hat im allgemeinen nichts mit der Gleichung zu tun. Daher ist es wichtig, die Diskretisierung passend zur zu lösenden Differentialgleichung zu wählen.
• Der Abbruchfehler des Upwind-Verfahrens hat die Form einer künstlichen Diffusion.
∂ϕ∂t
+U ∂ϕ∂x
+U 2Δt2
−UΔx2
$
%&
'
()∂2ϕ∂x2
= 0
Abbruchfehler = numerische Diffusion
...−νN∂2ϕ∂x2
Abbruchfehler – Numerische Diffusion
Stefan Hickel - Angewandte Strömungssimulation 26
• Für eine Diffusionsgleichung, die bereits einen Term 2. Ableitung im Raum besitzt, wird bei entsprechender Vorgehensweise ebenfalls zur Grundgleichung eine numerische Dissipation hinzukommen, die zu einer effektiven Verringerung der Reynoldszahl im numerischen Verfahren führt.
• Die effektiv berechnete Lösung kann als die einer geringeren
Reynoldszahl entsprechende Strömung interpretiert werden.
...−ν ∂2u∂x2
−νN∂2u∂x2
Re = uLν
→ Reeffektiv =uL
ν +νN
Abbruchfehler – Numerische Diffusion
Stefan Hickel - Angewandte Strömungssimulation 27
• Numerische Diffusionskoeffizient (Viskosität)des UDS Verfahrens
• Ein numerisches Verfahren ist instabil, wenn die numerische Viskosität negativ ist.
- Stabilitätskriterium für maximal erlaubten Zeitschritt Δt ≤ Δx / U - Max. Genauigkeit bei max. Zeitschritt, Stabilität vorausgesetzt.
Courant-Friedrichs-Lewy Zahl
νN =UΔx2
−U 2Δt2
CFL = Δt UΔx
∂ϕ∂t
+U ∂ϕ∂x
+U 2Δt2
−UΔx2
$
%&
'
()∂2ϕ∂x2
= 0
Abbruchfehler – Numerische Diffusion
Stefan Hickel - Angewandte Strömungssimulation 28
Beispiel: Unterexpandierter Freistrahl mit CFX gerechnet
UDS 1ter Ordnung
„High Resolution“ Verfahren
Abbruchfehler – Numerische Diffusion
Stefan Hickel - Angewandte Strömungssimulation 29
Lineare Interpolation (Central Differencing Scheme - CDS): • Taylorreihenentwicklung
• Einsetzen und Umsortieren ergibt
Ψe =ΨEλe +ΨP (1−λe ) λe =xe − xPxE − xP
ΨE =ΨP + (xE − xP )∂Ψ∂x
$
%&
'
()P
+(xE − xP )
2
2∂2Ψ
∂x2$
%&
'
()P
+ ...
Ψe =ΨP + (xe − xP )∂Ψ∂x
$
%&
'
()P
+(xe − xP )
2
2∂2Ψ
∂x2$
%&
'
()P
+ ...
Ψe =ΨEλe +ΨP (1−λe )−(xe − xP )(xE − xe )
2∂2Ψ
∂x2$
%&
'
()P
+ ...
Abbruchfehler – Numerische Diffusion
Stefan Hickel - Angewandte Strömungssimulation 30
Lineare Interpolation (Central Differencing Scheme - CDS): • Abbruchfehler
• ECDS ist näherungsweise proportional zum Quadrat der Gitterweite (xe-xP)(xE-xe) , d.h. es handelt sich um ein Verfahren 2. Ordnung.
Ψe =ΨEλe +ΨP (1−λe )−(xe − xP )(xE − xe )
2∂2Ψ
∂x2$
%&
'
()P
+ ...
ECDS ∝(xe − xP )(xE − xe )
2∂2Ψ
∂x2!
"#
$
%&P
Numerische Diffusion und Dispersion
Stefan Hickel - Angewandte Strömungssimulation 31
UDS diffusiver Fehler
CDS dispersiver Fehler
Tafelanschrieb:
Abbruchfehler – Numerische Diffusion
Stefan Hickel - Angewandte Strömungssimulation 32
Gradienten an Zellflächen
Stefan Hickel - Angewandte Strömungssimulation 33
• Für die Auswertung der diffusiven Terme der zellgemittelten Erhaltungsgleichungen werden auch die Gradienten von φ an den Grenzen des FV benötigt. z.B. Reibung in Impulsgleichung:
• Ansatz entsprechend zentraler Differenzen (CDS) da dies dem ungerichteten Charakter der Diffusion am besten entspricht.
∂∂t
ρu dVV∫ + ρuu
SV∫ ⋅ n dS = − p dS
SV∫ + τ ⋅ n
SV∫ dS + ρFV dVV∫
∂ϕ∂x"
#$
%
&'e
≈ϕE −ϕP
xE − xP
Gradienten an Zellflächen
Stefan Hickel - Angewandte Strömungssimulation 34
• Aus der Taylor-Reihenentwicklung folgt: mit dem Abbruchfehler:
• Für nicht gleichförmige Gitter ist das Verfahren von
1. Ordnung, der dominante Fehlerterm ist proportional zu Δx
• Für gleichförmige Gitter wird die Fehlerordnung um eine Stufe erhöht, wir erhalten ein Verfahren 2. Ordnung!
∂Ψ∂x
#
$%
&
'(e
=ψE −ψP
xE − xP+E
E∝−xE − xe( )2 − xe − xP( )2
2 xE − xP( )∂2ψ∂x2!
"#
$
%&e
+xE − xe( )3 + xe − xP( )3
6 xE − xP( )∂3ψ∂x3!
"#
$
%&e
+…
~ Δx ~ Δx2
Gradienten an Zellflächen
Stefan Hickel - Angewandte Strömungssimulation 35
Was passiert bei der Verfeinerung eines nicht gleichförmigen Gitters?
1. Abstände zwischen den Zellen in gleich große Abschnitte teilen:
• Der erste Fehlerterm ist nur an den alten, aber nicht an den neuen Gitterpunkten vorhanden.
• Der globale Abbruchfehler wird nur geringfügig weniger reduziert, als in einem Verfahren 2. Ordnung.
grob verfeinert
Gradienten an Zellflächen
Stefan Hickel - Angewandte Strömungssimulation 36
Was passiert bei der Verfeinerung eines nicht gleichförmigen Gitters?
2. Abstände zwischen den Zellen in Abschnitte mit gleicher Streckung teilen:
• Die Streckung nimmt ab. • Der erste Abbruchfehlerterm verringert sich schneller als der
zweite Abbruchfehlerterm, wir erhalten somit ein Verfahren 2. Ordnung.
grob verfeinert
Gradienten an Zellflächen
Stefan Hickel - Angewandte Strömungssimulation 37
Was passiert bei der Verfeinerung eines nicht gleichförmigen Gitters?
• Systematische Gitterverfeinerung auf nicht gleichförmigen Gittern kann die selbe Rate der Reduzierung des Abbruchfehlers ergeben wie die Verfeinerung auf gleichförmigen Gittern !
• Gitter sollten möglichst glatt sein, bzw. gleichmäßige Streckung haben.
Anmerkungen zum Abbruchfehler
Stefan Hickel - Angewandte Strömungssimulation 38
Schemata höherer Ordnung • Versprechen genauere Lösungen bei gleichem Gitter • Führen zu sehr großen Differenzensternen
• Einfach für explizite Methoden auf strukturiertem Gitter • Aufwändig/teuer für implizite Methoden und unstrukturierte
Gitter
P
Anmerkungen zum Abbruchfehler
Stefan Hickel - Angewandte Strömungssimulation 39
• Ein Verfahren mit hoher Fehlerordnung bedeutet nicht automatisch auch kleiner Fehler!
• Die Angabe bezieht sich auf die Konvergenzrate bei Gitterverfeinerung:
log E
Anzahl Gitterpunkte
1. Ordnung Upwind
2. Ordnung Zentrale Differenz
geringe Auflösung hohe Auflösung große Gitterweite kleine Gitterweite
ANSYS-CFX
Stefan Hickel - Angewandte Strömungssimulation 40
Interpolation in ANSYS-CFX
• 1st order Upwind Differencing Scheme (UDS1) • 2nd order Upwind Differencing Scheme (UDS2)
• 2nd order Central Differencing Scheme (CDS)
• 1st-2nd order blend factor (UDS1 <-> UDS2, 0 ≤ β ≤1)
• „High-Resolution“ Scheme (automatische lokale Anpassung von β)
Instationäre Probleme
Stefan Hickel - Angewandte Strömungssimulation 42
• Bei instationären Problemen muss eine 4. Dimension diskretisiert werden -> ZEIT
• Instationäre Probleme sind parabolisch in der Zeit (kein „Rückwärts-Einfluss“)
• Instationäre Probleme sind Anfangsrandwertprobleme, also Anfangswertproblem (AWP) und Randwertproblem (RWP)
• Für die Zeitdiskretisierung können analoge Techniken wie für die Ortsdiskretisierung verwendet werden
Instationäre Probleme
Stefan Hickel - Angewandte Strömungssimulation 43
Der Einfachheit halber betrachten wir eine gewöhnliche Differentialgleichung
Aufgabe: • Finden der Lösung φ zum Zeitpunkt tn+1 einen kurzen Moment Δt
nach tn =t0+nΔt (oBdA)
• Um die Lösung bei tn zu berechnen wird iterativ die Lösung bei tn-1=tn-Δt als Startwert verwendet. Anfangswert bei t0 ist gegeben.
Lösung: • Integration
und umstellen
dϕ(t)dt
= f (t,ϕ(t)) ; ϕ(t0 ) =ϕ 0
dϕdttn
tn+1
∫ dt =ϕ n+1 −ϕ n = f (t,ϕ(t)tn
tn+1
∫ )dt
ϕ n+1 =ϕ n + f (t,ϕ(t)tn
tn+1
∫ )dt
Zeitintegration
Stefan Hickel - Angewandte Strömungssimulation 44
Numerische Lösung des Integrals
Euler vorwärts (explizit)
Euler rückwärts (implizit)
(1. Ordnung)
(1. Ordnung)
f
t0 t0+Δt t
f
t0 t0+Δt t
ϕ n+1 =ϕ n + f (tn,ϕn )Δt
ϕ n+1 =ϕ n + f (tn+1,ϕn+1)Δt
Zeitintegration
Stefan Hickel - Angewandte Strömungssimulation 45
Mittelpunktsregel (implizit)
Trapezregel (implizit)
(2. Ordnung)
(2. Ordnung)
f
t0 t0+Δt t
f
t0 t0+Δt t
ϕ n+1 =ϕ n + f (tn+1 2,ϕn+1 2 )Δt
ϕ n+1 =ϕ n +f (tn,ϕ
n )+ f (tn+1,ϕn+1)
2Δt
Zeitintegration
Stefan Hickel - Angewandte Strömungssimulation 46
• Alle Verfahren geben gute Lösungen bei kleinen Δt
• Ordnung des Verfahrens besagt, wie schnell der Abbruchfehler gegen null geht, wenn Δt klein genug ist
• „klein genug“ ist vom Verfahren und dem Problem abhängig
• Einen Anhaltspunkt liefert auch hier die Courant-Friedrichs-Lewy Zahl; für CFL < 0.5 sind Zeitintegrationsfehlter i.A. klein
• Implizite Verfahren ermöglichen jedoch wesentlich größere Zeitschrittweiten Δt.
CFL = Δt UΔx
Eigenschaften expliziter Verfahren
Stefan Hickel - Angewandte Strömungssimulation 47
• Effizient, keine weitere Iteration erforderlich • Einfach zu implementieren • Geringer Speicherbedarf
• Instabil bei größeren Zeitschritten • Zeitschrittbegrenzung abschätzten mittels Stabilitätsanalyse
Courant-Friedrichs-Lewy (CFL) , siehe auch GNSM Vorlesung
• Zeitschritt muss der Strömung und der räumlichen Gitterweite angepasst werden!
Eigenschaften impliziter Verfahren
Stefan Hickel - Angewandte Strömungssimulation 48
• Implizite Methoden: Werte von φ zur Zeit t>t0 werden benötigt • Integrale können nur iterativ gelöst werden • Schwieriger zu programmieren • Größerer Speicheraufwand
• i.d.R. unbedingt stabil bei größeren Zeitschrittweiten als explizite Verfahren.
• Der höhere numerische Aufwand pro Zeitschritt (Iterationen) wird oftmals durch die Möglichkeit größerer Zeitschritte ausgeglichen. Größere Zeitschritte führen aber auch zu hoher numerischer Diffusion.
• Vorsicht bei LES.
Prädiktor-Korrektor-Methoden
Stefan Hickel - Angewandte Strömungssimulation 49
• Kombination aus expliziten und impliziten Verfahren • Die Lösung zum neuen Zeitpunkt wird mit einem Eulerzeitschritt
vorhergesagt (Prädiktor):
• Lösung wird korrigiert durch Anwendung des Prädiktors in der impliziten Trapezregel (Korrektor):
• Limit: Fehlerordnung 2
• Vielzahl an verschiedenen Varianten
ϕ n+1* =ϕ n + f (tn,ϕ
n )Δt
ϕ n+1 =ϕ n +12f (tn,ϕ
n )+ f (tn+1,ϕ n+1* )!" #$ Δt
Verfahren höherer Ordnung
Stefan Hickel - Angewandte Strömungssimulation 50
• Zum Erzielen einer höheren Genauigkeitsordnung sind weitere Stützstellen erforderlich
• Diese können entweder Punkte sein, an denen die Lösung bereits zu früheren Zeitpunkten (tn-1, tn-2, tn-3 ,...) bestimmt wurde
• Oder zusätzliche Punkte zwischen tn und tn+1 , welche nur für numerische Zwecke verwendet werden
• Beliebige Genauigkeitsordnung ist erzielbar
Mehrschritt-Verfahren
Runge-Kutta-Verfahren
Mehrschrittverfahren
Stefan Hickel - Angewandte Strömungssimulation 51
Interpolationspolynom durch f(t,φ) zu verschiedenen Zeitpunkten
Adams-Bashforth-Verfahren (explizit)
Adams-Moulton-Verfahren (implizit)
Oder z.B. AB als Prädiktor und AM als Korrektor Problem: Andere Methoden sind erforderlich, um die Simulation zu starten
ϕ n+1 =ϕ n +Δt12
23 f (tn,ϕn )−16 f (tn−1,ϕ
n−1)+ 5 f (tn−2,ϕn−2 )#$ %&
(3. Ordnung)
ϕ n+1 =ϕ n +Δt12
5 f (tn+1,ϕn+1)+8 f (tn,ϕ
n )− f (tn−1,ϕn−1)#$ %&
(3. Ordnung)
Runge-Kutta-Verfahren
Stefan Hickel - Angewandte Strömungssimulation 52
Hilfspunkte zwischen tn und tn+1 werden verwendet Runge-Kutta-Verfahren 2. Ordnung:
Prädiktorschritt für n+1/2 (Euler explizit) Korrektorschritt (Mittelpunktsregel) • Leicht zu implementieren und selbststartend • RK mit 3 bis 5 Schritten am beliebtesten • Große Ähnlichkeit zu den Prädiktor-Korrektor-Methoden
ϕ * =ϕ n +Δt2f (tn,ϕ
n )
ϕ n+1 =ϕ n +Δt ⋅ f (tn+12
,ϕ *)
Lokaler vs. physikalischer Zeitschritt
Stefan Hickel - Angewandte Strömungssimulation 53
• Die maximal erlaubte Zeitschrittweite Δt ist physikalisch und numerisch limitiert.
• Die Lösung darf sich physikalisch nicht schneller ausbreiten als dies durch den Einflussbereich der numerischen Diskretisierung abgebildet werden kann, sonst Instabilität bei expliziten Zeitschrittverfahren: - Einflussbereich: h (≈ Gitterweite) - Signalgeschwindigkeit: S = |u| + c (kompressibel)
S = |u| (inkompressibel) - Max. Zeitschritt: Δt ≈ min( h / S )
• Der erlaubte (physikalische) Zeitschritt wird oft durch die kleinste Zelle
bestimmt und kann daher sehr klein sein.
• Konvergenz zu stationärer Lösung kann durch lokalen Zeitschritt enorm beschleunigt werden: Δt = Δt(x)
Lokaler vs. physikalischer Zeitschritt
Stefan Hickel - Angewandte Strömungssimulation 54
Lokaler Zeitschritt: Transiente ist unphysikalisch! Unzureichend konvergiertes Ergebnis kann komplett falsch sein!
Physikalischer Zeitschritt: Transiente ist physikalisch. Auskonvergiertes stationäres Ergebnis wird sehr langsam (oder nie) erreicht.
Lokaler vs. physikalischer Zeitschritt
Stefan Hickel - Angewandte Strömungssimulation 55
• Auskonvergiertes stationäres Ergebnis ist identisch, wird aber mit lokalem Zeitschritt viel schneller erreicht.
• Transiente ist bei lokalem Zeitschritt unphysikalisch, daher Vorsicht: Unzureichend konvergiertes Ergebnis ist ebenfalls unphysikalisch!
• Lokaler Zeitschritt wird oft als Beschleuniger für sogenannte Duale Zeitschrittverfahren eingesetzt.
Lokaler Zeitschritt Globaler (physikalischer) Zeitschritt
Duales Zeitschrittverfahren
Stefan Hickel - Angewandte Strömungssimulation 56
• Zeitableitung wird aufgeteilt in - einen impliziten Teil mit (großer) physikalischer
Zeitschrittweite Δt - einen Pseudozeit-Anteil mit lokaler Zeitschrittweite Δτ
• Diskretisierung mit Euler implizit für die physikalische Zeit und Euler explizit für die Pseudozeit:
Und Iteration bis Konvergenz zu stationärer Lösung in τ. Für die gesuchte Lösung gilt
∂ϕ∂t
= f (ϕ(t))⇔ limτ→∞
∂ϕ∂τ
= f (ϕ(τ ))− ∂ϕ∂t
%
&'
(
)*
limτ→∞
ϕ(τ +Δτ )−ϕ(τ )Δτ
= −ϕ(t +Δt)−ϕ(t)
Δt+ f (ϕ(τ ))
%
&'
(
)*
ϕ(t +Δt)→ϕ(τ +Δτ ) ≈ϕ(τ )f (ϕ(τ )) → f (ϕ(t +Δt))
Zeitdiskretisierung in ANSYS-CFX
Stefan Hickel - Angewandte Strömungssimulation 57
• Duales Zeitschrittverfahren
• 1st order Euler (implizit) • 2nd order Euler (implizit)
• Lokale oder physikalische (globaler) Zeitschrittweite
• Anwender kann (muss nicht) CFL Zahl vorgeben
è Es gibt (fast) immer ein Simulationsergebnis!
Beispiel
Stefan Hickel - Angewandte Strömungssimulation 59
1-D Diffusions-Advektions-Gleichung
• Euler implizit
• Zentrale-Differenzen-Approximation + äquidistantes Gitter
• Einsetzen liefert die diskretisierte Gleichung:
∂ϕ∂t
= −U ∂ϕ∂x
+Γ∂2ϕ∂x2
ϕ n+1 =ϕ n + f (tn+1,ϕn+1)Δt
∂ϕ∂x"
#$
%
&'i
=ϕi+1 −ϕi−1
2Δx∂2ϕ∂x2"
#$
%
&'i
=ϕi+1 − 2ϕi +ϕi−1
Δx2
ϕin+1 =ϕi
n +Δt −U ϕi+1n+1 −ϕi−1
n+1
2Δx+Γ
ϕi+1n+1 − 2ϕi
n+1 +ϕi−1n+1
Δx2$
%&
'
()
Beispiel
Stefan Hickel - Angewandte Strömungssimulation 60
• Diskretisierte Gleichung
• Umsortieren nach Index i
• Lösung ist implizit gegeben durch
ϕi−1n+1 −
ΔtU2Δx
−ΔtΓΔx2
$
%&
'
()+ϕi
n+1 1+ 2ΔtΓΔx2
#
$%
&
'(+ϕi+1
n+1 ΔtU2Δx
−ΔtΓΔx2
$
%&
'
()=ϕi
n
ϕ n+1
ϕin+1 =ϕi
n +Δt −U ϕi+1n+1 −ϕi−1
n+1
2Δx+Γ
ϕi+1n+1 − 2ϕi
n+1 +ϕi−1n+1
Δx2$
%&
'
()
Beispiel
Stefan Hickel - Angewandte Strömungssimulation 61
• Aufsummieren der konvektiven und diffusiven Flüsse liefert ein algebraisches Gleichungssystem:
- Für jede Transportgleichung folgt in jedem Gitterpunkt eine algebraische Gleichung
- Index P – Punkt an dem die DGL approximiert wurde
- Index i – läuft über alle Knoten die in die Formulierung eingebunden sind
- Koeffizienten Ai – abhängig von der Geometrie und den Fluideigenschaften
- Koeffizient QP – enthält alle bekannten Terme
Matrix Struktur für 5 x 5 Netz
APϕPn+1 + Aiϕi
n+1
i∑ =QP