2.5 Numerische L¨ osungsverfahren
Das in Abschnitt 2.1 aufgestellte System von partiellen Differentialgleichungen, die eine kompressible reaktive Str¨omung in drei Dimensionen beschreiben, soll numerisch gel¨ost werden. Zur L¨osung der partiellen Differentialgleichungen wird eine Diskretisierung der r¨aumlichen Ableitungen vorgenommen, wodurch ein Sy-stem von gew¨ohnlichen Differentialgleichungen entsteht. Dieses wird dann mit Hilfe eines geeigneten Verfahrens in der Zeit integriert.
Diskretisierung im Raum
F¨ur die r¨aumliche Diskretisierung von partiellen Differentialgleichungen stehen mehrere Verfahren zur Auswahl. Sie lassen sich in zwei Gruppen einteilen, n¨amlich in finite Differenzenverfahren (FD) und in Verfahren der gewichteten Residuen.
Finite Differenzenverfahren zeichnen sich im wesentlichen dadurch aus, daß sie erheblich einfacher zu realisieren sind und auf strukturierten Gittern eine einfa-che Verwaltung der Gitterpunkte erlauben. Zu den Residuen-Verfahren z¨ahlen die Methode der finiten Volumen (FV) und finiten Elemente (FE). Diese Methoden zeichnen sich durch ihre gute Verwendbarkeit bei unstrukturierten oder geome-trisch komplizierten Gittern aus [Rie98]. Die Verwaltung der Nachbarschaftsbe-ziehungen der Knoten und Elemente erfordert allerdings einen erheblichen logi-stischen Aufwand.
Diskretisierung in der Zeit
Zur zeitlichen Integration der nach der ¨ortlichen Diskretisierung resultierenden gew¨ohnlichen Differentialgleichungen k¨onnen sowohl implizite als auch explizite Verfahren verwendet werden. Explizite Verfahren sind einfacher und ben¨otigen pro Rechenschritt bei gleicher Fehlerordnung deutlich weniger Rechenzeit. Bei impliziten Verfahren muß in jedem Rechenschritt die Jacobi-Matrix des Systems berechnet und invertiert werden.
Die Ber¨ucksichtigung detaillierter Chemie verursacht eine hohe Steifigkeit des Sy-stems, da die Zeitskalen der Reaktionen sehr unterschiedlich sind. Hinzu kommen auch die exponentielle Temperaturabh¨angigkeit der Kinetik durch das Arrhenius-Gesetz. Deshalb wird eine implizite Zeitdiskretisierung gew¨ahlt. Faßt man den L¨osungsvektor in ~s und alle Terme, die keine Zeitableitung enthalten, in F~ zu-sammen, dann gilt
B∂~s
∂t =F~(~s, t) , (2.100)
wobei B eine Diagonalmatrix ist, deren i-tes Diagonalelement den Wert 1 bzw.
0 besitzt, wenn die i-te Gleichung des zu l¨osenden Gleichungssystems eine Diffe-rentialgleichung bzw. eine algebraische Gleichung darstellt. Mit Hilfe des
Euler-26 Kapitel 2: Grundlagen
impliziten Ansatzes wird dieses Gleichungssystem linearisiert, und man erh¨alt B~sn+1−~sn
∆t = F~(~s, t)n+1 ,
= F~(~s, t)n+F~s(~sn+1−~sn) , (2.101)
~sn+1−~sn
| {z }
∆~s
= [B−∆tF~s]
| {z } M
−1∆t ~F(~s, t)n
| {z }
~b
,
wobei die Jacobi-Matrix F~s = F~s~
n die partiellen Ableitungen der Funktion F~ abgeleitet nach~szum Zeitpunktn enth¨alt. Somit erh¨alt man das folgende lineare Gleichungssystem
∆~s=M−1·~b , (2.102)
welches mit dem impliziten Extrapolationsverfahren LIMEX [DHZ87, DN87]
gel¨ost wird. Die verwendete Methode zur Schrittweiten- und Ordnungskontrol-le ist speziell f¨ur steife Differential-Algebraische-Gleichungssysteme entwickelt worden. Zur Erh¨ohung der Ordnung wird jeder Zeitschritt ∆t in j Teilschrit-te mit ∆tj = ∆t/j zerlegt. F¨ur jede Unterteilung wird dann in j Stufen eine L¨osung durch Integration bestimmt. Die Extrapolation f¨ur eine unendlich fei-ne Unterteilung ∆tj → 0 ergibt dann die gesuchte L¨osung f¨ur den Zeitschritt
∆t. Die Fehlerordnung betr¨agt dabei h¨ochstens j. Innerhalb der Extrapolation werden Fehlerabsch¨atzungen durchgef¨uhrt. Darauf basieren dann die Berechnung einer optimalen neuen Schrittweite f¨ur den n¨achsten Integrationsschritt sowie ein Kriterium f¨ur die Neuberechnung der Jacobi-Matrix. Dieses Verfahren ist zwar zun¨achst sehr aufwendig, da beijmax Verfeinerungsschritten Pjimaxi Teill¨osungen bestimmt werden m¨ussen. Durch die optimale Wahl der Schrittweite und Ord-nung sowie die m¨oglichst h¨aufige Wiederverwendung der Jacobi-Matrix ist es jedoch sehr effektiv.
27
Kapitel 3
Modellerweiterungen f¨ ur ionisierte Gase
3.1 Das Plasma
Um von einem Plasma sprechen zu k¨onnen, muß die Materie zwei Vorausset-zungen erf¨ullen: Es muß zu elektromagnetischen Wechselwirkungen zwischen den geladenen Teilchen kommen, und die Anzahl der positiven und negativen La-dungstr¨ager je Volumeneinheit muß ausgeglichen sein (Quasineutralit¨at). Zur Erzeugung des Plasmazustandes wird i.a. eine Energiequelle ben¨otigt, welche die Ionisation bewirkt und diesen Zustand aufrechterh¨alt. Damit ist in vielen F¨allen auch eine Aufheizung des Gases auf hohe Temperaturen verbunden, doch darf dies nicht als typische Plasmaeigenschaft angesehen werden.
Plasma-Modelle
Man spricht von einem vollst¨andigen thermodynamischen Gleichgewicht (VTG, engl.: complete thermodynamic equilibrium (CTE)), wenn das Gas durch ma-kroskopische Zustandsgr¨oßen, insbesondere durch eine einheitliche Temperatur, beschrieben werden kann.
Im Unterschied zu gew¨ohnlichen Aggregatzust¨anden ist es i.d.R. bei einem Plas-ma nicht m¨oglich, der Materie eine Temperatur zuzuordnen. Zwischen den Ga-sen des Plasmas besteht kein thermisches Gleichgewicht, weil der Energieaus-tausch durch St¨oße aufgrund der großen Massenunterschiede zwischen Elektro-nen und Schwerteilchen sehr ineffizient ist [LL94]. Auch wenn sich ein solches Plasma nicht im VTG befindet, k¨onnen bei gen¨ugend hoher Elektronendich-te lokal, also in jedem beliebigen kleinen Volumenelement, sowohl die Elek-tronen als auch die schweren Teilchen in guter N¨aherung durch die Maxwell-Geschwindigkeitsverteilung mit den Temperaturen Te und T beschrieben wer-den. Die Besetzungszahldichte der angeregten Zust¨ande folgt einer Boltzmann-Verteilung. Die Spezies-Konzentrationen lassen sich bestimmen, wenn man che-misches Gleichgewicht annimmt. Ein solcher Plasmazustand heißt lokales
ther-28 Kapitel 3: Modellerweiterungen f¨ur ionisierte Gase
modynamisches Gleichgewicht (LTG, engl.: local thermodynamic equilibrium (LTE)) [SW92, BFP94, LM98].
Statt einer einheitliches Gastemperatur k¨onnen Atome bzw. Molek¨ule und Ionen auch unterschiedliche Temperaturen T und Tion gem¨aß der Maxwell-Verteilung ihrer Geschwindigkeiten haben, so daß man zum Drei-Temperatur-Modell [Spa90]
gelangt.
Bei den typischen Bedingungen der Plasma¨atzverfahren zur Halbleiterherstellung (geringer Druck und kleiner Ionisationsgrad) verhindert die kontinuierliche Ener-giezuf¨uhrung durch die Plasmaquelle sogar den Energieaustausch zwischen den Elektronen, so daß sich selbst diese nicht im Gleichgewicht befinden und deren Verteilungsfunktion (engl.: electron energy distribution function (EEDF)) von der einer Maxwell-Verteilung abweicht [ARW98]. Trotzdem kann man den Elek-tronen aufgrund ihrer mittleren Geschwindigkeit eine orts- und zeitabh¨angige Temperatur zuordnen.
Die Debyesche Abschirml¨ange
Typische Plasmaeigenschaften treten dann auf, wenn die Reichweite der elektro-statischen Wechselwirkung zwischen den Teilchen gr¨oßer ist als die mittlere freie Wegl¨ange. Die Debyesche Abschirml¨ange [SH53]
λ˜D= ist gegeben durch die Entfernung von einer einzelnen Punktladung, innerhalb der das Potential infolge der Raumladungseffekte benachbarter Ladungen auf ein e-tel (≈0,37) sinkt, und kann als Reichweite der elektrostatischen Wechselwirkung interpretiert werden. 0 ist die Dielektrizit¨atskonstante, e die Elementarladung und ne die Teilchendichte der Elektronen. Die effektive Ionenladungszahl
Zeff = 1
ist als mit den Teilchendichten ni gewichteter Mittelwert ¨uber das Quadrat der Ladungszahlen Zi der Schwerteilchen definiert. Die Summation geschieht zwar
¨
uber alle Gasphasenspezies Ng, die Ladungszahlen der Neutralteilchen sind je-doch Null. Die Quasineutralit¨atsbedingung des Plasmas ist erf¨ullt, wenn in jedem Volumenelement
Ng
X
i=1
niZi = 0 (3.3)
gilt. In vielen Anwendungen kommen nur einfach geladene positive Ionen vor, so daß aufgrund der Quasineutralit¨atsbedingung die effektive Ionenladungszahl Eins