next up previous contents
Next: 7.2.2 Abbruchbedingungen, Lösungsmethoden Up: 7.2 Konvergenzverhalten des Gleichungslösers Previous: 7.2 Konvergenzverhalten des Gleichungslösers

7.2.1 Physikalisches versus numerisches Lösungsverhalten  

Löst man die gekoppelten nichtlinearen Halbleitergleichungen, so geschieht dies, indem man das Gleichungssystem im momentanen Punkt, der durch die aktuellen Werte der Lösungsvariablen bestimmt ist, als linear annimmt. Dieses linear angenommene Gleichungssystem wird mit Hilfe eines NEWTON-Verfahrens gelöst [9,61]. Man erhält auf diese Weise eine errechnete Änderung der Lösungsvariablen (Korrekturvektor). Wegen der Nichtlinearität des Gesamtproblems ist es bei entsprechender Entfernung von der Lösung nicht möglich, die neuen Lösungsvariablen durch einfache Addition des alten Variablenvektors mit dem kompletten Korrekturvektor zu erhalten. Um ein Überschießen aufgrund der Linearisierung zu vermeiden ist es daher notwendig, den errechneten Korrekturvektor zu dämpfen. Der neue Lösungsvektor errechnet sich dann durch Addition des alten Lösungsvektors mit dem gedämpften Korrekturvektor. Die so errechneten neuen Lösungsvariablen bestimmen den neuen Punkt des Gleichungssystems, der in der Folge wieder als linear angenommen wird und abermals mit Hilfe des NEWTON-Verfahrens gelöst wird. Das einmalige Lösen mit dem NEWTON-Verfahren entspricht dabei einer Iteration. Das nichtlineare Gleichungssystem wird somit durch eine endliche Anzahl von Iterationen gelöst, wobei bei jeder Iteration ein NEWTON-Verfahren angewendet wird.

Ein wesentlicher Punkt bei diesem Lösungsverfahren stellt dabei die Berechnung der Dämpfung des Korrekturvektors dar. Wird der Korrekturvektor nicht gedämpft, so terminiert das Lösungsverfahren im allgemeinen Fall nicht, da das Gleichungssystem zu großen Schwankungen der Lösungsverfahren unterliegt. Wird zu stark gedämpft, dann sind zu viele Iterationen notwendig, um die erforderte Abbruchbedingung zu erreichen, was letztlich zu einer langen Simulationszeit führt,

Die Dämpfung des Korrekturvektors ist umso größer, je weiter man von der Gesamtlösung entfernt ist. Eine Abschätzung für die Entfernung zur Gesamtlösung kann mit Hilfe der ,,Größe`` des Korrekturvektors erfolgen. Diese ,,Größe`` wird dabei durch bestimmte Vektornormen berechnet. Da der Korrekturvektor aus den Teilkorrekturvektoren der Lösungsvariablen besteht, können Normen dieser Teilkorrekturvektoren berechnet werden. Man erhält auf diese Weise ein Bewertungskriterium für die Lösungsvariablen des Potentials, der Trägerdichten, etc.. In MINIMOS-NT hat sich im Drift-Diffusionsmodell und im hydrodynamischen Modell die Dämpfung des Korrekturvektors nach dem Potential am Besten bewährt. Als Bewertungsgröße wird dabei die Maximumnorm (Unendlichnorm) [72] des Potentialkorrekturvektors herangezogen. Durch diese Vorgangsweise wird der Korrekturvektor umso stärker gedämpft, je größer die relative Änderung des Potentials ist. Ist man nahe an der Gesamtlösung, so werden die berechneten Normen der Korrekturen entsprechend kleiner. Für das Beenden einer Simulation wird daher gefordert, daß die Maximumnorm des Potentialkorrekturvektors unter einer geforderten Abbruchbedingung liegen soll.

Das Lösungsverhalten des Gleichungssystems hängt vom betrachteten Arbeitspunkt des Bauteils ab. Dies soll anhand einer Diodensimulation in Fluß- und Sperrichtung erläutert werden.

Wird eine Diode in Sperrichtung betrieben, so werden auch die Änderungen der Trägerkonzentrationen entsprechend der Potentialdämpfung berechnet. Während die Potentialdämpfung immer geringer wird, da sich das Potential immer mehr der Lösung nähert, bleiben die Trägerdichten noch länger unbestimmt. Dies geht so weit, daß die Potentialänderung schon lange unter der geforderten Abbruchbedingung liegt, während sich die Ladungsträger noch weiter verändern. Der Grund dafür liegt in der Ausbildung der Raumladungszone, die bewirkt, daß die relative Norm der Korrekturvektoren der Trägerkonzentrationen weiterhin sehr hoch bleibt. Würde man die Trägerkonzentrationen in der Raumladungszone auf die Hälfte reduzieren, so ergäbe dies eine hohe Korrekturnorm der Trägerdichten, während das Potential davon kaum beeinflußt wird. Der Grund dafür ist die Ladung der Raumladungszone, die zum Großteil durch die aktiven Dopanden hervorgerufen wird, während der Anteil der freien Ladungsträger die Gesamtladung der Raumladungszone kaum beeinflußt.

In Flußrichtung ist die Situation gänzlich verschieden. Die im Bauteil großen Ladungsträgerkonzentrationen bewirken hohe Stromdichten, die das Potential stark beeinflussen. Aus diesem Grund nimmt die Norm des Potentialkorrekturvektors gleichmäßig mit den Korrekturnormen der Trägerkonzentrationen ab. Damit das Iterationsverfahren in Flußrichtung bis zur geforderten Abbruchkorrekturnorm konvergiert, muß sogar die Dämpfung des Potentials verstärkt werden (z.B. durch Erhöhung des Dämpfungsfaktors $\delta$ von Gleichung (9.16) in [18]).

Die Situation im hydrodynamischen Modell ist sehr ähnlich jener im Drift-Diffusionsmodell. Bei der Diode in Sperrichtung kann sich der Energiefluß in den Raumladungszonen nur dann einstellen, wenn sich das Potential nicht mehr ändert. Dies bedingt, daß sich die Trägertemperaturen zusammen mit den Trägerkonzentrationen einstellen.

Werden Bauteile mit hohen Dotierungsbereichen simuliert, in denen sehr kleine Stromdichten fließen, dann kann es vorkommen, daß die errechnete Stromverteilung Zacken aufweist. Diese Zacken sind physikalisch unrealistisch und resultieren aus der diskretisierten Beschreibung des Problems. Eine solche Situation liegt im Fall einer hochdotierten, in Sperrichtung betriebenen Diode vor, wo mögliche Generationsmodelle in der Raumladungszone abgeschaltet sind. Der Grund liegt in der starken Schwankung der Stromberechnung, die von der hohen lokalen Trägerkonzentration der Majoritätsträger abhängt. Die hohe Trägerkonzentration ist die Folge der hohen Dotierung, die von Majoritätsträgern kompensiert werden muß. Die Ursache der unstetigen Stromdichte liegt darin, daß nicht nach den Stromdichten gelöst wird, sondern nach den Ladungsträgerkonzentrationen. Verändert man die lokale Konzentration der Majoritätsladungsträger in der letzten gültigen Dezimalstelle, so kann sich die neu berechnete Stromdichte stark von der alten unterscheiden. Der Effekt ist rein numerisch bedingt und resultiert aus der endlichen Zahlendarstellung.

Diese Tatsache ist besonders bei der Simulation von Durchbrüchen mit dem Drift-Diffusionsmodell entscheidend, wo die Generationsrate proportional zur Teilchenstromdichte ist. Werden jedoch mögliche Generationsmodelle in der Raumladungszone mitberücksichtigt, dann ist der generierte Sperrstrom vor dem Durchbruch bereits so groß, daß die Zacken der Stromverteilung nicht mehr auftreten sollten. Der Effekt des Auftretens von Zacken in der Stromverteilung ist besonders bei Gittertemperaturen um die Raumtemperatur oder darunter gegeben. Bei höheren Temperaturen ist durch die Zunahme der intrinsischen Konzentration der Diffusionsstrom bereits so groß, daß der Effekt in den Hintergrund tritt. Dies ist dadurch begründet, da bei steigender Kristallgittertemperatur die Stromdichte im Bauteil zunimmt, während die Majoritätsträgerkonzentration nahezu gleich bleibt.

Dieses Verhalten zeigt, daß es möglich ist, daß die Norm des Residuums nicht unter einen gewissen Abbruchwert fällt, obwohl die Korrekturvektoren bereits im Bereich der numerischen Genauigkeit bestimmt sind. Umgekehrt ist es nicht immer notwendig, die Lösungsvariablen bis zu hoher Genauigkeit zu lösen, da oft nur der Kontaktstrom interessiert, der sich bei einer angelegten Spannung einstellt.


next up previous contents
Next: 7.2.2 Abbruchbedingungen, Lösungsmethoden Up: 7.2 Konvergenzverhalten des Gleichungslösers Previous: 7.2 Konvergenzverhalten des Gleichungslösers
Martin Knaipp
1998-10-09