DE69033893T2

DE69033893T2 - Verfahren zur Analyse des Betriebes eines Halbleiterbausteins, Verfahren zur Analyse eines spezifischen physischen Phänomenes und Gerät zur Durchführung dieser Verfahren

Info

Publication number: DE69033893T2
Application number: DE69033893T
Authority: DE
Inventors: Mamoru Kurata; Shin Nakamura
Original assignee: Toshiba Corp
Current assignee: Toshiba Corp
Priority date: 1989-09-29
Filing date: 1990-09-28
Publication date: 2002-09-12
Anticipated expiration: 2010-09-29
Also published as: EP0421684B1; DE69033893D1; EP0421684A3; US6041424A; EP0421684A2

Description

Die Erfindung bezieht sich auf das Modellieren von Halbleitervorrichtungen und spezifischen physikalischen Phänomenen.
Die auf das Modellieren von Halbleitervorrichtungen angewendeten Grundgleichungen werden allgemein in der folgenden Form dargestellt:
∂p/∂τ = (-1/q)(divJp) + G - U ... (1)
∂n/∂τ = (1/q)(divJn) + G - U ... (2)
Jp = - qDp(gradp) - qupp(gradψ) ... (3)
Jn = qDn(gradn) - qunn(gradψ) ... (4)
div(gradψ) = (-q/E)(Nd - Na + p - n) ... (5)
Diese fünf Gleichungen sind beispielsweise aus Muller & Kamins, "Device electronics for integrated circuits", Wiley 1977 bekannt; sie werden üblicherweise auf das herkömmliche Modellierungsverfahren und das Modellierungsverfahren gemäß dieser Erfindung angewendet. Die Gleichungen (1) und (2) sind eine Kontinuitätsgleichung für Löcher bzw. eine Kontinuitätsgleichung für Elektronen. Bei diesen Gleichungen ist τ die physikalische Zeit, und G und U bezeichnen die Erzeugung und Rekombination von Ladungsträgern (im folgenden häufig als Träger bezeichnet). Die Gleichungen (3) und (4) beschreiben explizit eine Lochstromdichte bzw. eine Elektronenstromdichte. Wie es aus Gleichung (3) ersichtlich ist, ist die Lochstromdichte Jp eine algebraische Summe einer Diffusionsstromkomponente (das erste Glied), die proportional dem Trägerdichtengradienten gradp ist, und einer Driftkomponente (das zweite Glied), die proportional dem Potentialgradienten gradψ und der Trägerdichte ist. Wie es aus Gleichung (4) ersichtlich ist, ist auch die Elektronenstromdichte Jn eine algebraische Summe einer Diffusionsstromkomponente (das erste Glied), die proportional dem Trägerdichtengradienten gradn ist, und einer Driftkomponente (das zweite Glied), die proportional dem Potentialgradienten grad ψ und der Trägerdichte ist. Gleichung (5) ist eine Poisson-Gleichung, die die Beziehung zwischen dem Potential ψ und der Raumladungsdichte darstellt.
Das herkömmliche Verfahren des Modellierens von Halbleitervorrichtungen wird zuerst erläutert, siehe ebenfalls J. P. Darling "Parallel algorithms for the solution of nonlinear Poisson equation of semiconductor device theory and its implementation on the MPP", IEEE, Proc. of the 2nd Symposium on the Frontiers of Massively Parallel Computation, 10. bis 12. Oktober 1988, Seiten 289-294. Bei dem ersten Schritt eines Erhaltens der nummerischen Lösungen von den Gleichungen (1) bis (5) werden diese Gleichungen vom kontinuierlichen System (d. h. Differentialgleichungen) in ein diskretes System (d. h. Differenzgleichungen) umgeschrieben. Um die Gleichungen (1) bis (5) derart umzuschreiben, ist es notwendig, den Vorrichtungsraum, der zu analysieren ist, in ein rechtwinkliges Gitter- oder Maschenpunkt-Array aufzuteilen. Die Maschenaufteilung für ein zweidimensionales Raummodell ist schematisch in Fig. 1 dargestellt. Es sei der Maschenpunktabstand hx in der x-Richtung und der Maschenpunktabstand hy in der y-Richtung festgelegt, wie es in Fig. 1 dargestellt ist; es sei ferner der Mittelpunktabstand hx' in der x-Richtung und der Mittelpunktabstand hy' in der y-Richtung wie folgt festgelegt:
hx'(M) = (1/2) [hx(M' - 1) + hx(M')] ... (6-1)
hy'(N) = (1/2) [hy(N' - 1) + hy(N')] ... (6-2)
Aus den Gleichungen (3), (6-1) und (6-2) kann die Lochstromdivergenz div Jp durch die folgende Näherung in dem diskreten System dargestellt werden:
divJp = ∂Jpx/∂x + ∂Jpy/∂y = [Jpx(M') - Jpx(M' - 1)]/hx'(M) + [Jpy(N') - Jpy(N' - 1)]/hy'(N) ... (7)
Hinsichtlich des Elektronenstroms kann die Divergenz div Jn auf die gleiche Art und Weise umgeschrieben werden. Daher werden die Gleichungen (1) und (2) umgeschrieben in:
(1/q) {[Jpx(M') - Jpx(M' - 1)]/hx'(M) + [Jpy(N') - Jpy(N' - 1)]/hy'(N)} = - U(M,N) ... (8-1)
(1/q) {[Jnx(M') - Jnx(M' - 1)]/hx'(M) + [Jny(N') - Jny(N' - 1)]/hy'(N))} = U(M,N) ... (8-2)
Bei den Gleichungen (8-1) und (8-2) wird das Trägererzeugerglied G zwecks Vereinfachung weggelassen, da dieses Glied beim Analysieren von Halbleitervorrichtungen unter gewöhnlichen Betriebsbedingungen nicht erforderlich ist. Ferner wird zwecks Vereinfachung die physikalischen Zeitdifferentialglieder ∂p/∂τ und ∂n/∂τ unter der Annahme gleich Null gesetzt, dass die Berechnungen durchgeführt werden, um eine stationäre Gleichstromlösung zu finden. Es sei jedoch bemerkt, dass Berechnungen im Wesentlichen auf die gleiche Art und Weise durchgeführt werden können, um eine nicht-stationäre Lösung zu erhalten, bei denen die physikalischen Zeitdifferentialglieder nicht Null sind.
Auf ähnliche Weise kann die Poisson-Gleichung (5) ausgedrückt werden wie folgt:
[1/hx'(M)] {[ψ(M + 1,N) - ψ(M,N)]/hx(M') - [ψ(M,N) - ψ(M - 1,N)]/hx(M' - 1)} + [1/hy'(N)] {(ψ(M,N + 1) - ψ(M,N)]/hy(N') - [ψ(M,N) - ψ(M,N - 1)]/hy(N' - 1)} = - (q/&epsi;) [Γ(M,N) + p(M,N) - n(M,N)] ... (9)
wobei zwecks Einfachheit die Dotierungsfunktion, d. h. die Differenz Γ zwischen der Donatorkonzentration und der Akzeptorkonzentration als Γ = Nd - Na gesetzt wird.
Hier wird das Trägerrekombinationsglied U(M,N) in der Shockely-Read-Hall-Form für die nachfolgende Entwicklung ausgedrückt. Das heißt:
U(M,N) = [p(M,N)n(M,N) - ni²(M,N)]/ {τn [p(M,N) + ni(M,N)] + τp [n(M,N) + ni(M,N)]} ... (10)
Als nächstes werden die Gleichungen (3) und (4) für Ströme durch Anwenden der Scharfetter-Gummel-Näherung diskretisiert, um, wie es in der Technik bekannt ist, die Stabilität bei der nummerischen Analyse zu gewährleisten. Als ein Ergebnis dessen werden die x-Komponente und die y- Komponente der Lochstromdichte und diejenigen der Elektronenstromdichte wie folgt gegeben:
Jpx(M') = [q/hx(M')] [λpx1(M')P(M N) + λpx2(M')p(M+1,N)] ... (11-1)
Jpy(N') = [q/hy(N')] [λpy1(N')P(M N) + λpy2(N')p(M,N+1)] ... (11-2)
Jnx(M') = [q/hx(M')] [λnx1(M')n(M N) + λnx2(M')n(M + 1,N)] ... (11-3)
Jny(N') = [q/hy(N')] [λny1(N')n(M N) + λny2(N')n(M,N + 1)] ... (11-4)
Die λ-Elemente in den Gleichungen (11-1) bis (11-4) werden durch die folgenden Gleichungen definiert:
λpx1(M') = [up(M')/θ(M')] [β(M')/(1 - e-β(M'))}]
λpx2(M') = [up(M')/θ(M')] [β(M')/(1 - eβ(M'))}]
λpy1(N') = [up(N')/θ(N')] [β(N')/(1 - e-β(N'))}]
λpy2(N') = [up(N')/θ(N')] [β(N')/(1 - eβ(N'))}]
λnx1(M') = [un(M')/up(M')] λpx2(M')
λnx2(M') = [un(M')/up(M')] λpx1(M')
lny1(N') = [un(N')/up(N')] λpy2(N')
Xny2(N') = [un(N')/up(N')] λpy1(N')
... (12)
Somit wurden die Grundgleichungen (1) bis (5) zum Modellieren von Halbleitervorrichtungen in diskrete Grundgleichungen umgeschrieben. Die Gleichungen, die wir zu lösen haben, sind: Gleichung (9), Gleichung (8-1) mit den Gleichungen (11-1) bis (11-4), Hilfsbeziehungen (12) und die darin eingesetzte Rekombinationsgleichung (10), sowie Gleichung (8-2) mit den Gleichungen (11-1) bis (11-4), Hilfsbeziehung (12) und die darin eingesetzte Rekombinationsgleichung (10). Diese drei Gleichungen enthalten drei unbekannte Größen p, n und ψ. Es ist jedoch schwierig, die Gleichungen (8-1) und (8-2) zu lösen, es sei denn, dass diese Gleichungen transformiert werden, da diese Gleichungen jeweils nichtlineare Beziehungen für die drei unbekannten Größen enthalten, Somit werden die Gleichungen (8-1) und (8-2) in Taylor-Reihen für die unbekannten Größen entwickelt, wobei quadratische Glieder und Glieder höherer Ordnung vernachlässigt werden, um dadurch lineare Gleichungen für die Änderungen δp, δn und δψ der unbekannten Größen p, n und ψ herzuleiten.
Genauer gesagt nehmen wir zuerst an, dass die unbekannten Größen jeweils die Summe eines Testwerts (mit hochgestelltem Zeichen 0 gekennzeichnet) und einer Änderung sind. Das heißt:
p(M,N) = p&sup0;(M,N) + δp(M,N)
n(M,N) = n&sup0;(M,N) + δn(M,N)
ψ(M,N) = ψ&sup0;(M,N) + δψ(M,N) ... (13)
Dann werden die Stromgleichungen (11-1) bis (11-4) in Taylor-Reihen entwickelt, wobei die quadratischen Glieder und Glieder höherer Ordnung vernachlässigt werden, womit die folgenden Ergebnisse erhalten werden:
Jpx(M') Jpx&sup0;(M') + [∂Jpx&sup0;(M')/∂p(M,N)] δp(M,N) + [∂Jpx&sup0;(M')/∂p(M+1,N)] δp(M+1,N) + [∂Jpx&sup0;(M')/∂ψ(M,N)] δψ(M,N) + [∂Jpx&sup0;(M')/δψ(M+1,N)] δψ(M+1,N) ... (14-1)
Jpy(N') Jpy&sup0;(N') + [∂Jpy&sup0;(N')/∂p(M,N)] δp(M,N) + [∂Jpy&sup0;(N')/∂p(M,N+1)] δp(M,N+1) + [∂Jpy&sup0;(N')/∂ψ(M,N)] δψ(M,N) + [∂Jpy&sup0;(N')/∂ψ((M,N+1)] δψ(M,N+1) ... (14-2)
Jnx(M') Jnx&sup0;(M') + [∂Jnx&sup0;(M')/∂n(M,N)] δn(M,N) + [∂Jnx&sup0;(M')/∂n(M+1,N)] δn(M+1,N) + [∂Jnx&sup0;(M')/∂ψ(M,N)] δψ(M,N) + [∂Jnx&sup0;(M')/∂ψ(M+1,N)] δψ(M+1,N) .... (14-3)
Jny(N') Jny&sup0;(N') + [∂Jny&sup0;(N')/∂n(M,N)] δn(M,N) + [∂Jny&sup0;(N')/∂n(M,N+1)] δn(M,N+1) + [∂Jny&sup0;(N')/∂ψ(M,N)) δψ(M,N) + [∂Jny&sup0;(N')/∂ψ(M,N+1)] δψ(M,N+1) ... (14-4)
wobei das hochgestellte Zeichen 0 Größen zugeordnet ist, die p&sup0;, n&sup0; und ψ&sup0; entsprechen.
Zusätzlich zu den Stromgleichungen ist das Rekombinationsglied U(M,N) ebenfalls nichtlinear. Das Rekombinationsglied U(M,N) wird daher in Taylor-Reihen entwickelt. Das heißt:
U(M,N) = U&sup0;(M,N) + Up&sup0;(M,N)δp(M,N) + Un&sup0;(M,N)δn(M,N)
wobei
Up(M,N) = ∂U(M,N)/∂p(M,N) = [n(M,N) - τnU(M,N)]/ {τn [p(M,N) + ni(M,N)]} + τp [n(M,N) + ni(M,N))}
Un(M,N) = ∂U(M,N)/∂n(M,N) = [p(M,N) - τpU(M,N)]/ {τn [p(M,N) + ni(M,N)] + τp [n(M,N) + ni(M,N)]} ... (15)
Die Gleichungen (13), (14-1), (14-2), (14-3), (14-4) und (15) werden in die Gleichungen (8-1), (8-2) und (9) eingesetzt, wobei die folgende Matrix-Vektorgleichung für das zweidimensionale Problem erhalten wird:
AT(M,N)ΘT(M,N - 1) + BT(M,N)ΘT(M,N) + CT(M,N)ΘT(M,N+1) + DT(M,N)ΘT(M-1,N) + ET(M,N)ΘT(M+1,N) = FT(M,N) ... (16)
wobei für die Definition von Matrizen und Vektoren gilt:
Die Matrixelemente und die Vektorelemente können explizit wie folgt dargestellt werden:
a&sub1;&sub1; = - [1/qhy'(N)] [∂Jpy&sup0;(N-1)/∂p(M,N-1)]
a&sub1;&sub2; = 0
a&sub1;&sub3; = - [1/qhy'(N)] [∂Jpy&sup0;(N-1)/∂ψ(M,N-1)]
a&sub2;&sub1; = 0
a&sub2;&sub2; = - [1/qhy'(N)] [∂Jny&sup0;(N-1)/∂n(M,N-1)]
a&sub2;&sub3; = - [1/qhy'(N)] [∂Jny&sup0;(N-1)/∂ψ(M,N-1)]
a&sub3;&sub1; = 0
a&sub3;&sub2; = 0
a&sub3;&sub3; = γ1(M,N) = 1/[hy'(N)hy(N-1)]
b&sub1;&sub1; = [1/hx'(M)] [∂Jpx&sup0;(M)/∂p(M,N) - ∂Jpx&sup0;(M-1)/∂p(M,N)] + [1/hy'(N)] [∂Jpy&sup0;(N)/∂p(M,N) - ∂Jpy&sup0;(N-1)/∂p(M,N)] + ∂U&sup0;(M,N)/∂p(M,N)
b&sub1;&sub2; = ∂U&sup0;(M,N)/∂n(M,N)
b&sub1;&sub3; = [1/qhx'(M)] [∂Jpx&sup0;(M)/∂ψ(M,N) - ∂Jpx&sup0;(M-1)/∂ψ(M,N)] + [1/qhx'(M)] [∂Jpy&sup0;(N)/∂ψ(M,N) - ∂Jpy&sup0;(N-1)/∂ψ(M,N)]
b&sub2;&sub1; = - ∂U&sup0;(M,N)/∂p(M,N)
b&sub2;&sub2; = [1/qhx'(M)] [∂Jnx&sup0;(M)/∂n(M,N) - ∂Jnx&sup0;(M-1)/∂n(M,N)] + [1/qhy'(N)] [∂Jny&sup0;(M)/∂n(M,N) - ∂Jny&sup0;(N-1)/∂n(M,N)] - ∂U&sup0;(M,N)/∂n(M,N)
b&sub2;&sub3; = [1/qhx'(M)] [∂Jny&sup0;(N)/∂ψ(M,N) - ∂Jnx&sup0;(M-1)/∂ψ(M,N)] + [1/qhy'(N)] [∂Jny&sup0;(N)/∂ψ(M,N) - ∂Jny&sup0;(N-1)/∂ψ(M,N)]
b&sub3;&sub1; = q/&epsi;, b&sub3;&sub2; = - q/&epsi;
b&sub3;&sub3; = - [1/hx'(M)] [1/hx(M) + 1/hx(M-1)] - [1/hy'(N) + 1/hy(N-1)]
c&sub1;&sub1; = [1/qhy'(N)] [∂Jpy&sup0;(N)/∂p(M,N+1)]
c&sub1;&sub2; = 0
c&sub1;&sub3; = [1/qhy'(N)] [∂Jpy&sup0;(N)/∂ψ(M,N+1)]
c&sub2;&sub1; = 0
c&sub2;&sub2; = - [1/qhy'(N)] [∂Jny&sup0;(N)/∂n(M,N+1)]
c&sub2;&sub3; = [1/qhy'(N)] [∂Jny&sup0;(N)/∂ψ(M,N+1)]
c&sub3;&sub1; = 0, c&sub3;&sub2; = 0
c&sub3;&sub3; = 1/hy'(N)hy(N)
d&sub1;&sub1; = - [1/qhx'(M)] [∂Jpx&sup0;(M-1)/∂p(M-1,N)]
d&sub1;&sub2; = 0
d&sub1;&sub3; = - [1/qhX'(M)] [∂Jpx&sup0;(M-1)/∂ψ(M-1,N)]
d&sub2;&sub1; = 0
d&sub2;&sub2; = - [1/qhx'(M)] [∂Jnx&sup0;(M-1)/∂n(M-1,N)]
d&sub2;&sub3; = - [1/qhX'(M)] [∂Jnx&sup0;(M-1)/∂ψ(M-1,N)]
d&sub3;&sub1; = 0 , d&sub3;&sub2; = 0
d&sub3;&sub3; = 1/hx'(M)hx(M-1)
e&sub1;&sub1; = [1/qhx'(M)] [∂Jpx&sup0;(M)/∂p(M+1,N)]
e&sub1;&sub2; = 0
e&sub1;&sub3; = [1/qhX'(M)] [∂Jpx&sup0;(M)/∂ψ(M+1,N)]
e&sub2;&sub1; = 0
e&sub2;&sub2; = [1/qhx'(M)] [∂JnX&sup0;(M)/∂n(M+1,N)]
e&sub2;&sub3; = [1/qhx'(M)] [∂Jnx&sup0;(M)/∂ψ(M+1,N)]
e&sub3;&sub1; = 0, e&sub3;&sub2; = 0
e&sub3;&sub3; = 1/hx'(M)hx(M)
f&sub1; = - U&sup0;(M,N) + Up&sup0;(M,N)p&sup0;(M,N) + Un&sup0;(M,N)n&sup0;(M,N) - [1/qhx'(M)] [∂Jp&sup0;(M-1)/∂ψ(M-1,N)] · [y&sup0;(M-1,N) - ψ&sup0;(M,N)] + [1/qhx'(M)] [∂Jp&sup0;(M)/∂ψ(M,N)] · [ψ&sup0;(M,N) - ψ&sup0;(M+1,N)] - [1/qhy'(N)] [∂Jp&sup0;(N-1)/∂ψ(M,N-1)] · [ψ&sup0;(M,N-1) - ψ&sup0;(M,N)] + [1/ghy'(N)] [∂Jp&sup0;(N)/∂ψ(M,N)] · [y&sup0;(M,N) - ψ&sup0;(M,N+1)]
f&sub2; = U&sup0;(M,N) - Up&sup0;(M,N)p&sup0;(M,N) - Un&sup0;(M,N)n&sup0;(M,N) - [1/qhx'(M)] [∂Jn&sup0;(M-1)/∂ψ(M-1,N)] · [ψ&sup0;(M-1,N) - ψ&sup0;(M,N)] + [1/qhx'(M)] [∂Jn&sup0;(M)/∂ψ(M,N)] · [ψ&sup0;(M,N) - ψ&sup0;(M+1,N)] - [1/qhy'(N)] [∂Jn&sup0;(N-1)/∂ψ(M,N-1)] · [ψ&sup0;(M,N-1) - y&sup0;(M,N)] + [1/qhy'(N)] [∂Jn&sup0;(N)/∂ψ(M,N)] · [ψ&sup0;(M,N) - ψ&sup0;(M,N+1)]
f&sub3; = - (q/&epsi;)Γ(M,N) ... (18)
Wie es aus dem obigen ersichtlich ist, kann die Matrix- Vektorgleichung (16) durch Linearisierung der diskreten Gleichungen (8-1), (8-2) und (9) hergeleitet werden.
Die Gleichung (16) entspricht einem der vielen in der x- Richtung und der y-Richtung angeordneten Maschenpunkten. Somit müssen, um eine Lösung für die gesamte Vorrichtungsebene zu erhalten, die Gleichungen (16) für die Anzahl von Maschenpunkten gleichzeitig gelöst werden. In dem 1 Fall, in dem die Anzahlen von Maschenpunkten, die die gesamte Vorrichtungsebene bilden, beispielsweise in den Bereichen von 1 ≤ M < K und 1 ≤ N ≤ L sind, werden die gesamten Gleichungen (16) wie in Fig. 2 dargestellt. Zwecks Einfachheit wurden AT (M;N), ΘT (M,N) und FT (M,N) in AMN, ΘMN bzw. FMN umgeschrieben.
In Fig. 2 weist jede Matrix, beispielsweise B&sub1;&sub1;, 3 · 3 Dimensionen auf. Daher umfasst das die gesamte Vorrichtungsebene definierende Matrixsystem (3KL) · (3KL) Dimensionen. In dem Fall eines Standardvorrichtungsmodells ist L = 30 und M = 40. Somit weist das Matrixsystem für das Standardvorrichtungsmodell die Dimensionen von 3.600 · 3.600 = 12.960.000 = 1,296 · 10&sup7; auf.
Verschiedene Verfahren sind verfügbar, die verwendet werden, um die in Fig. 2 schematisch dargestellte Matrix- Vektorgleichung zu lösen. Das am häufigsten verwendete Verfahren ist die Gauss'sche Eliminierung. Das Gauss'sche Eliminierungsverfahren, das ein Matrixsystem der oben spezifizierten Größe beinhaltet, kann mittels eines Hochgeschwindigkeits-Computers mit großem Speicher durchgeführt werden, der nun handelsüblich verfügbar ist. Weitere alternative Verfahren, die verwendet werden können, um die Matrix-Vektorgleichung von Fig. 2 zu lösen, sind Integrationsverfahren, die insbesondere für die numerische Analyse von sehr großen Matrixsystemen wirksam sind.
Ein Matrixlösungsverfahren wird nicht nur bei der Vorrichtungsmodellierung, sondern ebenfalls beim Lösen von physikalischen und technischen Problemen verwendet. Da das Verfahren ein in der Technik bekanntes Berechnungsverfahren ist, sei hier angenommen, dass die Gleichung (16) für den gesamten Raum der Vorrichtung gelöst wurde. Dies bedeutet offensichtlich, dass die Lösungen der linearisierten Gleichung (16), die aus den nichtlinearen Gleichungen (8-1), (8-2) und (9) herrührt, erhalten wurden. Wenn diese Lösungen in die nichtlinearen Gleichungen (8-1), (8-2) und (9) eingesetzt werden, werden in den meisten Fällen beide Seiten der Gleichungen in Folge des nichtlinearen Effekts nicht gleich sein. Es sei daher die Lösung [p(M,N), n(M,N), ψ(M,N), 1 ≤ M ≤ K, 1 ≤ N ≤ L] als ein Satz von Probewerten [p&sup0;(M,N), n&sup0;(M,N), ψ&sup0;(M,N), 1 ≤ M ≤ K, 1 ≤ N ≤ L] eingesetzt und dann die gleichen Berechnungsprozeduren wie oben erwähnt ausgeführt. Wenn dieses Berechnungsverfahren eine notwendige Anzahl von Malen wiederholt wird, dann können wir Lösungen erhalten, die die Gleichungen (8-1), (8-2) und (9) vollständig erfüllen.
Das oben beschriebene Matrixlösungsverfahren kann erreicht werden, um sehr große Aufgabenstellungen mittels des Computers, der nun verfügbar ist, zu lösen, der einen großen Massenspeicher aufweist und Vorgänge mit hoher Geschwindigkeit durchführen kann. Je höher die Leistung des Computers ist, desto größer wird jedoch eine Aufgabenstellung, die Leute mittels des Computers zu lösen wünschen. Im Hinblick darauf ist es immer noch eine der größten Schwierigkeiten, Matrixaufgabenstellungen mit ausreichend hoher Geschwindigkeit zu lösen, um ein erfolgreiches Modellieren von Halbleitervorrichtungen zu erreichen.
Um diese Gleichungen mit hoher Geschwindigkeit zu lösen, wurden Verfahren und ein Gerät vorgeschlagen, die bestimmte Gleichungen ohne Anwendung von Matrixlösungsverfahren lösen können, wodurch die Funktionsweisen von Halbleitervorrichtungen und spezifischen physikalischen Phänomenen analysiert werden können. Die Erfindung bezieht sich auf ein verbessertes Verfahren zum Bestimmen von Empfindlichkeitskoeffizienten oder der Verstärkungsfaktoren, die helfen, den Wirkungsgrad der Analyse der Funktionsweise von Halbleitervorrichtungen und der spezifischen physikalischen Phänomene zu verbessern.
Eine Aufgabe der Erfindung ist es, ein Verfahren und ein Gerät zum Lösen von bestimmten Gleichungen ohne Anwenden von Matrixlösungsschemata bereitzustellen, um dadurch die Funktionsweisen von Halbleitervorrichtungen und spezifische physikalische Phänomene zu analysieren.
Diese Aufgabe wird mit der Erfindung gelöst, wie es in dem beigefügten Verfahren nach Anspruch 1 und dem Gerät nach Anspruch 20 dargelegt ist.
Erfindungsgemäß werden verschiedene Verfahren zum Lösen der Funktionsweisen bzw. des Betriebs von Halbleitervorrichtungen oder spezifische physikalische Phänomene und ebenfalls verschiedene Geräte zum Analysieren der Funktionsweisen von Halbleitervorrichtungen oder spezifische physikalische Phänomene bereitgestellt.
Das erste Verfahren, das ausgestaltet ist, um einen Empfindlichkeitskoeffizienten zu bestimmen, um die Funktionsweise einer Halbleitervorrichtung zu analysieren und bei dem simultane Gleichungen, die aus Elektronen- und Loch- Transportgleichungen und der Poisson-Gleichung bestehen, mittels eines Rechners gelöst werden, um die Modellierung der Halbleitervorrichtung zu erreichen, umfasst folgende Schritte:
Umschreiben der simultanen Gleichungen fp = 0, fn = 0 und fψ = 0 in die folgenden Gleichungen (19-1), (19-2), (19- 3), die künstliche Zeitdifferentialglieder dp/dt, dn/dt, dψ/dt und Empfindlichkeitskoeffizienten λp, λn und λψ enthalten:
dp/dt = - λpfp ... (19-1)
dn/dt = λnfn ... (19-2)
dψ/dt = λψfψ ... (19-3)
Bestimmen von Gitter- oder Maschenpunkten (M,N) der Halbleitervorrichtung, wodurch die Gleichungen (19-1), (19-2) und (19-3) in die folgenden transformiert werden:
dp(M,N)/dt = - λp(M,N)fp(M,N) ... (20-1)
dn(M,N)/dt = λn(M,N)fn(M,N) ... (20-2)
dψ(M,N)/dt = λψ(M,N)fψ(M,N) ... (20-3)
wobei fp, fn und fψ wie folgt definiert sind:
fp(M,N) = (1/q) {[Jpx(M) - Jpx(M-1)]/hx'(M)} + (1/q) {[Jpy(N) - Jpy(N-1)]/hy'(M)) + U(M,N) ... (21-1)
fn(M,N) = (1/q) {[Jnx(M) - Jnx(M-1)]/hx'(M)} + (1/q) {[Jny(N) - Jny(N-1)]/hy'(M)} - U(M,N) ... (21-2)
fψ(M,N) = [1/hx'(M)] {[ψ(M+1,N) - ψ(M,N)]/hx(M) - [ψ(M,N) - ψ(M-1,N)]/hx(M-1)} + [1/hy'(N)] {[ψ(M,N+1) - ψ(M,N)]/hy(N) - [ψ(M,N) - ψ(M,N-1)]/hy(N-1)} + (q/&epsi;) [Γ(M,N) + p(M,N) - n(M,N)] ... (21-3)
Zeitintegrieren der Gleichungen (20-1), (20-2) und (20-3), wodurch die Lösungen der simultanen Gleichungen erhalten werden, wobei die Empfindlichkeitskoeffizienten λp, λn und λψ so sind, dass die Eigenwerte der Fehlerfortpflanzungsmatrizen nicht größer als 1 sind.
Bei der obigen Beschreibung des erfindungsgemäßen ersten Verfahrens sei zwecks Einfachheit angenommen, dass die durch dieses Verfahren zu lösende Aufgabenstellung ein stationäres Gleichstromproblem ist, das ursprüngliche Funktionen fp und fn beinhaltet, die keine Differtialwerte für die physikalische Zeit τ enthalten. Um eine beliebige andere Aufgabenstellung mit ursprünglichen Funktionen zu lösen, die Differentialausdrücke der Zeit τ enthalten, ist es ausreichend, die Glieder ∂p/∂τ und ∂n/∂τ, die in den Gleichungen (8-1) und (8-2) enthalten sind, in die folgenden Differenzgleichungen umzuschreiben:
[p(M,N) - p&sub0;(M,N)]/δτ
[n(M,N) - n&sub0;(M,N)]/δτ
Wenn die Glieder ∂p/∂τ und ∂n/∂τ derart umgeschrieben sind, wird es möglich, die simultanen Gleichungen ohne die Notwendigkeit eines Änderns des Wesens des erfindungsgemäßen ersten Verfahrens zu lösen.
Das zweite Verfahren, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung zu analysieren, und bei dem simultane Gleichungen, die aus Elektronen- und Loch-Transportgleichungen und der Poisson-Gleichung bestehen, mittels eines Computers gelöst werden, um die Modellierung der Halbleitervorrichtung zu erreichen, umfasst folgende Schritte:
Umschreiben der simultanen Gleichungen in die folgenden Gleichungen (19-1), (19-2), (19-3), die künstliche Zeitdifferentialglieder dp/dt, dn/dt, dψ/dt, und Empfindlichkeitskoeffizienten λp, λn und λψ enthalten:
dp/dt = - λpfp ... (19-1)
dn/dt = λnfn ... (19-2)
dψ/dt = λψfψ ... (19-3)
wobei
fp = (1/q)divJp - (G-U) ... (30-1)
fn = (1/q)divJn + (G-U) ... (30-2)
fψ = divgradψ + (q/&epsi;)(Γ+p-n) ... (30-3)
n = niexp [θ(ψ-φn)] ... (30-4)
p = niexp [θ(φp-ψ)] ... (30-5)
θ = q/(kT) ... (30-6)
Anwenden der Gleichungen (30-4) und (30-5), womit p und n in den linken Seiten der Gleichungen (19-1) und (19-2) eliminiert werden, und Umschreiben von Gleichungen (19-1) und (19-2) in die folgenden Gleichungen (33-1) und (33-2) ergibt:
dφp/dt = - λpfp/(θp) + λψfψ ... (33-1)
dφn/dt = - λnfn/(θn) + λψfψ ... (33-2)
Bestimmen von Maschenpunkten (M,N) der Halbleitervorrichtung und Bestimmen der Empfindlichkeitskoeffizienten λp, λn und λψ derart, dass die Eigenwerte der Fehlerfortpflanzungsmatrix nicht größer als 1 sind, wodurch die simultanen Gleichungen gelöst werden.
Das dritte Verfahren, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung zu analysieren, und bei dem simultane Gleichungen, die aus Elektronen- und Loch-Transportgleichungen und der Poisson-Gleichung bestehen, mittels eines Computers gelöst werden, um die Modellierung der Halbleitervorrichtung zu erreichen, umfasst folgende Schritte:
Umschreiben der simultanen Gleichungen in die folgenden Gleichungen (19-1), (19-2), (19-3), die künstliche Zeitdifferentialglieder dp/dt, dn/dt, dψ/dt und Empfindlichkeitskoeffizienten λp, λn und λψ enthalten:
dp/dt = - λpfp ... (19-1)
dn/dt = λnfn ... (19-2)
dψ/dt = λψfψ ... (19-3)
wobei
fp = (1/q)divJp - (G-U) ... (30-1)
fn = (1/q)divJn + (G-U) ... (30-2)
fφ = divgradψ + (q/&epsi;)(Γ+p-n) ... (30-3)
n = niexp [θ(ψ-φn)] ... (30-4)
p = niexp [θ(φp-ψ)] ...(30-5)
θ = q/(kT) ... (30-6)
Anwenden von Gleichungen (30-4) und (30-5), wodurch p und n in den linken Seiten der Gleichungen (19-1) und (19-2) eliminiert werden, und Umschreiben von Gleichungen (19-1) und (19-2) in die folgenden Gleichungen (33-1) und (33-2):
dφp/dt = - λpfp/(θp) + λψfψ ... (33-1)
dφn/dt = - λnfn/(θn) + λψfψ ...(33-2)
Das Eliminieren des zweiten Glieds der Gleichung (33-1) und ebenfalls des zweiten Glieds der Gleichung (33-2), wodurch die Gleichungen (33-1) und (33-2) in die folgenden Gleichungen (34-1) und (34-2) transformiert werden:
dφp/dt = - λpffp/(θp) ... (34-1)
dφn/dt = - λnfn/(θn) ... (34-2)
Bestimmen von Maschenpunkten (M,N) der Halbleitervorrichtung und Bestimmen der Empfindlichkeitskoeffizienten λp, λn und λψ derart, dass die Eigenwerte der Fehlerfortpflanzungsmatrix nicht größer als 1 sind, wodurch die simultanen Gleichungen gelöst werden.
Das vierte Verfahren, das ausgestaltet ist, um die Funktionsweise der Halbleitervorrichtung zu analysieren, ist dadurch gekennzeichnet, dass das zweite oder dritte Verfahren während einer ersten Zeitspanne verwendet wird, und das andere der zweiten und dritten Verfahren während einer zweiten Zeitspanne verwendet wird, um dadurch die Zeitintegration zum Lösen der simultanen Gleichungen durchzuführen.
Das fünfte Verfahren, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung zu analysieren, und bei dem die simultanen Gleichungen, die aus Elektronen- und Loch-Transportgleichungen und der Poisson-Gleichung bestehen, mittels eines Computers gelöst werden, um die Modellierung der Halbleitervorrichtung durchzuführen, umfasst folgende Schritte:
Umschreiben der simultanen Gleichungen fp = 0, fn = 0 und fψ = 0 in die folgenden Gleichungen (19-1), (19-2), (19-3), die künstliche Zeitdifferentialglieder dp/dt, dn/dt, dt/dt und Empfindlichkoeffizienten λp, λn und λψ enthalten:
dp/dt = - λpfp ... (19-1.)
dn/dt = λnfn ... (19-2)
dψ/dt = λψfψ ... (19-3)
wobei
fp = (1/q)divJp - (G-U) ... (30-1)
fn = (1/q)divJn + (G-U) ... (30-2)
fψ = divgradψ + (q/&epsi;)(Γ+p-n) ... (30-3)
n = niexp [θ(ψ-φn)] ... (30-4)
p = niexp [θ(φp-y)] ... (30-5)
θ = q/(kT) ... (30-6)
Bestimmen von Maschenpunkten (M,N) der Halbleitervorrichtung und Übertragen von räumlicher Positionsabhängigkeit, die für die physikalischen Eigenschaften der Halbleitervorrichtung geeignet ist, auf die Empfindlichkeitskoeffizienten λP, λn und λψ, wodurch Gleichungen (19-1), (19-2) und (19-3) in die folgenden transformiert werden:
dp(M,N)/dt = - λp(M,N)fp(M,N) ... (20-1)
dn(M,N)/dt = λn(M,N)fn(M,N) ... (20-2)
dψ(M,N)/dt = λψ(M,N)fψ(M,N) ... (20-3)
wobei fp, fn und fψ wie folgt definiert sind:
λp(M,N)δt = 2ωp/ {( ∂fp(M,N)/∂p(K,L)) + ∂fp(M,N)/∂p(M,N)} ... (32-1)
λn(M,N)δt = 2ωn/ {( ∂fn(M,N)/∂n(K,L)) - ∂fn(M,N)/∂n(M,N)} ... (32-2)
λψ(M,N)δτ = 2ωψ/ {( ∂fψ(M,N)/∂ψ(K,L)) - ∂fψ(M,N)/∂ψ(M,N)} ... (32-3)
wobei δt die diskreten Zeitintervalle auf der Zeitachse, (MN) die Nummer von beliebigen Maschenpunkten im diskreten Raum, und (K,L) die Nummer jedes beliebigen benachbarten Punkt, der um den Raummaschenpunkt (M,N) existiert, ist.
Zeitintegrieren der Gleichungen (20-1), (20-2) und (20- 3), bis diese Gleichungen einen stationären Zustand erreichen, wodurch die Lösungen der simultanen Gleichungen erhalten werden, wobei die Glieder fp(M,N), fn(M,N) und fφ (M,N) auf den rechten Seiten der Gleichungen (20-1), (20-2) und (20-3) wie folgt festgelegt sind:
fp(M,N) = (1/q) {[Jpx(M) - Jpx(M-1)]/hx'(M)} + (1/q) {[Jpy(N) - Jpy(N-1)]/hy'(M)} + U(M,N) ... (26-1)
fn(M,N) = (1/q) {[Jnx(M) - Jnx(M-1)]/hx'(M)} + (1/q) {[Jny(N) - Jny(N-1)]/hy'(M)} - U(M,N) ... (26-2)
fψ(M,N) = [1/hx'(M)]{[ψ(M+1,N) - ψ(M,N)]/hx(M) - [ψ(M,N) - ψ(M-1,N)]/hx(M-1)} + [1/hy'(N)] {[ψ(M,N+1) - ψ(M,N)]/hy(N) - [ψ(M,N) - ψ(M,N-1)]/hy(N-1)} + (q/&epsi;) [Γ(M,N) + p(M,N) - n(M,N)] ... (26-3)
Bei der obigen Beschreibung des erfindungsgemäßen fünften Verfahrens wird angenommen, dass die durch dieses Verfahren zu lösende Aufgabenstellung zwecks Einfachheit ein stationäres Gleichstromproblem ist, das ursprüngliche Funktionen fp und fn enthält, die keine Differentialgröße für die physikalische Zeit τ enthalten. Um eine beliebige andere Aufgabenstellung zu lösen, die die ursprünglichen Funktionen einschließt, die Differentialgröße für die Zeit τ enthalten, ist es ausreichend, die in den Gleichungen (8-1) und (8-1) enthaltenden Glieder ∂p/∂τ und ∂n/∂t in die folgenden Differenzausdrücke umzuschreiben:
[p(M,N) - p&sub0;(M,N)]/δτ
[n(M,N) - n&sub0;(M,N)]/δτ
Wenn die Glieder ∂p/∂τ und ∂n/∂τ so umgeschrieben sind, wird es möglich, die simultanen Gleichungen ohne die Notwendigkeit des Änderns des Wesens des erfindungsgemäßen ersten Verfahrens zu lösen.
Das sechste Verfahren, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung zu analysieren, ist mit den zweiten, dritten und vierten Verfahren mit der Ausnahme identisch, dass die Empfindlichkeitskoeffizienten λp, λn und λψ aus den folgenden Gleichungen erhalten werden:
Das erste Gerät, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung zu analysieren, und bei der die simultanen Gleichungen, die aus Elektronen- und Loch-Transportgleichungen und der Poisson-Gleichung bestehen, in Differentialgleichungen dp/dt = -λpfp, dn/dt = λnfn und dψ/dt = λψfψ umgeschrieben werden, die die Zeitdifferentialglieder und Empfindlichkeitskoeffizienten enthalten, und wobei die Differentialgleichungen gelöst werden, um die Modellierung der Halbleitervorrichtung zu erreichen, umfasst:
eine logische Zelle oder m logischen Zellen;
m Verstärker mit einem einstellbaren Verstärkungsfaktor λi oder einem den m Verstärkern äquivalenten Multiplizierer, die mit den logischen Zellen verbunden ist und
Detektornetzwerkmittel, die mit den logischen Zellen verbunden sind, zum Erfassen, dass die logischen Zellen eine Einschwingstabilität (transient stability) erreichen,
wobei der Verstärkungsfaktor λi automatisch bestimmt wird.
Das zweite Gerät, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung durch das zweite Verfahren zu analysieren, wobei simultane Gleichungen, die aus Elektronen- und Loch-Transportgleichungen und der Poisson-Gleichung bestehen, in die Differentialgleichungen dφp/dt = -λpfp/(θp) + λψfψ, dφn/dt = λnfn/(θn) + λψfψ, δψ/dt = λψfψ umgeschrieben werden, die Zeitdifferentialglieder und Empfindlichkeitskoeffizienten enthalten, umfasst:
erste Mittel zum Kennzeichnen der Adressen von Maschenpunkten der Halbleitervorrichtung;
zweite Mittel zum Durchführen der Integration in Übereinstimmung mit φp, φn, ψ, λp, λn und λψ entsprechend einem beliebigen Maschenpunkt, dessen Adresse durch das erste Mittel gekennzeichnet wurde, wodurch φp, φn und ψ an dem Maschenpunkt hergeleitet werden;
dritte Mittel zum Speichern von so hergeleiteten φp, φn und ψ,; und
vierte Mittel zum Bestimmen, ob durch das zweite Mittel hergeleitete φp, φn und ψ innerhalb erlaubter Fehlergrenzen fallen oder nicht, zum Zurückführen von φp, φn und ψ an das zweite Mittel, so dass das zweite Mittel die Integration erneut durchführt, wenn φp, φn und ψ außerhalb der erlaubten Fehlergrenzen fallen, und zum Ausgeben von φp, φn und ψ, wenn φp, φn und ψ innerhalb der erlaubten Fehlergrenzen fallen.
Das dritte Gerät, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung durch das dritte Verfahren zu analysieren, wobei simultane Gleichungen, die aus Elektronen- und Loch-Transportgleichungen und der Poisson-Gleichung bestehen, in Differentialgleichungen dφp/dt = -λpfp/(θp), dφn/dt = -λnfn/(θn), dψ/dt = λψfψ umgeschrieben werden, die Zeitdifferentialglieder und Empfindlichkeitskoeffizienten enthalten, umfasst:
erste Mittel zum Kennzeichnen der Adressen Von Maschenpunkten der Halbleitervorrichtung;
zweite Mitteln zum Durchführen der Integration in Übereinstimmung mit φp, φn, ψ, λp, λn und λψ, entsprechend einem beliebigen Maschenpunkt, dessen Adresse durch das erste Mittel gekennzeichnet wurde, wodurch φp, φn und ψ an dem Maschenpunkt hergeleitet werden;
dritte Mittel zum Speichern der so hergeleiteten φp, φn und ψ; und
vierte Mittel zum Bestimmen, ob durch das zweite Mittel hergeleitete φp, φn und ψ innerhalb erlaubter Fehlergrenzen fallen oder nicht, zum Zurückführen von φp, φn und ψ zu dem zweiten Mittel, so dass das zweite Mittel die Integration erneut durchführt, wenn φp, φn und ψ außerhalb der erlaubten Fehlergrenzen fallen, und zum Ausgeben von φp, φn und ψ, wenn φp, φn und ψ innerhalb der erlaubten Fehlergrenzen fallen.
Das vierte Gerät, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung durch das vierte Verfahren zu analysieren, wobei simultane Gleichungen, die aus Elektronen- und Loch-Transportgleichungen und der Poisson-Gleichung bestehen, in Differentialgleichungen dφp/dt = -λpfp/(θp) + λψfψ, dφn/dt = -λnfn/(θn) + λψfψ, dψ/dt = λψfψ, die Zeitdifferenzglieder und Empfindlichkeitskoeffizienten enthalten, oder in Differentialgleichungen dφp/dt = - λpfp/(θp), dφn/dt = -λnfn/(θn), dψ/dt = λψfψ, die Zeitdifferentialglieder und Empfindlichkeitskoeffizienten enthalten, umgeschrieben werden, umfasst:
erste Mittel zum Kennzeichnen der Adressen von Maschenpunkten der Halbleitervorrichtung;
zweite Mittel zum Durchführen der Integration in Übereinstimmung mit φp, φn, ψ, λp, λn, und λψ, die jedem Maschenpunkt entsprechen, dessen Adresse durch das erste Mittel gekennzeichnet wurde, durch Verwenden eines der neunten und zehnten Verfahren während einer ersten Zeitspanne, und durch Verwenden des anderen der neunten und zehnten Verfahren während einer zweiten Zeitspanne, wodurch φp, φn und ψ an dem Maschenpunkt hergeleitet werden;
dritte Mittel zum Speichern der somit hergeleiteten φp, φn und ψ; und
vierte Mittel zum Bestimmen, ob durch das zweite Mittel hergeleitete φp, φn und ψ innerhalb erlaubter Fehlergrenzen fallen oder nicht, zum Zurückführen von φp, φn und ψ an das zweite Mittel, so dass das zweite Mittel die Integration erneut durchführt, wenn φp, φn und ψ außerhalb der erlaubten Fehlergrenzen fallen, und zum Ausgeben von φp, φn und ψ, wenn φp, φn und ψ innerhalb der erlaubten Fehlergrenzen fallen.
Das fünfte Gerät, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung durch das fünfte Verfahren zu analysieren, wobei aus Elektronen- und Loch- Transportgleichungen und der Poisson-Gleichung bestehenden simultanen Gleichungen in Differentialgleichungen dp/dt = -λpfp, dn/dt = λnfn, dψ/dt = λψfψ umgeschrieben werden, die Zeitdifferentialglieder und Empfindlichkeitskoeffizienten enthalten, umfasst:
erste Mittel zum Kennzeichnen der Adressen von Maschenpunkten der Halbleitervorrichtung;
zweite Mittel zum Durchführen der Integration in Übereinstimmung mit p, n, ψ, λp, λn und λψ, entsprechend jedem Maschenpunkt, dessen Adresse durch das erste Mittel gekennzeichnet wurde, wodurch p, n und ψ an dem Maschenpunkt hergeleitet werden; und
dritte Mittel zum Bestimmen, ob durch das zweite Mittel hergeleiteten p, n und ψ innerhalb erlaubter Fehlergrenzen fallen, zum Zurückführen von p, n und ψ an das zweite Mittel, so dass das zweite Mittel die Integration erneut durchführt, wenn p, n und ψ außerhalb der erlaubten Fehlergrenzen fallen, und zum Ausgeben von p, n und ψ, wenn p, n und ψ innerhalb der erlaubten Fehlergrenzen fallen.
Das sechste Gerät, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung durch das sechste Verfahren zu analysieren, wobei aus Elektronen- und Loch- Transportgleichungen und der Poisson-Gleichung bestehende simultane Gleichungen entweder in Differentialgleichungen dφp/dt = -λpfp/(θp) + λψfψ, dφn/dt = -λnfn/(θn) + λψfψ, dψ/dt = λψfψ, die Zeitdifferentialglieder und Empfindlichkeitskoeffizienten enthalten, oder in Differentialgleichungen dφp/dt = -λpfp/(θp), dφn/dt = - λnfn/(θn), dψ/dt = λψfψ, die Zeitdifferentialglieder und Empfindlichkeitskoeffizienten enthalten, oder in beide Sätze von simultanen Gleichungen umgeschrieben werden, umfasst:
erste Mittel zum Kennzeichnen der Adressen von Maschenpunkten der Halbleitervorrichtung;
zweite Mittel zum Durchführen der Integration in Übereinstimmung mit φp, φn, ψ, λp, λn und λψ, entsprechend einem beliebigem Maschenpunkt, dessen Adresse durch das erste Mittel gekennzeichnet wurde, wodurch φp, φn und ψ an dem Maschenpunkt hergeleitet werden;
dritte Mittel zum Speichern der so hergeleiteten φp, φn und ψ; und
vierte Mittel zum Bestimmen, ob durch das zweite Mittel hergeleitete φp, φn und ψ innerhalb erlaubter Fehlergrenzen fallen oder nicht, zum Zurückführen von φp, φn und ψ an das zweite Mittel, so dass das zweite Mittel die Integration erneut durchführt, wenn φp, φn und ψ außerhalb der erlaubten Fehlergrenzen fallen, und zum Ausgeben von φp, φn und ψ, wenn φp, φn und ψ innerhalb der erlaubten Fehlergrenzen fallen.
Wie es beschrieben wurde, werden bei dem ersten bis sechsten erfindungsgemäßen Verfahren die aus Elektronen- und Loch-Transportgleichungen und der Poisson-Gleichung bestehende simultane Gleichungen in Differentialgleichungen umgeschrieben, die Differentialglieder und Empfindlichkeitskoeffizienten λ enthalten, und die Empfindlichkeitskoeffizienten λ werden gemäß der Gleichungen (32-1), (32-2) und (32-3) oder (35-1), (35-2) und (35-3) eingestellt. Folglich können diese Verfahren die Matrixaufgabenstellung mit hoher Geschwindigkeit ohne externe Steuerung sogar dann lösen, wenn die beteiligten Gleichungen in dem Problem "massiv" bzw. sehr zahlreich sind.
Die ersten bis sechsten erfindungsgemäßen Geräte bauen Schaltungsnetzwerke in Übereinstimmung mit den ersten bis sechsten Analyseverfahren auf. Die so aufgebauten Schaltungsnetzwerke werden betrieben, bis die simultanen Gleichungen stationäre Zustände erreichen, wodurch die Lösungen der Gleichungen erhalten werden. Da die Empfindlichkeitskoeffizienten λ jeweils gemäß spezifischer Gleichungen bestimmt werden können, können die ersten bis sechsten Geräte die Lösungen der simultanen Gleichungen mit hoher Geschwindigkeit erhalten.
Das siebente Verfahren, das ausgestaltet ist, um einen Empfindlichkeitskoeffizienten zu bestimmen, um die Funktionsweise einer Halbleitervorrichtung zu analysieren, und bei dem eine die Beziehung zwischen dem elektrischen Potential und der Raumladung darstellende Poisson-Gleichung mittels eines Computers gelöst wird, um die Modellierung der Halbleitervorrichtung zu erreichen, umfasst folgende Schritte:
Umschreiben der Poisson-Gleichung in die folgende Gleichung (22), die Zeitdifferentialglieder ∂ψ/∂t und einen Empfindlichkeitskoeffizienten λ aufweisen,
∂ψ/∂t = λfψ ... (22)
Bestimmen von Maschenpunkten (M,N) der Halbleitervorrichtung, wodurch die Gleichung (22) in die folgende Gleichung (23) transformiert wird:
∂ψ(M,N)/∂t = λ(M,N)fψ(M,N) ... (23)
Zeitintegrierung der Gleichung (23), um die Poisson- Gleichung zu lösen, wobei der Empfindlichkeitskoeffizient A so bestimmt wird, dass die Eigenwerte der Fehlerfortpflanzungsmatrix nicht mehr als 1 sind.
Das achte Verfahren zum Analysieren einer die Beziehung zwischen elektrischem Potential und Raumladung darstellenden Poisson-Gleichung mittels eines Computers, um dadurch die Modellierung einer Halbleitervorrichtung zu erreichen, umfasst folgende Schritte:
Umschreiben der Poisson-Gleichung in die folgende Gleichung (22), die Zeitdifferentialglieder ∂ψ/∂t und einen Empfindlichkeitskoeffizienten λ enthält.
∂ψ/∂t = λψfψ ... (22)
Bestimmen von Maschenpunkten (M,N) der Halbleitervorrichtung und Auferlegen einer Raumpositionsabhängigkeit, die für die physikalischen Eigenschaften der Halbleitervorrichtung geeignet ist, wodurch die Gleichung (22) in die folgende Gleichung (23) transformiert wird:
∂ψ(M,N)/∂t = λ(M,N)fψ(M,N) ... (23)
Bestimmen der Empfindlichkeitskoeffizienten λ(M,N) durch die folgende Gleichung (54-11):
λ(M,N)δt ≤ 2/(2H + α) ... (54-11)
wobei δt das diskrete Zeitintervall auf der Zeitachse und (M,N) die Nummer eines beliebigen diskreten Raummaschenpunkts ist;
Zeitintegrieren der Gleichungen (23), bis diese Gleichungen einen stationären Zustand erreichen, wodurch die Lösungen der Poisson-Gleichung erhalten werden.
Das siebente Gerät, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung zu analysieren, und bei dem eine die Beziehung zwischen dem elektrischen Potential und der Raumladung darstellender Poisson-Gleichung in eine Differentialgleichung dψ/dt = λfψ umgeschrieben wird, und wobei die Differentialgleichung gelöst wird, um die Modellierung der Halbleitervorrichtung zu erreichen, umfasst:
eine logischen Zelle bzw. eine Logikzelle oder m logische Zellen;
m Verstärker mit einem einstellbaren Verstärkungskoeffizienten λi oder einen den m Verstärkern äquivalenten Multiplizierer, der mit den logischen Zellen verbunden ist und
Detektornetzwerkmittel, die mit den logischen Zellen verbunden sind, zum Erfassen, dass die logischen Zellen eine Einschwingstabilität (transient stability) erreichen,
wobei der Verstärkungskoeffizient λi automatisch bestimmt wird.
Das achte Gerät zum Analysieren der Poisson-Gleichung durch das bei dem achten Verfahren definierte Verfahren durch Umschreiben der Poisson-Gleichung in Differentialgleichungen, ∂ψ/∂t = λψfψ, die ein Zeitdifferentialglied und einen Empfindlichkeitskoeffizienten enthält, und dann Lösen der Differentialgleichungen, um dadurch die Modellierung der Halbleitervorrichtungen zu erreichen, umfasst:
erste Mittel zum Kennzeichnen der Adressen von Maschenpunkte der Halbleitervorrichtung;
zweite Mittel zum Durchführen der Integration in Übereinstimmung mit ψ und λψ, die einem beliebigen Maschenpunkt entsprechen, dessen Adresse durch das erste Mittel gekennzeichnet wurde, wodurch ψ an dem Maschenpunkt hergeleitet wird;
dritte Mittel zum Speichern des durch das zweite Mittel hergeleiteten ψ; und
vierte Mittel zum Bestimmen, ob das durch das zweite Mittel hergeleitete ψ innerhalb einer erlaubten Fehlergrenze fällt, zum Zurückführen von ψ an das zweite Mittel, so dass das zweite Mittel die Integration erneut durchführt, wenn ψ außerhalb der erlaubten Fehlergrenze fällt, und zum Ausgeben von ψ, wenn ψ innerhalb der erlaubten Fehlergrenze fällt.
Das neunte Verfahren, das ausgestaltet ist, um einen Empfindlichkeitskoeffizienten zu bestimmen, um ein spezifisches physikalisches Phänomen zu analysieren, und bei dem entweder eine einzelne Gleichung oder simultane Gleichungen f(x) = 0 mittels eines Computers gelöst werden, umfasst folgende Schritte:
Umschreiben der Gleichung oder der simultanen Gleichungen in die folgende Gleichung (25), die ein Zeitdifferentialglied dx/dt und einen Empfindlichkeitskoeffizienten λ enthält:
dx/dt = λf(x) ... (25)
Zeitdifferenzieren der Gleichung (25), wodurch die Lösung der Gleichung oder der simultanen Gleichungen erhalten wird, wobei der Empfindlichkeitskoeffizient λ der Art bestimmt wird, dass die Eigenwerte der Fehlerfortpflanzungsmatrizen nicht größer als 1 sind.
Bei dem neunten Verfahren kann die Vektorfunktion f Operatoren für eine unbekannte Größe x und die Funktion von x, wie beispielsweise ein Differential-Operator, einen Integral-Operator und weitere Operatoren enthalten. Alternativ kann die Vektorfunktion f Operatoren für die physikalische Zeit τ und die Funktion von τ enthalten, wie beispielsweise einen Differential-Operator, einen Integral- Operator und weitere Operatoren.
Das neunte Gerät, das ausgestaltet ist, um eine bestimmtes physikalisches Phänomen zu analysieren, und bei der entweder eine einzelne Gleichung oder simultane Gleichungen f(x) = 0 in eine Differentialgleichung dx/dt = λf(x) umgeschrieben werden und die Differentialgleichung gelöst wird, umfasst:
einen analogen oder digitalen Berechnungsabschnitt zum Durchführen eines Integrationsverfahrens auf der Zeitachse; und
Mittel zum Einstellen und Optimieren der Größe einer Koeffizientenmatrix λ, während der Berechnungsabschnitt den Integrationsprozess durchführt,
wobei der Verstärkungskoeffizient λi automatisch bestimmt wird.
Das zehnte Verfahren, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung zu analysieren, und bei dem eine Aufgabenstellung eines Bestimmens von äquivalenten Schaltungskonstanten der Halbleitervorrichtung durch das Verfahren des kleinsten Fehlerquadrats mittels eines Computers gelöst wird, umfasst folgende Schritte:
Umschreiben des Problems des Bestimmens der äquivalenten Schaltungskonstanten in die folgende Gleichung (24), die ein Zeitdifferentialglied dq/dt und einen Empfindlichkeitskoeffizienten λ enthält:
dq/dt = - λ(Aq - b) ... (24)
Auferlegen von unterschiedlichen Werten, die für die durch die äquivalenten Schaltungskonstanten der Halbleitervorrichtung dargestellten physikalischen Eigenschaften geeignet sind, an den Empfindlichkeitskoeffizienten λ; und
Zeitintegrieren der Gleichung (24), wodurch die Lösung der Aufgabenstellung erhalten wird.
Das zehnte Gerät, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung zu analysieren, und bei dem eine Aufgabenstellung des Bestimmens der äquivalenten Schaltungskonstanten der Halbleitervorrichtung in eine Differentialgleichung dq/dt = - λ(Aq - b) umgeschrieben wird, und die Differentialgleichung durch das Verfahren des kleinsten Fehlerquadrats gelöst wird, umfasst:
eine logische Zelle oder m logische Zellen;
m Verstärker mit einem einstellbaren Verstärkungskoeffizienten λi oder einem den m Verstärkern äquivalenten Multiplizierer, der mit den logischen Zellen verbunden ist; und
Detektornetzwerkmittel, die mit den logischen Zellen verbunden sind, zum Erfassen, dass die logischen Zellen eine Einschwingstabilität (transient stability) erreichen.
Wie es beschrieben wurde, werden bei den ersten bis zehnten erfindungsgemäßen Verfahren die Gleichungen, die gelöst werden müssen (d. h. die simultanen Gleichungen, die aus Elektronen- und Loch-Transportgleichungen und der Poisson-Gleichung aufgebaut sind, wobei die Poisson-Gleichung die Beziehung zwischen dem elektrischen Potential und der Raumladung darstellt, und die Aufgabenstellung des Bestimmens der äquivalenten Schaltungskonstanten der Halbleitervorrichtung, und dem Problem des Analysierens eines spezifischen Phänomens) in eine Differentialgleichung umgeschrieben, die Differentialglieder und Empfindlichkeitskoeffizienten λ enthält, und die Empfindlichkeitskoeffizienten λ werden auf Werte eingestellt, die für die physikalischen Eigenschaften der zu lösenden Gleichungen geeignet sind. Folglich können diese Verfahren die Matrixaufgabenstellungen mit hoher Geschwindigkeit sogar dann lösen, wenn die beteiligten Gleichungen in den Problemen massiv bzw. sehr zahlreich sind.
Bei den ersten bis zehnten erfindungsgemäßen Geräten wird die unabhängige Variable t in der Gleichung (25) als Zeit angesehen, und die unbekannte Größe x ist die an einen bestimmten Anschluss eines elektrischen Schaltungsnetzwerks anliegenden Spannung, die Gleichung (25) erfüllt. Das Einschwingphänomen (transient phenomenon) des elektrischen Schaltungsnetzwerks entwickelt sich mit der Zeit, so dass es stabil wird. Wenn das Phänomen stabil wird, ist das Gerät insgesamt so aufgebaut, dass es die Lösung der Gleichung (25) liefert. Ferner kann das Gerät, um die ersten bis zehnten Verfahren zu erreichen, durch eine digitale Schaltung aufgebaut sein.
Um die Matrixgleichung Ax = b durch Durchführen einer digitalen Berechnung bei den herkömmlichen Verfahren zu lösen, wird die Gleichung in x = A&supmin;¹b umgeschrieben, falls eine Eliminierungstechnik angewendet wird. Diese Berechnung umfasst eine sogenannte Matrixinversionsberechnung (inverse matrix calculation) oder eine im Wesentlichen dazu äquivalente Matrixberechnung. Ein großer Zeitbetrag muss aufgewendet werden, um eine Matrixinversionsberechnung oder eine dazu äquivalente Berechnung zu erreichen.
Bei dem oben beschriebenen Berechnungsverfahren, ob x einen stationären Wert erreichen kann oder nicht, hängt davon ab, ob die ausgewählte diagonale Koeffizientenmatrix λ geeignet ausgewählt ist oder nicht. Im Allgemeinen ist, wenn λ zu klein ist, der Zeitdifferentialbetrag dx/dt ebenfalls zu klein. In diesem Fall ändert sich der Anfangsprobewert x = x&sub0; wenig, und ein großer Zeitbetrag verstreicht unvermeidbar, bevor x den stationären Wert erreicht. Wenn λ umgekehrt zu groß ist, ist der Zeitdifferentialbetrag dx/dt ebenfalls zu groß. In diesem Fall ist das durch die Differentialgleichung (25) dargestellte Verhalten des Systems (oder des Schaltungsnetzwerks) instabil. Bei einem tatsächlichen Beispiel, bei dem die numerische Analyse mittels einer digitalen Simulation durchgeführt wurde, änderte sich die unbekannte Größe x übermäßig, wenn λ zu groß war; bei Verstreichen einer infinitesimalen Zeitspanne begann sich die Größe x in der umgekehrten Richtung zu einem stationären Wert hin zu ändern.
Es ist offensichtlich, dass sich die unbekannte Größe x wiederholt in den Vorwärts- und Rückwärtsrichtungen mit dem Ablauf der Zeit ändert. Mit anderen Worten fluktuiert die Größe x auf der Zeitachse zu und weg von dem stationären Wert, und kann sich nicht auf dem stationären Wert innerhalb einer kurzen Zeitspanne einschwingen. Wie es in der Technik bekannt ist, ist, je mehr sich die Größe x ändert, desto größer ist die Amplitude ihrer Fluktuation. In einigen Fällen kann sich die Größe x unbestimmt in einer Richtung (entweder der positiven Richtung oder der negativen Richtung) ändern, sie kann unbeschränkt divergieren.
Um das obige zusammenzufassen ist es bedeutsam, um das Gerät der Erfindung gemäß dem Verfahren der Erfindung zu betreiben, die Größe der Diagonal-Koeffizientenmatrix λ auf einen solchen Wert einzustellen, der die unbekannte Größe x veranlasst, sich moderat zu ändern.
Es sei bemerkt, dass, da λ eine Diagonal- Koeffizientenmatrix ist, n Koeffizienten (λ&sub1;, λ&sub2;, ..., λn) n Elementen (x&sub1;, x&sub2;, ..., xn) entsprechen, die die Komponenten des Vektors x sind. Somit weist das Verfahren und das Gerät erfindungsgemäß einen Freiheitsgrad auf, dass n Koeffizienten auf Werte eingestellt werden können, die für die n unbekannten Größen-Elemente geeignet sind. Daher kann das Verfahren dieser Erfindung einen stationären Wert der unbekannten Größe zuweisen, und es ist beträchtlich effizienter als das herkömmliche Verfahren, d. h. das Iterationsverfahren, bei dem in den meisten Fällen sogenannte "Relaxationskoeffizienten", jeweils für eine Matrixgleichung, verwendet werden.
Es wird nun erläutert, wie λ zu optimieren ist. Wie es herausgestellt wurde, ist λ eine Diagonal-Matrix mit n Elementen λ&sub1;, λ&sub2;, ..., λn. Sie kann ausgedrückt werden als:
λ = diag(λ&sub1;, λ&sub2;, ... λn) ... (26)
Gleichung (25) muss eine maximale Zeitdifferentialgröße enthalten, deren absoluter Wert dxi/dt max ist, der es dem Gerät ermöglicht, stabil für alle unbekannten Größen xi zu arbeiten, wobei i = 1, 2, ..., n ist. Folglich darf keine Komponente von λi die durch die maximale Zeitdifferentialgröße bestimmten Grenzen überschreiten. Diese Bedingung wird wie folgt dargestellt:
λi ≤ (1/ aijxj - bi ) dxi/dt max
i = 1, 2, ..., n
... (27)
Die Betrag oder die Magnitude der rechten Seite der Gleichung (27) hängt von der Art der zu lösenden Aufgabenstellung ab. Gemäß der Matrixalgebra bezieht sich die Bedingung, die auf das Bestimmen des Betrags auf der rechten Seite anwendbar ist, auf den Eigenwert der Gegenstandsmatrix. Das Gerät dieser Erfindung umfasst Verstärker für m logische Zellen. Wie es später ausführlich beschrieben wird, umfasst jeder der Verstärker einen Verstärkungskoeffizienten λi, der groß genug ist, um tatsächliche Aufgabenstellungen zu lösen, und der mittels einer Wahlvorrichtung (Dial) oder eines Controllers eingestellt werden kann. Es ist ferner wünschenswert, dass m logische Zellen mit einem jeweiligen Multiplizierer verbunden sind, dessen Eingangssignal mit λi durch den Multiplizierer multipliziert wird.
Bei den Verfahren wird die Parameterfunktion λ aus den Stabilitätsbedingungen bestimmt. Definitionsgemäß ermöglichen es die Stabilitätsbedingungen, den Akkumulierungseffekt von Abrundungsfehlern zu vernachlässigen, die während einer numerischen Analyse durchgeführt wurden. (Siehe G. D. Smith, Solutions of Partial Differential Equations by Means of Computers, Science Inc., 1979, Seite 63.) Um den Parameter A zu bestimmen, wird x&sup0; + δx für die Variable x in der Gleichung (25) eingesetzt, die Gleichung für δx linearisiert und das nicht δx enthaltende Glied wird als eine Konstante eliminiert, wodurch die Gleichung (25) in die folgende Gleichung (26) transformiert wird:
d(δx)/dt = λ(δf)(δx) ... (26)
Ein Diskretisieren der Gleichung (26) ergibt:
δxt + dt = λdt(δf)δxt + δxt ... (27)
Die den Übergang von δxt zu δxt + dt darstellende Gleichung kann in die folgende Matrixgleichung (28) umgeschrieben werden:
δxt+dt = Aδxt ... (28)
wobei A eine als "Fehlerfortpflanzungsmatrix" bekannte Koeffizientenmatrix ist. Die Funktion λ erscheint in der Gleichung in Form des Produkts λdt (dt: Zeitinkrement). Somit wird das Bestimmen des Parameters λ durch die Art und Weise des Bestimmens des Zeitinkrements dt beeinflusst. Um die Stabilität eines beliebigen Systems zu erläutern, müssen wir das Gerschgorin'sche Kreistheorem betrachten. Dieses mathematische Theorem ist gemäß G. D. Smith, Solutions of Partial Differential Equations by Means of Computers, Science Inc., 1979, Seite 68, wie folgt definiert.
Es sei angenommen, dass ein beliebiger Eigenwert einer Matrix n-ter Ordnung A = (aij) ein Element des folgenden geschlossenen Schaltungssatzes Gi ist, vorausgesetzt A = G + F, und F = (fij), G = diag (a&sub1;&sub1;, a&sub2;&sub2;, ..., ann):
u Gi
wobei Gi = {z (komplexe Zahl); z - aii ≤ ri ist.
ri = Σ fij ... (29)
Dieses Theorem definiert die Verteilung der Eigenwerte einer Matrix in einer komplexen Ebene, wie es in Fig. 3 dargestellt ist. Diese Theorem wird auf die Matrixgleichung (28) angewendet. Die Bedingung für die Stabilität, d. h. der Betrag des Eigenwerts sollte 1 oder kleiner sein, wird erfüllt, wenn die durch die Gleichung (29) definierten Kreise innerhalb des Einheitskreises in der komplexen Ebene existieren. Die Gleichung zum Bestimmen von λ wird in eine positive Form umgeschrieben, um die Stabilitätsbedingung zu erfüllen, und ein Schaltungsnetzwerkgerät wird aufgebaut. Was hier bedeutsam ist, ist, dass dieses Schaltungsnetzwerkgerät den Verstärkungsfaktor aufweist, der λ entspricht. Wenn die Gleichung nicht in eine positive Form umgeschrieben werden kann, kann λ durch einen numerischen Prozess bestimmt werden. Das Gerät kann einen Multiplizierer aufweisen, der einem Verstärker äquivalent ist.
Um das Gerät dieser Erfindung zusammenzubauen, werden m Verstärker für m logische Zellen vorgesehen. Jeder der Verstärker ist ausgestaltet, um einen Verstärkungskoeffizienten λi aufzuweisen, der eingestellt werden kann, und ist so groß, dass das Gerät die tatsächlichen einzelnen Probleme lösen kann. Um den Verstärkungskoeffizienten λi einzustellen, weist das Gerät eine Wahlvorrichtung (Dial) oder einen Controller auf. Es ist ferner wünschenswert, dass m logische Zellen mit einem jeweiligen Multiplizierer verbunden sind, dessen Eingangssignal durch den Multiplizierer mit λi multipliziert wird.
Diese Erfindung kann vollständiger aus der folgenden ausführlichen Beschreibung in Verbindung mit den beigefügten Zeichnungen verstanden werden, in denen zeigt bzw. zeigen:
Fig. 1 ein Diagramm, das die Maschenaufteilung des Raums einer Halbleitervorrichtung darstellt;
Fig. 2 ein Diagramm, das schematisch eine Matrix- Vektorgleichung darstellt, die aus Elektronen- und Loch- Transportgleichungen und der Poisson-Gleichung aufgebaut ist;
Fig. 3 ein Diagramm zum Erläutern des Gerschgorin'schen Kreistheorems;
Fig. 4 ein Ablaufdiagramm, das ein Verfahren zum Analysieren der Funktionsweise einer Halbleitervorrichtung gemäß einer ersten Ausführungsform dieser Erfindung erläutert;
Fig. 5 ein Ablaufdiagramm, das ein Verfahren zum Analysieren der Funktionsweise einer Halbleitervorrichtung gemäß einer zweiten Ausführungsform dieser Erfindung erläutert;
Fig. 6 ein Beispiel von Schaltbildern, die ein Gerät zum Durchführen des Verfahrens gemäß der zweiten Ausführungsform der Erfindung zeigen;
Fig. 7A ein Diagramm, das eine Schaltung zum Bestimmen der Äquivalent-Schaltungskonstante einer Halbleitervorrichtung darstellt;
Fig. 7B ein Äquivalent-Schaltbild der in Fig. 7A dargestellten Schaltung;
Fig. 8 ein Ablaufdiagramm, das ein Verfahren zum Bestimmen der Äquivalent-Schaltungskonstanten einer Halbleitervorrichtung gemäß einer dritten Ausführungsform der Erfindung erläutert;
Fig. 9A bis 9D grafische Darstellungen, die vier in komplexen Ebenen aufgetragene y-Parameter zeigen, wobei die y-Parameter den wahren und den Probewerten entsprechen;
Fig. 10 ein Beispiel von Schaltbildern, die ein Gerät zum Durchführen des Verfahrens gemäß der dritten Ausführungsform der Erfindung darstellen;
Fig. 11 ein Ablaufdiagramm, das ein Verfahren zum Analysieren der Funktionsweise einer Halbleitervorrichtung gemäß einer vierten Ausführungsform dieser Erfindung erläutert;
Fig. 12A bis 12C Diagramme, die ein Objekt darstellen, das von dem Verfahren gemäß der vierten Ausführungsform der Erfindung analysiert wird;
Fig. 13A bis 13F Diagramm, die darstellen, wie eine numerische Konvergenz stattfindet, während das in den Fig. 12A bis 12C dargestellte Objekt durch das Verfahren gemäß der vierten Ausführungsform der Erfindung analysiert wird;
Fig. 14 ein Diagramm, das zeigt, wie eine numerische Konvergenz auftritt, während ein umgekehrt vorgespanntes Objekt durch das Verfahren gemäß der vierten Ausführungsform dieser Erfindung berechnet wird;
Fig. 15 ein Ablaufdiagramm, das ein Verfahren zum Analysieren der Funktionsweise einer Halbleitervorrichtung gemäß einer fünften Ausführungsform der Erfindung erläutert;
Fig. 16 ein Ablaufdiagramm, das ein Verfahren zum Analysieren der Funktionsweise einer Halbleitervorrichtung gemäß einer sechsten Ausführungsform der Erfindung erläutert;
Fig. 17 ein Blockdiagramm, das ein Gerät zum Durchführen des funktions-analysierenden Verfahrens gemäß der sechsten Ausführungsform der Erfindung zeigt;
Fig. 18 ein Blockdiagramm, das die in dem in Fig. 17 gezeigten Gerät enthaltene Maschenpunkteinheit darstellt;
Fig. 19 ein Ablaufdiagramm, das ein Verfahren zum Analysieren der Funktionsweise einer Halbleitervorrichtung gemäß einer siebenten Ausführungsform der Erfindung erläutert;
Fig. 20 ein Blockdiagramm, das ein Gerät darstellt, das ausgestaltet ist, um das funktions-analysierende Verfahren gemäß der siebenten Ausführungsform durchzuführen;
Fig. 21 ein Ablaufdiagramm, das ein Verfahren zum Analysieren der Funktionsweise einer Halbleitervorrichtung gemäß einer achten Ausführungsform der Erfindung erläutert;
Fig. 22 ein Blockdiagramm, das die in dem in Fig. 21 gezeigten Gerät enthaltene Maschenpunkteinheit darstellt;
Fig. 23A und 23B Ablaufdiagramme, die ein Verfahren zum Analysieren der Funktionsweise einer Halbleitervorrichtung gemäß einer neunten Ausführungsform der Erfindung erläuterten;
Fig. 24 ein Blockdiagramm der in einem Gerät zum Durchführen des funktions-analysierenden Verfahrens enthaltenen Maschenpunkteinheit gemäß der neunten Ausführungsform dieser Erfindung;
Fig. 25A eine Querschnittsansicht der von den siebenten bis neunten Verfahren der Erfindung analysierten Halbleitervorrichtungsmodelle;
Fig. 25B ein Diagramm, das die Verunreinigungsverteilungen in den in Fig. 25A gezeigten Vorrichtungsmodellen darstellt; und
Fig. 26A eine grafische Darstellung, die die Ergebnisse der von den siebenten bis neunten erfindungsgemäßen Verfahren durchgeführten Funktionsanalysen zeigt.
Die Erfindung wird nun ausführlich mit Bezug auf die beigefügten Zeichnungen beschrieben.
Fig. 4 ist ein Ablaufdiagramm, das eine Ausführungsform eines ersten Verfahrens um Analysieren der Betriebs einer erfindungsgemäßen Halbleitervorrichtung erläutert. Bei dem Verfahren der ersten Ausführungsform werden, um die Gleichungen (8-1), (8-2) und (8-3) zu lösen, diese Gleichungen durch die folgenden Gleichungen ersetzt, die Zeitdifferentialglieder enthalten:
dp(M,N)/dt = - λp(M,N)fp(M,N) ... (20-1)
dn(M,N)/dt = λn(M,N)fn(M,N) ... (20-2)
dψ(M,N)/dt = λψ(M,N)fψ(M,N) ... (20-3)
wobei fp, fn und fψ wie folgt definiert sind:
fp(M,N) = (1/q) {[Jpx(M) - Jpx(M-1)]/hx'(M)} + (1/q) {[Jpy(N) - Jpy(N-1)]/hy'(M)} + U(M,N) ... (21-1)
fn(M,N) = (1/q) {[Jnx(M) - Jnx(M-1)]/hx'(M)} + (1/q) {[Jny(N) - Jny(N-1)]/hy'(M)} - U(M,N) ... (21-2)
fψ(M,N) = [1/hx'(M)] {[ψ(M+1,N) - ψ(M,N)]/hx(M) - [ψ(M,N) - ψ(M-1,N)]/hx(M-1)} + [1/hy'(N)] {[ψ(M,N+1) - ψ(M,N)]/hy(N) - [ψ(M,N) - ψ(M,N-1)]/hy(N-1)} + (g/&epsi;) [Γ(M,N) + p(M,N) - n(M,N)] ... (21-3)
Die Zeitdifferentialglieder auf den linken Seiten der Gleichungen (20-1), (20-2) und (20-3) werden durch Setzen einer endlichen Anzahl von Maschenpunkten auf der Zeitachse differenzgenähert, wodurch die folgenden Gleichungen erhalten werden:
p¹(M,N) = P&sup0;(M,N) - δtλp(M,N)fp&sup0;(M,N) ... (36-1)
n¹(M,N) = n&sup0;(M,N) + δtλn(M,N)fn&sup0;(M,N) ... (36-2)
ψ¹(M,N) = ψ&sup0;(M,N) + δtλψ(M,N)fφ&sup0;(M,N) ... (36-3)
wobei 1 ≤ M ≤ K, 1 ≤ N ≤ L ist. Bei den Gleichungen (36- 1), (36-2) und (36-3) bezeichnet die hochgestellte 0 einen alten Zeitwert und die hochgestellte 1 einen neuen Zeitwert.
Das Verfahren der ersten Ausführungsform zum Analysieren der Funktionsweise einer Halbleitervorrichtung wird nun mit Bezug auf das Ablaufdiagramm von Fig. 4 erläutert. Zuerst wird bei Schritt a1 die den Aufbau und die Verunreinigungskonzentrationen der Halbleitervorrichtung darstellenden Datenelemente in den Speicher eines Computers eingegeben. Als nächstes werden bei Schritt b1 die Positionen der Maschenpunkte der Vorrichtung bestimmt. Bei Schritt c1 werden die eine Vorspannung V anzeigenden Daten in den Speicher des Computers eingegeben. Der Ablauf geht dann zu Schritt d1, bei dem die anfänglichen Probewerte p&sup0;(M,N), n&sup0;(M,N), ψ&sup0;(M,N) für die grundvariablen Größen p(M,N), n(M,N) und ψ(M,N) angewendet werden. Ferner wird bei Schritt e1 eine Integrationsstartzeit t¹ eingestellt.
Als nächstes werden bei Schritt f1 die Werte auf der rechten Seite der Gleichungen (36-1), (36-2) und (36-3) für alle Maschenpunkte 1 ≤ M ≤ K, 1 ≤ N ≤ L erhalten. Die somit erhaltenen Werte werden als gleich p¹(M,N), n¹(M,N), ψ¹(M,N) betrachtet; wobei 1 ≤ M ≤ K, 1 ≤ N ≤ L. Somit schreitet die Integration auf der Zeitachse einen Schritt vorwärts.
Bei dem nächsten Schritt, d. h. Schritt g1, werden die Absolutwerte der Änderungen der variablen Größen p¹(M,N) - p&sup0; (M,N) , n¹(M,N) - n&sup0;(M,N) , ψ¹(M,N) - ψ&sup0;(M,N) für alle Maschenpunkte erhalten. Ferner wird bei Schritt g1 bestimmt, ob die Absolutwerte dieser Änderungen für alle Maschenpunkte in vorbestimmte Fehlergrenzen fallen oder nicht.
Wenn es auch nur ein NEIN bei einem Maschenpunkt in dem Vorrichtungsraum bei Schritt g1 gibt, geht der Ablauf zu Schritt h1, wobei die korrigierten Werte p¹(M,N), n¹(M,N), ψ¹(M,N) als neue Probewerte verwendet werden, und somit wird t¹ gleich t&sup0; eingestellt und die Integration auf der Zeitachse durchgeführt.
Umgekehrt wird angenommen, wenn bei Schritt g1 bei allen Maschenpunkten JA ist, dass stationäre Lösungen erhalten wurden. Diese Lösungen sind von einer solchen Art, dass die Zeitdifferentialglieder auf den linken Seiten der Gleichungen (20-1), (20-2) und (20-3) Null werden, d. h. fp(M,N) = 0, fn(M,N) = 0, fψ(M,N) = 0; 1 ≤ M ≤ K, 1 ≤ N ≤ L. Dies bedeutet, dass die Lösungen der Gleichungen (8-1), (8-2) und (8-3) erhalten wurden.
Was bis jetzt beschrieben wurde, ist das Grundkonzept des Verfahrens der ersten Ausführungsform der Erfindung. Im Vergleich mit dem herkömmlichen Verfahren ist das Verfahren der ersten Ausführungsform in der folgenden Hinsicht gekennzeichnet.
1. Das Verfahren umfasst keine Prozesse zum Linearisieren von nichtlinearen Gleichungen.
2. Keine massive Matrix-Vektorgleichung, wie beispielsweise Gleichung (16) oder die aus Gleichung (16) hergeleitete Gleichung, die schematisch in Fig. 2 gezeigt ist, muss gelöst werden und folglich müssen keine Berechnungen durchgeführt werden, die so komplex sind, dass große Speicher erforderlich sind. Es ist ausreichend, die Werte der rechten Seiten der Gleichungen (36-1), (36-2) und (36-3) zu finden und die Größen der Grundvariablen zu korrigieren.
3. Da die als Einstellfaktoren verwendeten Koeffizienten λp (M,N), λn(M,N) und λψ(M,N) (nachstehend "Empfindlichkeitskoeffizienten" genannt) unterschiedliche Werte für die Maschenpunkte des Vorrichtungsraums annehmen können, können die Lösungen der Gleichungen wirksam durch Setzen der Empfindlichkeitskoeffizienten auf Werte konvergieren, die für die physikalischen Eigenschaften der Vorrichtung geeignet sind.
4. Der Prozess des Lösens der Gleichungen (20-1) bis (20-3) ist nichts anderes als eine numerische Integration auf der Zeitachse. Die Integration kann natürlich mittels eines großen digitalen Rechners ausgeführt werden, wobei sie jedoch ebenfalls durch ein Gerät mit elektronischen Schaltungen gelöst werden kann, dessen Zeitantwort (d. h. Einschwingphänomen) einen stationären Zustand erreicht, wodurch die spezifischen Gleichungen gelöst werden, wie es später in Verbindung mit der Beschreibung der erfindungsgemäßen Geräte erläutert wird.
Wie es aus dem obigen offensichtlich ist, ist das Verfahren der Ausführungsform dieser Erfindung neuartig und unterscheidet sich beträchtlich von dem herkömmlichen Verfahren, das auf dem Lösen von Matrix-Vektorproblemen aufgebaut ist, und öffnet möglicherweise ein neues Gebiet der Vorrichtungs-Modellierungstechnologie.
Das Gerät, das das Verfahren der ersten Ausführungsform durchführt, wird beschrieben. Das Gerät umfasst, NT logische Zellen bzw. Logikzellen und NT Verstärker, die einen Verstärkungsfaktor λi aufweisen und jeweils mit den logischen Zellen verbunden sind, wobei NT die Gesamtzahl der festgelegten Maschenpunkte ist, wenn die ursprüngliche Kontinuitätsgleichung eines kontinuierlichen Systems in ein diskretes System mittels Differenznäherung oder Finite- Elemente-Näherung transformiert wird, um von einem Computer gelöst zu werden. Die Verstärker sind ausgestaltet, um eine Verstärkung aufzuweisen, die groß genug ist, um tatsächliche Probleme zu bewältigen. Der Verstärkungskoeffizient wird von der Eingabe mittels eines Addierers, eines Multiplizierers und eines Dividierers bestimmt, wie es aus Gleichung (54-12) ersichtlich ist. Das Gerät umfasst ferner eine Wahlvorrichtung (Dial) oder einen Controller zum Einstellen des in Gleichung (54-12) enthaltenen Parameters ω, der dem gesamten Schaltungsnetzwerk zugeordnet ist.
Ein Teil der Gleichung (54-12) zum Berechnen des Verstärkungskoeffizienten kann ein von dem Wert von α oder H abhängiges Näherungsglied sein. Ferner kann der Verstärkungskoeffizient aus der Eingabe in die Schaltung bestimmt werden.
Die Koeffizienten des Verstärkers können abhängig von der Größe der zu lösenden Aufgabenstellung kleiner als 1 sein, und die Verstärker können als Dämpfungsglieder arbeiten.
Ein Multiplizierer, die eine Eingabe mit einen Multiplizierungskoeffizienten λi multiplizieren, kann anstelle der Verstärker verwendet werden.
Das Verfahren der zweiten Ausführungsform zum Analysieren der Funktionsweise einer Halbleitervorrichtung wird nun mit Bezug auf das Ablaufdiagramm von Fig. 5 erläutert.
Wenn die p-Schicht und die n-Schicht, beispielsweise einer pn-Diode, jeweils auf Null-Potential bzw. ein positives Potential eingestellt sind, wird eine Sperrschicht nahe des pn-Übergangs oder der Grenzfläche zwischen der p-Schicht und der n-Schicht ausgebildet. Wenn das Potential V gleich oder kleiner als die Durchbruchspannung der Diode ist, ist der durch die Diode fließende Strom vernachlässigbar klein. In diesem Fall kann die Verteilung des Potentials ψ in der gesamten Vorrichtung durch Lösen nur der Poisson-Gleichung (5) genau bestimmt werden, wie es in der Technik bekannt ist (siehe Kurata, Numerical Analysis for Semiconductor Devices, Heath, 1982, Kapitel 8, und Kurata, Operation Theory of Bipolar Transistors, Kindai Kagaku Co., 1980, Abschnitte 3.2 und 5.1). Die Lochdichte p und die Elektronendichte n werden wie folgt gegeben:
p = nie-θψ ... (37-1)
n = nieθ(ψ-V) ... (37-2)
Durch Einsetzen der Gleichungen (37-1) und (37-2) in die Gleichung (5) wird die folgende Gleichung für einen zweidimensionalen Raum erhalten:
(∂²ψ/∂x²) + (∂²ψ/∂y²) = - (q/&epsi;) [Γ(x,y) + nie-θψ - nieθ(ψ-V)] ... (38)
Die Gleichung (38) wird in ein diskretes System auf die gleiche Art und Weise transformiert, wie es beschrieben wurde. Das diskrete System ist nichts anderes als die Gleichung (9), bei der die Gleichungen (37-1) und (37-2) auf der rechten Seite eingesetzt sind, wie es nachstehend gezeigt ist:
[1/hx'(M)] {[ψ(M+1,N) - ψ(M,N)]/hx(M) - [ψ(M,N) - ψ(M-1,N)]/hx(M-1)} + [1/hy'(N)] {[ψ(M,N+1) - ψ(M,N)]/hy(N) - [ψ(M,N) - ψ(M,N-1)]/hy(N-1)} = - ρ(M,N)/&epsi; ... (39)
wobei
ρ(M,N) = q [Γ(M,N) + nie-θψ(M,N) - nieθ(ψ(M,N)-V)] ... (40)
Wenn das herkömmliche Verfahren benutzt wird, um Gleichung (39) zu lösen, ist es notwendig, eine Matrixgleichung zu lösen und ein Integrationsverfahren durchzuführen, wie es oben beschrieben wurde. Im Gegensatz dazu besteht der erste Schritt des erfindungsgemäßen Verfahren der zweiten Ausführungsform darin, die folgende Gleichung (41) zu lösen, die Zeitdifferentialglieder und einen Empfindlichkeitskoeffizienten enthält:
∂ψ(M,N)/∂t = λ(M,N) · < < [1/hx'(M)]{[ψ(M+1,N) - ψ(M,N)]/hx(M') - [ψ(M,N) - ψ(M-1,N)]/hx(M' - 1)} + [1/hy'(N)]{[ψ(M,N + 1) - ψ(M,N)]/hy(N') - [ψ(M,N) - ψ(M,N-1)]/hy(M'-1)} + ρ(M,N)/&epsi; > > ... (41)
Diese Gleichung wird auf der Zeitachse auf die gleiche Art und Weise wie bei dem Verfahren der ersten Ausführungsform der Erfindung integriert, wie es mit Bezug auf das Ablaufdiagramm von Fig. 5 erläutert wird. Bei dem ersten Schritt a2 werden die Datenelemente, die den Aufbau und die Verunreinigungskonzentrationen der Halbleitervorrichtung darstellen, in den Speicher eines Computers eingegeben. Als nächstes werden bei Schritt b2 die Positionen von Maschenpunkten der Vorrichtung bestimmt. Bei Schritt c2 werden die Daten, die eine Vorspannung V zeigen, in den Speicher des Computers eingegeben. Bei dem nächsten Schritt, Schritt d2, wird der Anfangs-Probewert ψ&sup0; des Potentials in den Speicher eingegeben. Dann wird bei Schritt e2 eine Integrationsstartzeit t¹ eingestellt.
Danach geht der Ablauf zu Schritt f2, bei dem ψ¹ berechnet wird, das wie folgt ausgedrückt wird:
ψ¹(M,N) = ψ&sup0;(M M) + δt·λ(M,N) < < [1/hx'(M)] {[ψ&sup0;(M+1,N) - ψ&sup0;(M,N)]/hx(M') - [ψ&sup0;(M,N) - ψ&sup0;(M-1,N)]/hx(M-1)} + [1/hy'(N)] {[ψ&sup0;(M,N+1) - ψ&sup0;(M,N)]/hy(N') - [ψ&sup0;(M,N) - ψ&sup0;(M,N-1)]/hy(M'-1)} + ρ(M,N)/&epsi; > > ... (42)
Der Wert der rechten Seite der Gleichung (42) wird für alle Maschenpunkte 1 ≤ M ≤ K, 1 ≤ N ≤ L gefunden. Die somit erhaltenen Werte werden als gleich zu ψ¹(M,N) betrachtet; 1 ≤ M ≤ K, 1 ≤ N ≤ L. Somit schreitet die Integration auf der Zeitachse einen Schritt vorwärts.
Bei dem nächsten Schritt, d. h. Schritt g2, werden die Absolutwerte der Änderungen in der variablen Größe ψ¹(M,N) - ψ&sup0;(M,N) für alle Maschenpunkte erhalten. Ferner wird bei Schritt g2 bestimmt, ob die Absolutwerte dieser Änderungen für alle Maschenpunkte in vorbestimmte Fehlergrenzen fallen oder nicht.
Wenn auch nur ein NEIN bei einem Maschenpunkt bei Schritt g1 gibt, geht der Ablauf zu Schritt h2, bei dem die korrigierten Werte ψ¹(M,N) als neue Probewerte verwendet werden, und somit als gleich ψ&sup0;(M,N) angesehen werden. Dann wird die Zeit t¹ auf t&sup0; eingestellt, und die Integration wird an der Zeitachse auf die gleiche Art und Weise, wie bei dem Verfahren der ersten Ausführungsform durchgeführt.
Umgekehrt wird angenommen, wenn bei Schritt g2 bei allen Maschenpunkten JA ist, dass stationäre Lösungen erhalten wurden. Diese Lösungen sind von einer solchen Art, dass die Zeitdifferentialglieder auf der linken Seite von Gleichung (20-3) Null sind, d. h. fψ(M,N) = 0; 1 ≤ M ≤ K, 1 ≤ N ≤ L. Dies bedeutet, dass die Lösung der Gleichung (9) oder (39) erhalten wurde.
Wenn die Absolutwerte der Änderungen in der unbekannten variablen Größe ψ(M,N), d. h. ψ¹(M,N) - ψ&sup0;(M,N) , zu klein sind, dann kann die Integration auch nach der Zeit t keinen wahren Wert von ψ(M,N) in einer kurzer Zeit ergeben, obgleich der Integrationsprozess stationär ist. Wenn andererseits die Absolutwerte der Änderungen zu groß sind, ist der Integrationsprozess instabil und ψ(M,N) fluktuiert oder divergiert, wobei es misslingt, den wahren Wert von ψ(M,N) zu erhalten.
Der Empfindlichkeitskoeffizient λ(M,N) ist ein Einstellfaktor, um die variable Größe zu veranlassen, sich geeignet zu ändern. Die folgende Gleichung kann als eine Gleichung angewendet werden, die den Empfindlichkeitskoeffizienten λ(M,N) festlegt, um die Poisson-Gleichung (42) zu lösen, die ein zweidimensionales Vorrichtungsmodell darstellt. Zuerst werden verschiedene Größen für die Gleichung (42) wie folgt definiert:
H = 1/[hx(M'-1)hx'(M)] + 1/[hx(M')hx'(M)] + 1/[hy(N'-1)hy'(N)] + 1/[hy(N')hy'(N)] ... (43-1)
α = (qθni/&epsi;) [exp{-θψ&sup0;(M,N)} + exp{θ(ψ&sup0;(M,N)-V)}] ... (43-2)
Daher:
λ(M,N)·δt ≤ 2/(2H+α) ... (43-3)
Die Bedeutung der Gleichung (43-3) wird qualitativ beschrieben. Im Allgemeinen ist entweder die Elektronendichte n oder die Lochdichte p nahezu gleich der Verunreinigungskonzentration in dem Hochverunreinigungsbereich einer Halbleitervorrichtung, und der neutrale Zustand der Raumladung ist in diesem Bereich erfüllt. In diesem Fall ist, da α weitaus größer als H und der Zähler 2 ist, λ(M,N)δt extrem klein. Folglich sind die Änderungen ψ¹(M,N) - ψ&sup0;(M,N) ebenfalls extrem klein. In einem Bereich, der eine relativ niedrige Verunreinigungskonzentration aufweist, wo es wahrscheinlicher ist, dass eine Sperrschicht gebildet wird, ist α nahezu gleich H oder vernachlässigbar klein verglichen mit H und λ(M,N)δt ist relativ groß. In diesem Fall werden die Änderungen ψ¹(M,N) - ψ&sup0; (M,N) relativ groß sein.
Aus den oben angegebenen Gründen ändert sich die variable Änderung ψ¹(M,N) - ψ&sup0;(M,N) in dem Niedrig- Verunreinigungsbereich stark, wobei ein Wert geliefert wird, der den neutralen Zustand der Raumladung erfüllt, für die variable Größe ψ(M,N) gegeben wird. Als ein Ergebnis wird die Bildung einer Sperrschicht gefördert, und die variable Größe ψ(M,N) konvergiert zu der wahren Lösung.
Das Gerät, das das Verfahren der zweiten Ausführungsform durchführt, wird mit Bezug auf Fig. 6 beschrieben. Wie es in Fig. 6 gezeigt ist, umfasst das Gerät NT logische Zellen und NT Verstärker ai, die einen Verstärkungskoeffizienten λi aufweisen bzw. mit den logischen Zellen verbunden sind, wobei NT die Gesamtzahl der definierten Maschenpunkte ist, wenn die ursprüngliche Kontinuitätsgleichung in eine diskrete Gleichung mittels einer Differenznäherung oder Finite- Elemente-Näherung transformiert wird, um durch einen Computer gelöst zu werden. Die Verstärker sind ausgestaltet, damit die Verstärkungskoeffizienten λi groß genug sind, um tatsächliche Probleme und Variablen zu bewältigen. Die Verstärker umfassen ferner eine Wahlvorrichtung oder einen Controller zum Variieren des Verstärkungskoeffizienten jedes Verstärkers. Die Koeffizienten der Verstärker können abhängig von der Größe des zu lösenden Problems kleiner als 1 sein. In diesem Fall können die Verstärker als Dämpfungsglieder arbeiten. Die Verstärker können ebenfalls durch Multiplizierer, die einem Verstärker äquivalent sind, jeweils zum Erhalten des Produkts einer Eingabe und eines Verstärkungskoeffizienten ersetzt werden.
Genauer gesagt umfasst das in Fig. 6 gezeigte Gerät zwölf logische Zellen 1 bis 12 und zwölf Verstärker a1 bis a12, d. h. NT = 12. Die logischen Zellen 1 bis 12 weisen eine nichtlineare Eingangs-/Ausgangs-Antwort auf. Jeder Verstärker ist mit dem Eingangsanschluss der zugeordneten logischen Zelle verbunden. Das Gerät umfasst ferner Resonanzschaltungen b1 bis b12, die jeweils aus einem Widerstand und einem Kondensator aufgebaut sind und mit dem Anschlußpunkt der zugeordneten logischen Zelle und des Verstärkers verbunden sind. Die variable Verstärkung jedes Verstärkers entspricht λ.
Ströme werden von den Stromquellen c1 bis c12 an die logischen Zellen 1 bzw. 12 geliefert. Stromerfassungsmittel 11, 21, 22, 32, ... 1112 und 1212, die für die logischen Zellen vorgesehen sind und ein Netzwerk bilden, erfassen, dass die Eingangsströme in dem Gerät stationär werden. Eine Adresseneinheit 100 ist zum Kennzeichnen der Stromerfassungsmittel 11, 21, ... 1212 vorgesehen. Der von einem beliebigen, durch eine Adresseneinheit 100 gekennzeichneten Stromerfassungsmittel erfasste Strom wird in ein digitales Signal mittels eines A/D-Wandlers 101 umgewandelt. Das digitale Signal wird an eine Schaltung 103 geliefert, die ein Differenzierglied, einen Multiplizierer und einen Addierer umfasst. Das Differenzierglied, der Multiplizierer und der Addierer führen Operationen an dem Signal durch, wobei Messdaten geliefert werden. Die Messdaten werden in eine Schaltung 104 eingegeben.
Das in Fig. 6 dargestellte Gerät kann einen Sweeper aufweisen, wobei in diesem Fall die Anzahl erforderlicher logischer Zellen auf m verringert wird, wobei m < NT ist.
Das erfindungsgemäße Verfahren der dritten Ausführungsform, d. h. ein Verfahren des kleinsten Quadrats (least square method) zum Bestimmen der Äquivalent- Schaltungskonstanten einer Halbleitervorrichtung wird mit Bezug auf Fig. 7A, 7B und 8 beschrieben.
Es sei angenommen, dass vier y-Parameter y11, y12, y21 und y22 - alle komplexen weggelassen - bei N Frequenzen fk (k = 1, 2, .., N) gemessen werden. Um das Anpassen der vier y- Parameter mittels der in Fig. 7A und 7B gezeigten Äquivalenzschaltung zu erreichen, die M äquivalente Schaltungskonstanten q&sub1;, q&sub2;, ..., qM aufweist, ist es notwendig, die folgenden Matrixgleichungen zu lösen:
Aq = b, A = (aij), b = (bi) ... (44-1)
aij = [α*lmkiαlmkj + αlmkiα*lmkj]/ ηlmk ² ... (44-2)
bi = - [β*lmkαlmki + βlmkα*lmki]/ ηlmk ² ... (44-3)
wobei
αlmki = ∂ylmk/∂qi ... (45-1)
βlmk = ylmk&sup0; - nlmk - Σαlmkqi&sup0; ... (45-2)
In den Gleichungen (44-2) und (44-3) gibt das Zeichen * eine konjugiert komplexe Zahl an, und 1 und m sind entweder 1 oder 2.
Die Matrixgleichungen werden in dq/dt = - λ(Aq - b) umgeschrieben, was gelöst wird, wie es mit Bezug auf das Ablaufdiagramm von Fig. 8 erläutert wird. Zuerst werden bei Schritt a3 die gemessenen y-Parameter η11k, η12k, η21k, η22k (k = 1, 2, ... N) gelesen. Als nächstes werden bei Schritt b3 anfängliche Probewerte der unbekannten Größe q&sub1;&sup0;, q&sub2;&sup0;, ..., qM&sup0; gegeben. Bei Schritt c3 werden eine Anfangszeit und ein Zeitintervall eingestellt, und die Zeitintegration wird gestartet. Dann werden bei Schritt d3 αlmki und βlmk durch Anwenden des theoretischen Werts ylmk und des Messwerts ηlmk des y-Parameters erhalten. Ferner werden bei Schritt e3 aij und bij erhalten.
Der Ablauf geht zu Schritt f3, bei dem qi¹, das durch die folgende Gleichung (46) ausgedrückt wird, für i = 1, 2, ..., M berechnet wird:
qi¹ = qj&sup0; - δt·λi(Ag&sub0;-b) ... (46)
Folglich schreitet die Integration auf der Zeitachse um einen Schritt vorwärts.
Danach werden bei Schritt g3 die Absolutwerte der Änderungen der variablen Größe qi¹ - qi&sup0; für i = 1, 2, ..., M erhalten. Ferner wird bei Schritt g3 bestimmt, ob die Absolutwerte dieser Änderungen für alle Konstanten der Äquivalent-Schaltung q1, q2, ..., qM in vorbestimmte Fehlergrenzen fallen.
Bei NEIN bei Schritt g3, geht der Ablauf zu Schritt h3, bei dem die korrigierten Werte qi¹ als neue Probewerte verwendet werden, und somit als gleich qi&sup0; betrachtet werden. Dann wird die Zeit t¹ auf t&sup0; eingestellt, und die Integration auf der Zeitachse wird auf die gleiche Art und Weise wie bei Verfahren der ersten und zweiten Ausführungsformen durchgeführt.
Umgekehrt, bei JA bei Schritt g3, wird angenommen, dass die Konstanten der Äquivalent-Schaltung q&sub1;, q&sub2;, ..., qM erhalten wurden. Diese Konstanten sind von einer solchen Art, dass das Zeitdifferentialglied auf der linken Seite der Gleichung dq/dt = - λ(Aq - b) Null wird. Somit wurde die Lösung der Gleichung Aq = b erhalten.
Eine numerische Analyse, die tatsächlich durchgeführt wurde, um eine Aufgabenstellung eines Bestimmens der Konstanten der äquivalenten Schaltung eines Transistors zu lösen, wird nun mit Bezug auf Fig. 7A und 7B erläutert. Fig. 7A ist eine schematische Darstellung eines bipolaren Transistors. Es sei angenommen, dass die Emitter-geerdete Äquivalent-Schaltung des bipolaren Transistors aus drei Netzwerken mit vier Anschlüssen aufgebaut ist, wie es in Fig. 7B dargestellt ist. Dann können die y-Parameter der Äquivalent-Schaltung ohne weiteres bestimmt werden. Genauer gesagt sind die somit erhaltenen y-Parameter yi des wesentlichen Teils 1000 der Schaltung:
y&sub1;&sub1;i = {gb [gegee(1-α) - ω²cecci] + jωgb(cegee+ccige+ccigee)}/Δyi ... (47-1)
y&sub1;&sub2;i = - [jωccigb(ge+gee+jωce)]/Δyi ... (47-2)
y&sub2;&sub1;i = [gb(αgegee+ω²cecci) - jωccigb(ge+gee)]/Δyi ... (47-3)
y&sub2;&sub2;i = [-ω2cecci(ge+gee) + jωcci(gbge+gegee+gee+geegb)}/Δyi ... (47-4)
Δyi = gbge + gbgee + gegee(1-α) - ω2cecci + jω [ce(gee+gb) + cci(ge+gee)] ...(47-5)
In Übereinstimmung mit der Regel zum Synthetisieren von Netzwerken mit vier Anschlüssen, werden die y-Parameter y* des Abschnitts 100 und die externe Kollektorkapazität 1001, die zusammen kombiniert sind, wie folgt definiert:
Y&sub1;&sub1;* = y&sub1;&sub1;i + jωcco ... (48-1)
Y&sub1;&sub2;* = y&sub1;&sub2;i - jωcco ... (48-2)
Y&sub2;&sub1;* = y&sub2;&sub1;i - jωcco ... (48-3)
y&sub2;&sub2;* = y&sub2;&sub2;i + jωcco ... (48-4)
Ferner werden in Übereinstimmung mit der Regel des Synthetisierens von Netzwerken mit vier Anschlüssen die y- Parameter des gesamten bipolaren Transistors berechnet, d. h. Abschnitte 1000, 1001 und 1002, die die Kollektor-Konduktanz gc umfassen, was die folgenden Ergebnisse ergibt:
y&sub1;&sub1; = (gcy&sub1;&sub1;*+Δy*)/(gc+y&sub2;&sub2;*) ... (49-1)
y&sub1;&sub2; = (gcy&sub1;&sub2;*)/(gc+y&sub2;&sub2;*) ... (49-2)
y&sub2;&sub1; = (gcy&sub2;&sub1;*)/(gc+y&sub2;&sub2;*) ... (49-3)
y&sub2;&sub2; = (gcy&sub2;&sub2;*)/(gc+y&sub2;&sub2;*) ... (49-4)
Δy* = y&sub1;&sub1;*y&sub2;&sub2;* - y&sub1;&sub2;*y&sub2;&sub1;* ... (49-5)
Es gibt acht Konstanten der in Fig. 7A und 7B dargestellten Schaltung, die bestimmt werden müssen. Diese sind: ge, gee, gb, gc, α, ce, cci und cco. Dle Objektfunktion, die die Kennlinie der y-Parameteranpassung ist, wird durch die folgende Gleichung dargestellt:
E = ylmk - ηlmk ²/ηlmk ² ... (50)
Dies stellt einen relativen Fehler der Anpassung zwischen dem theoretischen Wert ylmk des y-Parameters und dem Messwert ηlmk des y-Parameters dar. Hier sind 1 und n entweder 1 oder 2, und k ist die Anzahl der verwendeten Messfrequenzen, die 1, 2, ... oder N ist.
Um die tatsächliche Berechnung zu ermöglichen, werden acht unbekannte Größen p&sub1; bis p&sub8; und acht Normierungseinheiten r&sub1; bis r&sub8; für diese Größen wie folgt definiert. Es sei bemerkt, dass qi = ripi, i = 1, 2, ..., 8 durch Definition der Normierung ist.
p&sub1; = ge(1-α), r&sub1; = 1
p&sub2; = age, r&sub2; = 1
p&sub3; = gee, r&sub3; = 1
p&sub4; = gb, r&sub4; = 1
p&sub5; = gc, r&sub5; = 1
p&sub6; = ce, r&sub6; = 10&supmin;¹²
p&sub7; = cci, r&sub7; = 10&supmin;¹²
p&sub8; = cco, r&sub8; = 10&supmin;¹²
Um die Anpassung der y-Parameter zu erreichen, müssen Messdaten erfasst werden. Idealerweise sollten die Messdaten einer tatsächlichen Vorrichtung (d. h. eines bipolaren Transistors) für diesen Zweck verwendet werden. In Übereinstimmung mit einer Anforderung, dass die Daten mit bekannter Natur verwendet werden sollten, werden jedoch geeignete Werte (nachstehend als "wahre Werte" bezeichnet) den Schaltungskonstanten p&sub1; bis p&sub8; zugeordnet, um y-Parameter zu erzeugen, wodurch künstliche "Messdaten" geliefert werden. Die den Schaltungskonstanten zugeordneten wahren Werte sind:
ge = 0.6siemens p&sub1; = ge(1-α) = 0.012siemens
α = 0.98 p&sub2; = gea = 0.588siemens
gee = 0.08siemens = p&sub3;
gb = 0.04siemens = p&sub4;
gc = 0.15siemens = p&sub5;
ce = 200fF = p&sub6;
cci = 3fF = p&sub7;
cco = 50fF = p&sub8;
Vierzig "gemessene" Frequenzen werden bereitgestellt, die von 0,1 GHz bis 19,6 GHz reichen, die sich jeweils um 0,5 GHz von der direkt benachbarten unterscheiden. Bevor das Anpassen der y-Parameter gestartet wird, werden die folgenden anfänglichen Probewerte an die Schaltungskonstanten gegeben:
ge = 0.7siemens p&sub1; = 0.014siemens
α = 0.98 p&sub2; = 0.686siemens
gee = 0.1siemens = p&sub3;
gb = 0.2siemens = p&sub4;
gc = 0.2siemens = p&sub5;
ce = 300fF = p&sub6;
cci = 5fF = p&sub7;
cco = 60fF = p&sub8;
Vier y-Parameter, denen die wahren Werte und die Probewerte, alle oben spezifiziert, gegeben wurden, werden in einer komplexen Ebene aufgetragen. Die y-Parameter sind in der komplexen Ebene verteilt, wie es in den Fig. 9A bis 9D dargestellt ist. Da die Messdaten künstlich erfasst wurden, wie es beschrieben wurde, enthalten sie keine Rauschbestandteile, die den Messdaten der tatsächlichen bipolaren Transistoren inhärent sind. Nichtsdestoweniger, wie es die Fig. 9A bis 9D klar zeigen, unterscheiden sich der wahre Wert und der Probewert für die gleiche Frequenz beträchtlich voneinander. Um die Differenz zwischen dem wahren Wert und dem Probewert quantitativ zu analysieren, wird die Quadratwurzel der durch Gleichung (50) definierten Objektfunktion bestimmt, wodurch der relative Fehlervektor für alle y-Parameter ausgewertet wird. Das heißt:
R. E. = [ ylmk - ηlmk ²/ ηlmk ²]1/2 ... (51)
Bei dem in Fig. 9A bis 9D gezeigten Fall ist der relative Fehler R. E. = 3,19 oder 319%. Dieser Wert ist zu groß, da die Bedingung für eine wünschenswerte Anpassung der y-Parameter darin besteht, dass der relative Fehlervektor für alle y-Parameter etwa 1% oder weniger sein sollte.
Die Intervalle, mit denen die Integration auf der Zeitachse durchgeführt wird, werden auf 10&supmin;&sup9; sec eingestellt, d. h. 1 Nanosekunde. Die Ergebnisse der tatsächlich mittels eines Computers durchgeführten fünf Tests waren, wie es in Tabelle 1 gezeigt ist. Der erste Job S bestand aus 500 Zeitschritten einer Berechnung, wobei die in der zweitäußersten linken Spalte von Tabelle 1 spezifizierten λ&sub1; bis λ&sub8; verwendet wurden. Der bei dem ersten Job S erhaltene relative Fehler R. E. ist 4,12%. Dieser Wert von 4,12% ist tatsächlich weitaus kleiner als 319%, wobei er jedoch nicht als eine ausreichende Bedingung für eine wünschenswerte Anpassung der y-Parameter angesehen werden kann. Folglich wird der zweite Job W ausgeführt, wobei in Übereinstimmung mit den Zwischenergebnissen des ersten Jobs S p&sub6; (= ce) stark geändert wird, während p&sub7; wenig geändert wird. Die bei dem zweiten Job W erhaltenen Ergebnisse sind: R. E. = 2,47%. Offensichtlich wurde der relative Fehler verringert. Dann wird der dritte Job I ausgeführt, wobei die Ergebnisse des Jobs W als Anfangsbedingung verwendet, und 500 Berechnungszeitschritte durchgeführt werden, wobei λ&sub1; bis λ&sub8; identisch zu denjenigen verwendet werden, die bei dem zweiten Job W verwendet wurden. Der dritte Job I ergibt den relativen Fehler R. E. von 1,07%. Ferner wird der vierte Job Q auf eine ähnliche Art und Weise durchgeführt, wobei λ&sub1; bis λ&sub8; verwendet werden, die mit denjenigen bei dem zweiten Job W verwendeten mit der Ausnahme identisch sind, dass λ&sub7; = 5.E3 ist (der vorherige Wert: 5.E2). Der vierte Job Q ergibt den relativen Fehler R. E. von 0,210%, was weitaus geringer als die obere Grenze für eine wünschenswerte Anpassung der y- Parameter ist.
Diese Art von Integration umfasst den folgenden allgemeinen Trend. Das heißt, wenn ein relativ gutes Ergebnis, z. B. R. E. von etwa einigen Prozent, in einem Job erhalten wird, wird der nächste Job mit dem Ergebnis des vorherigen Jobs als Anfangswert ausgeführt, wobei λ identisch mit oder geringfügig größer als dasjenige bei dem vorherigen Job verwendet wird. Dann ergibt der nächste Job bessere Ergebnisse.
Der letzte Job U wird mit den Ergebnissen des vierten Jobs Q als Anfangswert durchgeführt. Wenn 192 Zeitschritte beendet sind, erreicht der relative Fehler R. E. einen Minimalwert von 0,1% oder den bestmöglichen Wert. Die somit festgelegten y-Parameter werden in der komplexen Ebene aufgetragen. Die Ergebnisse dieser Auftragung sind wie es in den Fig. 9A bis 9D gezeigt ist. Soweit es die grafischen Darstellungen von Fig. 9A bis 9D betrifft, erscheint es, dass die resultierenden y-Parameter vollständig mit denjenigen koinzidieren, die die wahren Werte von p&sub1; bis p&sub8; geben. Somit kann die Anpassung der y-Parameter, die durch die numerische Analyse erreicht wurde, als ausreichend genau angesehen werden.
Fig. 10 stellt ein Gerät dar, um die oben beschriebene numerische Analyse durchzuführen, die durchgeführt wurde, um die Konstanten der Äquivalent-Schaltung des bipolaren Transistors zu bestimmen. Wie Fig. 10 klar zeigt, umfasst dieses Gerät acht logische Zellen 1 bis 8, acht Verstärker a1 bis aß, die jeweils mit den Eingängen der logischen Zellen 1 bis 8 verbunden sind, acht Resonanzkreise, die jeweils aus einem Widerstand (R1, R2, ..., oder R8) und einem Kondensator (C1, C2, ..., oder C8) aufgebaut sind und mit den Eingängen der logischen Zellen 1 bis 8 verbunden sind. Die logischen Zellen 1 bis 8 entsprechen den acht Schaltungskonstanten p&sub1; bis p&sub8;, die jeweils eine nichtlineare Eingangs-/Ausgangs- Antwort V = g(v) aufweisen. Die Verstärker a1 bis a8 weisen variable Verstärkungen auf, die λ&sub1;, λ&sub2;, ..., λ&sub8; entsprechen. Wie es hervorgehoben wurde, ist es bedeutsam, die Verstärkungskonstanten der Verstärker a1 bis aß für eine erfolgreiche numerische Analyse der Schaltungskonstanten der bipolaren Transistoren einzustellen. Tabelle 1
Die Kondensatoren C1, C2, ..., C8 geben Ströme C (dv/dt), die den Zeitdifferentialgliedern entsprechen, wohingegen die Widerstände R1, R2, ..., R8 eine Relaxationszeitkonstante τ = CR definieren. Stromquellen S1, S2, ..., S8 sind jeweils mit den Eingängen der Verstärker a1. a2, ..., a8 verbunden. Diese Stromquellen S1 bis S8 entsprechen den konstanten Gliedern bi, die in der Gleichung (43-3) enthalten sind, wie es später beschrieben wird.
Das in Fig. 10 gezeigte Gerät umfasst ferner Stromquellen 11, 21, ..., 81, ..., 18, 28, ... 88. Jede dieser acht Stromquellen ist jeweils mit den Ausgängen der logischen Zellen 1 bis 8 verbunden. Die von diesen Quellen gelieferten Ströme werden als TijVj ausgedrückt, was bedeutet, dass der j-te Parameter pj mit einer Normierungseinheit rj normiert ist, wodurch eine Größe qj erhalten wird, die sich auf die Ausgangsspannung fj der j-ten logischen Zelle bezieht. Das in Fig. 10 gezeigte Gerät kann daher die folgenden Gleichungen (52) und (53) liefern, um die Einschwingantwort des in Fig. 10 gezeigten Schaltungsnetzwerks zu bestimmen:
Ci(dvi/dt) = - vi/Ri + λi ( Tijvj + Ii) ... (52)
Vi = g(vi) ... (53)
wobei i die Nummer einer beliebigen logischen Zelle, vi und Vi die Eingangs- und Ausgangsspannungen der i-ten logischen Zelle, Ri, Ci und Ii der Widerstand, der Kondensator bzw. die Stromquelle, Tij die Kopplungsstärke zwischen dem Eingang der i-ten logischen Zelle und dem Ausgang der j-ten logischen Zelle, und g die Eingangs-/Ausgangs-Antwortfunktion jeder logischen Zelle ist. Beide Gleichungen (52) und (53) reduzieren sich auf die Gleichung dq/dt = - λ(Aq - b), vorausgesetzt, dass die folgenden Bedingungen erfüllt sind.
(1) Tij/Ci entspricht aij.
(2) Ii/Ci entspricht der Konstanten bi.
(3) Die Funktion g wird durch die lineare Funktion Vi = vi ersetzt.
(4) Die Zeitkonstante CiRi, die zum Lösen von dq/dt = - λ(Aq - b) nicht benötigt wird, wird durch Ri = ∞ ersetzt.
Mit anderen Worten wird die Gleichung dq/dt = - λ(Aq - b) automatisch mit dem in Fig. 10 gezeigten Schaltungsnetzwerk basierend auf den obigen Bedingungen (1) bis (4) gelöst. Es ist jedoch notwendig, um die Stromquellen TijVj und 11, die von Spannungen abhängen, mittels analoger Techniken vorzusehen, die im rechten Abschnitt von Fig. 10 dargestellte Steuervorrichtung zu benutzen. Die Steuervorrichtung umfasst eine Adresseneinheit 100, einen D/A-Wandler 101 und eine Schaltung 103 mit einem Differenzierer, einem Multiplizierer und einem Addierer.
Das in Fig. 10 gezeigte Gerät benötigt von der Vorrichtung (nicht gezeigt), die die y-Parameter gemessen hat, gelieferte Messdaten 104. Die Messdaten 104 werden als ein Bezug zum Erreichen der Anpassung der y-Parameter verwendet. Die Messdaten 104 und die Ausgangsspannungen der logischen Zellen 1 bis 8 werden von der Schaltung 103 differenziert, multipliziert und addiert, wodurch Gleichungen (44-1), (44-2) und (44-3) gelöst werden. Die Daten und die Spannungen können entweder in der Form von digitalen Daten oder in der Form von analogen Daten verarbeitet werden. Mit Blick auf die Genauigkeit der Berechnung erscheint es jedoch besser, sie in der Form von analogen Daten zu verarbeiten.
In dem Fall, bei dem die die Glieder der Gleichungen (44-1), (44-2) und (44-3) darstellenden Daten digital sind, werden sie in analoge Signale mittels des D/A-Wandlers umgewandelt, und die so bereitgestellten analogen Signale werden an die Adressen (i, j) übertragen, wobei i und j = 1, 2, ..., 8 ist, die durch die Adresseneinheit 100 gekennzeichnet wurden.
Nun wird ein Verfahren zum Bestimmen der Empfindlichkeitskoeffizienten beschrieben, was das erfindungsgemäße Verfahren der vierten Ausführungsform ist und ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung zu analysieren. Ähnlich wie das oben beschriebene Verfahren der zweiten Ausführungsform kann dieses Verfahren λ bestimmen, um dadurch die Poisson- Gleichung (5) der Halbleitervorrichtung zu lösen.
Beim Modellieren der Vorrichtung wird, wenn beispielsweise die p-Schicht und die n-Schicht einer pn-Diode jeweils auf Null-Potential bzw. ein positives Potential gesetzt sind, eine Sperrschicht nahe dem p-n-Übergang oder der Grenzfläche zwischen der p-Schicht und der n-Schicht gebildet. Wenn das Potential V gleich oder kleiner als die Durchbruchspannung der Diode ist, ist der durch die Diode fließende Strom vernachlässigbar klein. In diesem Fall kann die Verteilung des Potentials der gesamten Vorrichtung genau bestimmt werden, in dem nur die Poisson-Gleichung (5) gelöst wird, d. h. eine der fünf Grundgleichungen (1) bis (5), wie es allgemein in der Technik bekannt ist (siehe Kurata, Numerical Analysis for Semiconductor Devices, Heath, 1982, Kapital 8, und Kurata, Operation Theory of Bipolar Transistors, Kindai Kaguku Co., 1980, Abschnitte 3.2 und 5.1). Die Lochdichte p und die Elektronendichte n werden wie folgt gegeben:
P = nie-θψ ... (37-1)
n = nieθ(ψ-V) ... (37-2)
Durch Einsetzen der Gleichungen (37-1) und (37-2) in Gleichung (5) wird die folgende Gleichung für einen zweidimensionalen Raum erhalten:
(∂²ψ/∂x²) + (∂²ψ/∂y²) = - (q/&epsi;) [Γ(x,y) + nie-θψ - nieθ(ψ-V)] ... (38)
Gleichung (38) wird in ein diskretes System auf die gleiche Art und Weise transformiert, wie es beschrieben wurde. Das diskrete System ist nichts anderes als Gleichung (9), wobei die Gleichungen (37-1) und (37-2) auf der rechten Seite eingesetzt werden, wie es nachstehend gezeigt ist:
[1/hx'(M)] {[ψ(M+1,N) - ψ(M,N)/hx(M) - [ψ(M,N) - ψ(M-1,N)]/hx(M-1)} + [1/hy'(M)] {[ψ(M,N+1) - ψ(M,N)]/hy(N) - [ψ(M,N) - ψ(M,N-1)]/hy(N-1)} = - ρ(M,N)/&epsi; ... (39)
wobei
ρ(M,N) = q[Γ(M,N) + nie - θψ(M,N) - nieθ(ψ(M,N)-V)] ... (40)
Wenn das herkömmliche Verfahren verwendet wird, um Gleichung (39) zu lösen, wird es notwendig sein, eine Matrix- Vektorgleichung zu lösen und ein Iterationsverfahren durchzuführen, wie es oben beschrieben wurde. Im Gegensatz dazu besteht der erste Schritt des erfindungsgemäßen Verfahrens der vierten Ausführungsform darin, die folgende Gleichung (41) zu lösen, die Zeitdifferentialglieder und einen Empfindlichkeitskoeffizienten enthält:
∂ψ(M,N)/∂t = λ(M,N) · < < [1/hx'(M)] · {[ψ(M+1,N) - ψ(M,N)]/hx(M') - [ψ(M,N) - ψ(M-1,N)/hx(M'-1)} + [1/hy'(N)] · {[ψ(M,N+1) - ψ(M,N)/hy(N') - [ψ(M,N) - ψ(M,N-1)/hy(M'-1)} + ρ(M,N)/&epsi;> > ... (41)
Diese Gleichung wird auf der Zeitachse integriert, wie es mit Bezug auf das Ablaufdiagramm von Fig. 11 erläutert wird. Bei dem ersten Schritt a4 werden die Datenelemente, die die Struktur und Verunreinigungskonzentrationen der Halbleitervorrichtung darstellen, in den Speicher eines Computers eingegeben. Als nächstes werden bei Schritt b4 die Positionen der Maschenpunkte der Vorrichtung bestimmt. Bei Schritt c4 werden die Daten, die eine Vorspannung V zeigen, in den Speicher des Computers eingegeben. Bei dem nächsten Schritt, Schritt d4, wird der Anfangsprobewert ψ&sub0; des Potentials in den Speicher eingegeben. Dann wird bei Schritt e4 die Integrationsstartzeit t¹ eingestellt.
Danach geht der Ablauf zu Schritt f4, bei dem ψ¹ berechnet wird, was wie folgt ausgedrückt wird:
ψ¹(M,N) = ψ&sup0;(M,N) + δt·λ(M,N) · < < [1/hx'(M)] · {[ψ&sup0;(M+1,N) - ψ&sup0;(M,N)]/hx(M') - [ψ&sup0;(M,N) - ψ&sup0;(M-1,N)]/hx(M'-1)} + [1/hy'(N)] · {[ψ&sup0;(M,N+1) - ψ&sup0;(M,N)]/hy(N') - [ψ&sup0;(M,N) - ψ&sup0;(M,N-1)]/hy(M'-1)} + ρ(M,N)/e> > ... (42)
Der Wert der rechten Seite der Gleichung (42) wird für alle Maschenpunkte 1 ≤ M ≤ K, 1 ≤ N ≤ L gefunden. Die so erhaltenen Werte werden als gleich ψ¹(M,N); 1 ≤ M ≤ K, 1 ≤ N ≤ L angesehen. Somit rückt die Integration einen Schritt auf der Zeitachse vorwärts.
Bei dem nächsten Schritt, d. h. Schritt g4, werden die Absolutwerte der Änderungen in der variablen Größe ψ¹(M,N) - ψ&sup0; (M,N) für alle Maschenpunkte erhalten. Ferner wird bei Schritt g4 bestimmt, ob die Absolutwerte dieser Änderungen für alle Maschenpunkte in vorbestimmte Fehlergrenzen fallen oder nicht.
Wenn es auch nur ein NEIN bei einem Maschenpunkt bei Schritt g4 gibt, geht der Ablauf zu Schritt h4, bei dem die korrigierten Werte ψ¹ (M,N) als neue Probewerte verwendet und folglich als ψ&sup0;(M,N) betrachtet werden. Dann wird die Zeit t¹ auf t&sup0; eingestellt, und die Integration wird auf der Zeitachse auf die gleiche Art und Weise wie bei dem Verfahren der ersten Ausführungsform durchgeführt.
Umgekehrt wird angenommen, wenn bei Schritt g4 bei allen Maschenpunkten JA ist, wird angenommen, dass stationäre Lösungen erhalten wurden. Diese Lösungen sind von einer solchen Art, dass die Zeitdifferentialglieder auf der linken Seite der Gleichung (20-3) Null werden, d. h. fψ(M,N) = 0; 1 ≤ M ≤ K, 1 ≤ N ≤ L. Dies bedeutet, dass die Lösung der Gleichung (39) erhalten wurde.
Wenn die Absolutwerte der Änderungen in der unbekannten variablen Größe ψ(M,N), d. h. ψ¹(M,N) - ψ&sup0;(M,N) zu klein sind, kann auch nach der Zeit t die Integration keinen wahren Wert von ψ(M,N) innerhalb einer kurzen Zeit liefern, obgleich der Integrationsprozess stationär ist. Wenn andererseits die Absolutwerte der Änderungen zu groß sind, ist der Integrationsprozess instabil, und ψ(M,N) fluktuiert oder divergiert, wobei es misslingt, den wahren Wert von ψ(M,N) zu erhalten.
Der Empfindlichkeitskoeffizient λ(M,N) ist ein Einstellfaktor zum Veranlassen, dass sich die variable Größe geeignet verändert. Der die obige Bedingung erfüllende Empfindlichkeitskoeffizient λ(M,N) wird wie folgt bestimmt, um die Poisson-Gleichung (42) zu lösen, die ein zweidimensionales Halbleitervorrichtungsmodell darstellt. Zuerst wird die Variable ψ in Gleichung (41) durch ψ&sup0; + δψ&sup0; ersetzt, und die Gleichung (41) wird für δψ linearisiert. Dann wird die Stabilität der linearisierten Gleichung betrachtet. Wenn jedes Glied, das nicht δψ enthält, als eine Konstante weggelassen wird, erhält man:
δψ¹(M,N) = δt·λ(M,N) · {Aδψ&sup0;(M,N-1) + Bδψ&sup0;(M-1,N) + Cδψ&sup0;(M,N+1) + Dδψ&sup0;(M+1,N) + Hδψ&sup0;(M,N) - αδψ&sup0;(M,N)} + δψ&sup0;(M,N) ... (54-1)
A = 2/[hy(M'-1)·hy'(N)] ... (54-2)
B = 2/[hx(M'-1)·hx'(M)] ... (54-3)
C = 2/[hy(N')·hy'(N)] ... (54-4)
D = 2/[hx(M')·hx'(M)] ... (54-5)
H = 1/[hx(M'-1)hx'(M)] + 1/[hx(M')hx'(M)] + 1/[hy(N'-1)hy'(N)] + 1/[hy(N')hy'(N)] ... (54-6)
α = (qθni/&epsi;) · [exp{-θψ&sup0;(M,N)] + exp{θ(ψ&sup0;(M,N)-V)}] ... (54-7)
Diese Gleichungen werden in die folgende Vektorgleichung (54-8) kombiniert:
ψ1 = ψ ... (54-8)
Gemäß dem Gerschgorin'schen Theorem sind die Koordinaten der durch das Theorem definierten Mitte von Platten reelle Zahlen, wenn die Bedingung zum Stabilisieren der Gleichung (42) erfüllt ist, d. h. wenn ein Eigenwert der Fehlerfortpflanzungsmatrix einen Absolutwert von 1 oder weniger aufweist. In diesem Fall ist die Mitte der Platten auf der reellen Achse. Somit ist es ausreichend, die folgenden zwei Ungleichungen zu erfüllen:
1 - λ(M,N)·δt(H+α) - λ(M,N)·δt·H > - 1 ... (54-9)
1 - λ(M,N)·δt(H+α) + λ(M,N)·δt·H < 1 ... (54-10)
In der Ungleichung (54-10) ist -λ(M,N)δtα offensichtlich kleiner als 0. Daher ist es die Ungleichung (54-9), die die Bedingung der Stabilität für die Gleichung (42) gewährleistet. Mit anderen Worten ist es ausreichend, um die Gleichung (42) zu stabilisieren, das Zeitinkrement Δt und die Funktion λ, die Funktionen des Raumintervalls H bzw. des nichtlinearen Glieds α sind, durch die folgende Gleichung (54-11) zu bestimmen.
λ(M,N)·δt ≤ 2/(2H+α) ... (54-11)
Das Raumintervall H ist eine Funktion, die sich nur auf die Position bezieht und überhaupt nicht von der Zeit abhängt, wohingegen das nichtlineare Glied α, der ψ enthält, eine Funktion der Position und Zeit ist. Daher ist λ doch eine Funktion von sowohl Position als auch Zeit.
Wenn die obige Erläuterung auf ein eindimensionales System oder ein zweidimensionales System angewendet wird, ist es nur das Raumintervall H, das sich ändert. Das Verfahren des Bestimmens von λ kann nicht nur auf die Berechnung der finiten Differenz angewendet werden, sondern auch auf weitere diskretisierende Verfahren (beispielsweise wie Finite- Elemente-Verfahren und Finite-Volumen-Verfahren), ohne modifiziert zu werden.
Wenn dieses Verfahren in der Praxis angewendet wird, wird λδt auf einen Wert geringfügig kleiner als die rechte Seite der Gleichung (54-11) eingestellt. Genauer gesagt wird λδt mit dem Parameter ω (eine positive Zahl kleiner als und nahezu gleich 1) wie folgt definiert:
λ·δt = 2ω/(2H+α) ... (54-12)
Die Bedeutung der Gleichung (54-12) wird qualitativ beschrieben. Im allgemeinen ist entweder die Elektronendichte n oder die Lochdichte p nahezu gleich der Verunreinigungskonzentration in dem Hochverunreinigungsbereich einer Halbleitervorrichtung, und ein neutraler Zustand der Raumladung wird in diesem Bereich beibehalten. In diesem Fall ist, da α weitaus größer als H und der Zähler 2 ist, λ(M,N)δt extrem klein. Folglich sind die Änderungen ψ¹(M,N) - ψ&sup0;(M,N) ebenfalls extrem klein. In einem Bereich, der eine relativ niedrige Verunreinigungskonzentration aufweist, wo es wahrscheinlicher ist, dass eine Sperrschicht gebildet wird, ist α nahezu gleich H oder vernachlässigbar klein im Vergleich mit H, und λ(M,N)δt ist relativ groß. In diesem Fall sind die Änderungen ψ¹(M,N) - ψ&sup0;(M,N) relativ groß.
Aus den oben angegebenen Gründen ändert sich die Variable ψ¹(M,N) - ψ&sup0; (M,N) über die Zeit in den Niedrigverunreinigungsbereichen, vorausgesetzt, dass ein Wert, der den neutralen Zustand für die Raumladung erfüllt, der variablen Größe ψ(M,N) als eine Anfangsbedingung zum Starten der Berechnen gegeben wurde. Als ein Ergebnis dieses wird die Bildung einer Sperrschicht gefördert, und die variable Größe ψ(M,N) konvergiert auf die wahre Lösung.
Diese Berechnung kann natürlich mittels eines großen digitalen Computers durchgeführt werden, wobei sie jedoch ebenfalls durch ein Gerät mit elektronischen Schaltungen bzw. fester Verdrahtung gelöst werden kann, dessen Zeitantwort (d. h. Einschwingphänomen) einen stationären Zustand erreicht, womit spezifische Gleichungen gelöst werden, wie es später in Verbindung mit der Beschreibung der erfindungsgemäßen Geräte erläutert wird.
Das Gerät, das das Verfahren der vierten Ausführungsform durchführt, wird beschrieben. Das Gerät umfasst NT logische Zellen und NT Verstärker mit einem Verstärkungsfaktor λi, die jeweils mit den logischen Zellen verbunden sind, wobei NT die Gesamtzahl festgelegter Maschenpunkte ist, wenn die ursprüngliche Kontinuitätsgleichung eines kontinuierlichen Systems in ein diskretes System mittels einer Differenznäherung oder Finite-Elemente-Näherung transformiert wird, um durch einen Computer gelöst zu werden. Die Verstärker sind ausgestaltet, um einen Verstärkungsfaktor aufzuweisen, der groß genug ist, um tatsächliche Probleme zu bewältigen. Der Verstärkungskoeffizient wird aus der Eingabe mittels eines Addierers, eines Multiplizierers und eines Dividierers bestimmt, wie es aus Gleichung (54-12) ersichtlich ist. Das Gerät umfasst ferner eine Wahlvorrichtung (Dial) oder einen Controller zum Einstellen des Parameters ω, der in der Gleichung (54-12) enthalten und dem gesamten Schaltungsnetzwerk zugeordnet ist.
Teil der Gleichung (54-12) zum Berechnen des Verstärkungskoeffizienten kann ein Näherungsausdruck abhängig von dem Wert von α oder H sein. Ferner kann der Verstärkungskoeffizient aus der Eingabe in die Schaltung bestimmt werden.
Die Verstärkungskoeffizienten des Verstärkers können abhängig von der Größe der zu lösenden Aufgabenstellung kleiner als 1 sein, und die Verstärker können als Dämpfungsglieder fungieren.
Multiplizierer, die eine Eingabe mit einem Multiplizierer λi multiplizieren, können anstelle der Verstärker verwendet werden.
Eine numerische Analyse wird erläutert, die tatsächlich durchgeführt wurde, um die Poisson-Gleichung, die die Beziehung zwischen dem Potential und der Raumladung darstellt, mittels eines Computers zu lösen, womit die Modellierung einer Halbleitervorrichtung durchgeführt wurde, wobei ein pn-Übergang die in den Fig. 12A bis 12C spezifizierten Abmessungen und Verunreinigungsverteilung aufweist. Um diese Vorrichtung als ein dreidimensionales Problem zu analysieren, wird die Vorrichtung in 21840 Maschen (= 20 · 26 · 42) aufgeteilt. Um die Vorrichtung als ein zweidimensionales Problem zu analysieren, wird die x-z-Ebene in 840 Maschen (= 20 · 42) aufgeteilt. Mit Null Vorspannung wird der Probewert wie folgt eingestellt:
Tabelle 2 und Fig. 13A bis 13E zeigen, wie eine numerische Konvergenz stattfindet und welche Beziehung die positive Iterationzahl und der Parameter ω bei einem dreidimensionalen System aufweisen, auf das die Empfindlichkeitskoeffizienten λ durch das oben beschriebene Verfahren angewendet wurden. Hier wird der Fehler w definiert als:
w = δψ/(ψ&sup0;+δψ) ∞
Die Iteration wird durchgeführt, bis der Fehler w den Wert von 1.E-8 annimmt. Tabelle 2 und die Fig. 13A bis 13E zeigen die folgenden Tatsachen:
1) Wenn ω klein ist, muss die Iterationzahl groß sein, um eine Fehlerkonvergenz zu erreichen. Je größer der Parameter ω ist, desto schneller ist die Konvergenz. Wenn jedoch ω > 0,9 ist, wird die Konvergenzgeschwindigkeit sehr stark verringert.
2) Wenn ω < 1,0 ist, konvergiert der Fehler monoton, obgleich der Fehler w am Anfang der Iteration stark fluktuiert.
3) Wenn ω = 1,0 ist, vermindert sich der Fehler w, wobei jedoch keine Konvergenz auftritt.
4) Wenn ω > 1,0 ist, tritt keine Konvergenz auf; der Fehler w fluktuiert unregelmäßig.
Wie es aus Tabelle 2 ersichtlich ist, konvergiert der Fehler w am schnellsten, wenn ω = 0,9 ist. Tabelle 2 Anzahl erforderlicher Iterationen (dreidimensionales Modell, Null-Vorspannung)
*Anmerkung: I. E. steht für die Iterationzahl
Tabelle 3 und Fig. 13F zeigen, wie eine numerische Konvergenz stattfindet und welche Beziehung die positive Iterationzahl und der Parameter ω bei einem zweidimensionalen System aufweisen, auf das die Empfindlichkeitskoeffizienten λ von dem oben beschriebene Verfahren angewendet wurden. Der Fehler w konvergiert auf die gleiche Art und Weise wie in dem Fall des dreidimensionalen Modells (Null Vorspannung), wobei er jedoch am schnellsten konvergiert, wenn ω = 0,88 ist. Tabelle 3 zeigt, dass der Wert des Parameters ω die Berechnungszeit beeinflusst, wenn eine positive Iteration durchgeführt wird, und dass sich die Iterationzahl nur ein wenig ändert, vorausgesetzt, dass der Parameter ω auf etwa 0,9 eingestellt ist. Tabelle 3 Anzahl erforderlicher Iterationen (zweidimensionales Modell, Null Vorspannung)
*Anmerkung: I. E. steht für die Iterationzahl
Die Zeit, die erforderlich ist, um zwei- und dreidimensionale Modelle bei Null-Vorspannung durch das Verfahren der Erfindung zu analysieren, wird mit der Zeit verglichen, die erforderlich ist, um dasselbe durch das herkömmliche Verfähren zu analysieren, das negative Iteration beinhaltet. Das herkömmliche Verfahren besteht aus dem Newton'schen Verfahren zum Lösen nichtlinearer Gleichungen; wobei die Gauss'sche Eliminierung, das BCG-Verfahren und das CR-Verfahren angewendet werden, um massive simultane Gleichungen zu lösen, was die in Tabelle 4 spezifizierten Zeitspannen erfordert. In Tabelle 4 sind ebenfalls die Zeitspannen spezifiziert, die das Verfahren dieser Erfindung erfordert, um zwei- und dreidimensionale Modelle zu analysieren. Wie es aus Tabelle 4 offensichtlich ist, kann das erfindungsgemäße Verfahren nichtlineare Gleichungen zwei- bis viermal schneller als das herkömmliche Verfahren lösen. Bei dem Newton'schen Verfahren ist der Fehler zuerst extrem groß und konvergiert plötzlich einige Zeit später. Bei dem Verfahren der Erfindung verringert sich der Fehler allmählich, um eine Lösung zu ergeben, und der erlaubte Fehler w, der ein Kriterium des Bestimmens der numerischen Konvergenz ist, kann so groß wie 1.E-4 sein. Im Hinblick darauf ist die positive Iteration offensichtlich vorteilhaft gegenüber der bei dem herkömmlichen Verfahren verwendeten negativen Iteration. Tabelle 4 - Berechnungszeit
Die durch ein erfindungsgemäßes Berechnungsverfahren erhaltenen Ergebnisse mit einer umgekehrten Vorspannung, werden nun mit Bezug auf Fig. 14 beschrieben. Die umgekehrte Vorspannung VR wird schrittweise von 1 V auf 200 V erhöht. Die Probewerte, die angewendet wurden, um die Lösung zu erhalten, werden wie folgt basierend auf der Lösung bestimmt, die durch Anlegen einer Vorspannung V1 bei dem direkt vorhergehenden Schritt und durch die Differenz zwischen der angelegten Vorspannung V2 und der Vorspannung V1, d. h. V2 - V1, erhalten wurde:
Fig. 14 zeigt ferner die Beziehung zwischen der Differenz, die während des Berechnens des Falls der umgekehrten Vorspannung erhalten wurde, und der positiven Iterationzahlen. Obgleich die Iterationzahl, die erforderlich ist, um die Differenz zu konvergieren, ansteigt, wenn die angelegte Spannung ansteigt, konvergiert die Differenz in jedem Fall, wodurch sich eine Lösung ergibt. Dies zeigt, dass das Verfahren mit der Erfindung verwendet werden kann, um einen Fall einer umgekehrten Vorspannung zu berechnen.
Ein Verfahren zum Bestimmen der Empfindlichkeitskoeffizienten wird nun beschrieben, das das Verfahren der erfindungsgemäßen fünften Ausführungsform ist und ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung zu analysieren. Bei dem Verfahren der fünften Ausführungsform der Erfindung werden, um die Gleichungen (8-1) (8-2) und (9) zu lösen, diese Gleichungen in die folgenden Gleichungen transformiert, die jeweils ein Zeitdifferentialglied enthalten:
dp(M,N)/dt = -λp(M,N)fp(M,N) ... (20-1)
dn(M,N)/dt = λn(M,N)fn(M,N) ... (20-2)
dψ(M,N)dt = -λψ(M,N)fψ(M,N) ... (20-3)
wobei fp, fn und fψ, die in den rechten Seiten dieser Gleichungen sind, wie folgt definiert sind:
fp(M,N) = (1/q) {[Jpx(M) - Jpx(M-1)]/hx'(M)} + (1/q){[Jpy(N) - Jpy(N-1)]/hy'(M)} + U(M,N) ... (21-1)
fn(M,N) = (1/q) {[Jnx(M) - Jnx(M-1)]/hx'(M)} + (1/q) {[Jny(N) - Jny(N-1)]/hy'(M)} + U(M,N) ... (21-2)
fψ(M,N) = [1/hx'(M)] · {[ψ(M,N) - ψ(M,N)]hx(M) - [ψ(M,N) - ψ(M-1,N)]/hx(M-1)} + [1/hy'(N)] · {[ψ(M,N+1) - ψ(M,N)]/hy(M) - [ψ(M,N) - ψ(M,N-1)]/hy(N-1)} + (q/&epsi;) [Γ(M,N) + p(M,N) - n(M,N)] ... (21-3)
Die Zeitdifferentialglieder auf den linken Seiten der Gleichungen (20-1), (20-2) und (20-3) werden durch Setzen einer endlichen Anzahl von Maschenpunkten auf der Zeitachse differenzgenähert, wodurch die folgenden Gleichungen erhalten werden:
p¹(M,N) = p&sup0;(M,N) - δtλp(M,N)fp&sup0;(M,N) ... (36-1)
n¹(M,N) = n&sup0;(M,N) - δtλn(M,N)fn&sup0;(M,N) ... (36-2)
ψ¹(M,N) = ψ&sup0;(M,N) - δtλψ(M,N)fψ&sup0;(M,N) ... (36-3)
wobei 1 ≤ M ≤ K, 1 ≤ N ≤ L ist. In den Gleichungen (36- 1), (36-2) und (36-3) bezeichnet die hochgestellte 0 einen alten Zeitpunkt und die hochgestellte 1 einen neuen Zeitpunkt.
Das Verfahren der fünften Ausführungsform zum Analysieren der Funktionsweise einer Halbleitervorrichtung wird nun mit Bezug auf das Ablaufdiagramm von Fig. 15 erläutert. Zuerst werden bei Schritt a5 die Datenelemente, die die Struktur und die Verunreinigungskonzentrationen der Halbleitervorrichtung darstellen, in den Speicher eines Rechners eingegeben. Als nächstes werden bei Schritt b5 die Positionen der Maschenpunkte der Vorrichtung bestimmt. Bei Schritt c5 werden die Daten, die eine Vorspannung V zeigen, in den Speicher des Computers eingegeben. Der Ablauf geht zu Schritt d5, bei dem anfängliche Probewerte p&sup0;(M,N), n&sup0;(M,N), ψ&sup0;(M,N) für die variablen Grundgrößen p(M,N), n(M,N) und ψ(M,N) angewendet werden. Ferner wird bei Schritt e5 eine Integrationsstartzeit t¹ eingestellt.
Als nächstes werden bei Schritt f5 die Werte der rechten Seiten der Gleichungen (36-1), (36-2) und (36-3) für alle Maschenpunkte 1 ≤ M ≤ K, 1 ≤ N ≤ L erhalten. Die somit erhaltenen Werte werden als gleich p¹(M,N), n¹(M,N), ψ¹(M,N); 1 ≤ M ≤ K, 1 ≤ N ≤ L betrachtet. Somit rückt die Integration auf der Zeitachse einen Schritt vorwärts.
Bei dem nächsten Schritt, d. h. Schritt g5, werden die Absolutwerte der Änderungen in den variablen Größen p¹(M,N) - p&sup0;(M,N) , n¹(M,N) - n&sup0; (M,N) , ψ¹(M,N) - ψ&sup0;(M,N) für alle Maschenpunkte erhalten. Ferner wird bei Schritt g5 bestimmt, ob die Absolutwerte dieser Änderungen für alle Maschenpunkte in vorbestimmte Fehlergrenzen fallen oder nicht.
Wenn es auch nur ein NEIN bei einem Maschenpunkt bei Schritt g5 gibt, geht der Fluss zu Schritt h5, wobei die korrigierten Werte p¹(M,N), n¹(M,N), ψ¹(M,N) als neue Probewerte verwendet werden, und folglich als gleich als p&sup0;(M,N), n&sup0;(M,N), ψ&sup0;(M,N) angesehen werden. Dann wird die Zeit t¹ auf t&sup0; eingestellt, und die Integration wird auf der Zeitachse durchgeführt.
Umgekehrt wird angenommen, wenn bei Schritt g5 bei allen Maschenpunkten JA ist, dass stationäre Lösungen erhalten wurden. Diese Lösungen sind von einer solchen Art, dass die Zeitdifferentialglieder auf den linken Seiten der Gleichungen (20-1), (20-2) und (20-3) Null werden, d. h. fp(M,N) = 0, fn(M,N) = 0, fψ(M,N) = 0; 1 ≤ M ≤ K, 1 ≤ N ≤ L. Dies bedeutet, dass die Lösungen in den Gleichungen (8-1), (8-2) und (9) erhalten wurden.
Bei dem Verfahren der fünften Ausführungsform können die Empfindlichkeitskoeffizienten λp(M,N), λn(M,N) und λψ(M,N) auf die gleiche Art und Weise wie bei dem Verfahren der vierten Ausführungsform bestimmt werden. Somit werden sie durch die folgende Gleichung (55-4) dargestellt. Da die Eigenwerte der Fehlerfortpflanzungsmatrix A in Einheitskreisen in einer komplexen Ebene existieren, sind λp, λn und λψ bestimmt.
Tatsächliche Gleichungen, die die Empfindlichkeitskoeffizienten λp(M,N), λn(M,N) und λψ(M,N) ausdrücken, die für die Erfindung bedeutsam sind, werden erläutert. Um das Verfahren der Erfindung erfolgreich durchzuführen, ist es von größter Bedeutung diesen Koeffizienten geeignete Werte zuzuordnen. Zu diesem Zweck wird das Gerschgorin'sche Kreistheorem angewendet.
Es wird nun erläutert, wie die Gleichungen (19-1), (19- 2) und (19-3), d. h. die Grundgleichungen der Erfindung, verarbeitet werden, um dadurch die Gleichungen (32-1), (32-2) und (32-3) herzuleiten, die die Bedingungen der sechsten Ausführungsform eines erfindungsgemäßen Verfahrens darstellen, das ausgestaltet ist, um die Betriebsbedingung einer Halbleitervorrichtung zu analysieren. Bei diesem Prozess wird das Gerschgorin'sche Kreistheorem angewendet. Genauer gesagt werden Gleichungen (19-1), (19-2) und (19-3) in die folgende allgemeine Form umgeschrieben:
Dann erhalten wir:
dx/dt = λf(x) ... (25)
Wenn ferner angenommen wird, dass sich die Variable x um δx&sub1; während der Zeitspanne δt ändert, dann:
δx&sub1; = λδtf(x&sub0;) ... (56)
Es sei ferner angenommen, dass sich die Variable x durch δx&sub2; während der nächsten Zeitspanne δt ändert, dann:
δx&sub2; = λδtf(x&sub0; + δx&sub1;) = λδt [f(x&sub0;) + (δf/δx)δx&sub1;] = δx&sub1; + λδt(δf/δx)δx&sub1; ... (57)
Somit wird die Fehlerfortpflanzungsmatrix A durch [1 + λδt(δf/δx)] definiert. Das zweite Glied der Gleichung (57) wird in das folgende umgeschrieben:
δp(M,N) = -λp(M,N)δt · {fp(M,N) + (∂n(M,N)/∂n(M,N-1))δn(M,N-1) + (∂fp(M,N)/∂p(M,N-1))δp(M,N-1) + (∂fp(M,N)/∂ψ(M,N-1))δψ(M,N-1) + (∂fp(M,N)/∂n(M,N))δn(M,N) + (∂fp(M,N)/∂p(M,N))δp(M,N) + (∂fp(M,N)/∂ψ(M,N))δψ(M,N) + (∂fp(M,N)/∂n(M,N+1))δn(M,N+1) + (∂fp(M,N)/∂p(M,N+1))δp(M,N+1) + (∂fp(M,N)/∂ψ(M,N+1))δψ(M,N+1) + (∂fp(M,N)/∂n(M-1,N))δn(M-1,N) + (∂fp(M,N)/∂p(M-1,N))δp(M-1,N) + (∂fp(M,N)/∂ψ(M-1,N))∂ψ(M-1,N) + (∂fp(M,N)/∂n(M+1,N))δn(M+1,N) + (∂fp(M,N)/∂p(M+1,N))δp(M+1,N) + (∂fp(M,N)/∂ψ(M+1,N))δψ(M+1,N)} ... (58-1)
δn(M,N) = λn(M,N)δt · {fn(M,N) + (∂fn(M,N)/∂n(M,N-1))δn(M,N-1) + (∂fn(M,N)/∂p(M,N-1))δp(M,N-1) + (∂fn(M,N)/∂ψ(M,N-1))δψ(M,N-1) + (∂fn(M,N)/∂n(M,N)δn(M,N) + (∂fn(M,N)/∂p(M,N))∂p(M,N) + (∂fn(M,N)/∂ψ(M,N))δψ(M,N) + (∂fn(M,N)/∂n(M,N+1))δn(M,N+1) + (∂fn(M,N)/∂p(M,N+1))δp(M,N+1) + (∂fn(M,N)/∂ψ(M,N+1))δψ(M,N+1) + (∂fn(M,N)/∂n(M-1,N))δn(M-1,N) + (∂fn(M,N)/∂p(M-1,N))δp(M-1,N) + (∂fn(M,N)/∂ψ(M-1,N))δψ(M-1,N) + (∂fn(M,N)/∂n(M+1,N)δn(M+1,N) + (∂fn(M,N)/∂p(M+1,N))δp(M+1,N) + (∂fn(M,N)/∂ψ(M+1,N)δψ(M+1,N)} ... (58-2)
δψ(M,N) = λψ(M,N)δt · {fψ(M,N) + (∂fψ(M,N)/∂n(M,N-1))δn(M,N-1) + (∂fψ(M,N)/∂p(M,N-1))δp(M,N-1) + (∂fψ(M,N)/∂ψ(M,N-1))δψ(M,N-1) + (∂fψ(M,N)/∂n(M,Nj)δn(M,N) + (∂fψ(M,N)/∂p(M,N))δp(M,N) + (∂fψ(M,N)/∂ψ(M,N))δψ(M,N) + (∂fψ(M,N)/∂n(M,N+1))δn(M,N+1) + (∂fψ(M,N)/∂p(M,N+1))δp(M,N+1) + (∂fψ(M,N)/∂ψ(M,N+1))δψ(M,N+1) + (∂fψ(M,N)/∂p(M-1,N))bn(M-1,N) + (∂fψ(M,N)/∂p(M-1,N))δp(M-1,N) + (∂fψ(M,N)/∂ψ(M-1,N))δψ(M-1,N) + (∂fψ(M,N)/∂n(M+1,N))δn(M+1,N) + (∂fψ(M,N)/∂p(M+1,N))δp(M+1,N) + (∂fψ(M,N)/∂ψ(M+1,N))δψ(M+1,N)} ... (58-3)
Es sei nun Gleichung (58-1) für die Variable p betrachtet. Wenn δn (,) und δψ(,) 0 sind, dann:
δp(M,N) = -λp(M,N)δt { ∂fp(M,N)/∂p(K,L)δp(K,L) + ∂fp(M,N)/δp(M,N)δp(M,N)} ... (59)
wobei wie folgt definiert ist:
(∂fp(M,N)/∂p(x,L))δp(K,L) = (∂fp(M,N)/∂p(M-1,N))δp(M-1,N) + (∂fp(M,N)/∂p(M+1,N))δp(M+1,N) + (∂fp(M,N)/∂p(M,N-1))δp(M,N-1) + (∂fp(M,N)/∂p(M,N+1))δp(M,N+1) ... (60)
Die folgende Fehlerfortpflanzungsgleichung (61), die Gleichung (57) entspricht, wird aus Gleichung (59) hergeleitet, die der Gleichung (56) entspricht:
δp²(M,N) = δp¹(M,N) - λp(H,N)δt · (∂fp(M,N)/∂p(K,L))δp¹(K,L) + ∂fp(M,N)/∂p(M,N)∂p(M,N)δp¹(M,N)} ... (61)
Wenn Gleichung (61) in eine Matrixgleichung transformiert wird, wird ein diagonales Glied gefunden, das ausgedrückt wird als:
1 - λp(M,N)δt(∂fp(M,N)/∂p(M,N)) ... (62)
Die Summe der absoluten Werte von nicht-diagonalen Gliedern wird wie folgt gegeben:
+ λp(M,N)δt ( ∂fp(M,N)/∂p(K,L) ) ... (63)
Hier wird das Gerschgorin'sche Kreistheorem angewendet.
Gemäß Definition sind die Komponenten der Fehlerfortpflanzungsmatrix reelle Zahlen, wenn die Eigenwerte für die Matrix innerhalb von Einheitskreisen in einer komplexen Ebene existieren. Es folgt, dass die Zentren Gi der Einheitskreise auf einer reellen Achse verteilt sind. Daher müssen die folgenden zwei Ungleichungen gelten:
- 1 ≤ 1 - λp(M,N)δt(∂fp(M,N)/∂p(M,N)) - λp(M,N)δt ∂fp(M,N)/∂p(K,L) ... (64)
1 - λp(M,N)δt(∂fp(M,N)/∂p(M,N)) + λp(M,N)δt ( ∂fp(M,N)/∂p(K,L) ) ≤ 1 ... (65)
Ein Umschreiben der Gleichung (64) ergibt:
2/{∂fp(M,N)/∂p(M,N) + ∂fp(M,N)/∂p(K,L) } ≥ λp(M,N)δt ... (66)
Dabei ändert sich die Gleichung (65) in:
0 ≥ λp(M,N)δt{( ∂fp(MIN)/∂p(K,L) ) - ∂fp(MIN)/∂p(M,N)} ... (67)
Es ist offensichtlich, dass λp(M,N)δt entweder einen positiven Wert annimmt, wenn der Wert in den Klammern negativ ist. Alternativ nimmt λp(M,N)δt einen negativen Wert an, wenn der Wert in den Klammern positiv ist. Um physikalische Anforderungen zu erfüllen, muss λp(M,N)δt einen positiven Wert aufweisen. Somit wird die Gleichung (67) nicht verwendet; es wird nur Gleichung (66) verwendet. In der Praxis wird der Parameter ωp (0 ≤ wp ≤ 1) verwendet. Das heißt:
λp(M,N)δt = 2ωp/{( ∂fp(M,N)/∂p(K,L) ) + ∂fp(M,N)/∂p(M,N)} ... (68)
Auf ähnliche Weise können λn (M,N)δt und δψ(M,N)δt für die Variablen n und ψ bestimmt werden.
Ein weiteres Verfahren zum Analysieren der Funktionsweise einer Halbleitervorrichtung wird nun beschrieben, was das Verfahren der sechsten Ausführungsform der Erfindung ist. Bei dem Verfahren der sechsten Ausführungsform werden, um Gleichungen (8-1), (8-2) und (9) zu lösen, diese Gleichungen in die folgenden transformiert, die jeweils ein Zeitdifferentialglied enthalten:
dp(M,N)/dt = -λp(M,N)fp(M,.N) ... (20-1)
dn(M,N)/dt = λn(M,N)fn(M,N) ... (20-2)
dψ(M,N)/dt = λψ(M,N)fψ(M,N) ... (20-3)
wobei fp, fn, und fψ, die in den rechten Seiten dieser Gleichungen sind, wie folgt definiert sind:
fp(M,N) = (1/q) {[Jpx(M) - Jpx(M-1)]/hx'(M)} + (1/q) ([Jpy(N) - Jpy(N-1)]/hy'(M)} + U(M,N) ... (21-1)
fn(M,N) = (1/q) {[Jnx(M) - Jnx(M-1)]/hx'(M)} +(1/q) {[Jny(N) - Jny(N-1)J/hy'(M)} + U(M,N) ... (21-2)
fψ(M,N) = [1/hx'(M)] {[ψ(M+1,N) - ψ(M,N)]/hx(M) - [ψ(M,N) - ψ(M-1,N)]/hx(M-1)} + [1/hy'(N)] {[ψ(M,N+1) - ψ(M,N)]/hy(N) - [ψ(M,N) - ψ(M,N-1)]/hx(N-1)} + (q/&epsi;) [Γ(M,N) + P(M,N) - n(M,N)] ... (21-3)
Die Zeitdifferentialglieder auf den linken Seiten der Gleichungen (20-1), (20-2) und (20-3) sind differenzgenähert, indem eine endliche Anzahl von Maschenpunkten auf der Zeitachse eingestellt werden, wodurch die folgenden Gleichungen erhalten werden:
p¹(M,N) = p&sup0;(M,N) - δtλp(M,N)fp&sup0;(M,N) ... (36-1)
n¹(M,N) = n&sup0;(M,N) + δtλp(M,N)fn&sup0;(M,N) ... (36-2)
ψ¹(M,N) = ψ&sup0;(M,N) + δtλψ(M,N)fψ&sup0;(M,N) ... (36-3)
wobei 1 ≤ M ≤ K, 1 ≤ N ≤ L ist. In den Gleichungen (36- 1), (36-2) und (36-3) bezeichnet die hochgestellte 0 einen alten Zeitpunkt und die hochstellte 1 einen neuen Zeitpunkt.
Das Verfahren der sechsten Ausführungsform zum Analysieren der Funktionsweise einer Halbleitervorrichtung wird nun mit Bezug auf das Ablaufdiagramm von Fig. 16 erläutert. Zuerst werden bei Schritt a6 die Daten-Elemente, die die Struktur und Verunreinigungskonzentrationen der Halbleitervorrichtung darstellen, in den Speicher eines Computers eingegeben. Als nächstes werden bei Schritt b6 die Positionen der Maschenpunkte der Vorrichtung bestimmt. Bei Schritt c6 werden die Daten, die eine Vorspannung V aufweisen, in den Speicher des Computers eingegeben. Der Ablauf geht zu Schritt d6, bei dem die anfänglichen Probewerte p&sup0;(M,N), n&sup0;(M,N), ψ&sup0;(M,N) für die variablen Grundgrößen p(M,N), n(M,N) und ψ(M,N) angewendet werden. Ferner wird bei Schritt e6 eine Integrationsstartzeit t¹ eingestellt.
Als nächstes werden bei Schritt f6 Empfindlichkeitskoeffizienten λp(M,N), λn(M,N) und λψ für alle Maschenpunkte in Übereinstimmung mit den Gleichungen (32-1), (32-2) und (32-3) berechnet. Ferner werden bei Schritt g6 die Werte der rechten Seiten der Gleichungen (36- 1), (36-2) und (36-3) für alle Maschenpunkte 1 ≤ M ≤ K, 1 ≤ N ≤ L erhalten. Die somit erhaltenen Werte werden als gleich p¹ (M,N), n¹(M,N) und ψ¹(M,N); 1 ≤ M ≤ K, 1 ≤ N ≤ L angesehen. Somit schreitet die Integration um einen Schritt auf der Zeitachse vorwärts.
Bei Schritt h6 werden die Absolutwerte der Änderungen in den variablen Größen p¹(M,N) - p&sup0;(M,N) , n¹(M,N) - n&sup0;(M,N) , ψ¹(M,N) - ψ&sup0;(M,N) für alle Maschenpunkte erhalten. Ferner wird bei Schritt h6 bestimmt, ob die Absolutwerte dieser Änderungen für alle Maschenpunkte in vorbestimmte Fehlergrenzen fallen oder nicht.
Wenn es auch nur ein NEIN bei einem Maschenpunkt bei Schritt h6 gibt, geht der Ablauf zu Schritt 16, bei dem die korrigierten Werte p¹(M,N), n¹(M,N), ψ¹(M,N) als neue Probewerte verwendet werden können, und folglich als p&sup0;(M,N), n&sup0;(M,N) und ψ&sup0;(M,N) angesehen werden. Dann wird die Zeit t¹ auf t&sup0; eingestellt, und die Integration wird auf der Zeitachse auf die gleiche Art und Weise wie bei dem Verfahren der fünften Ausführungsform durchgeführt.
Umgekehrt wird angenommen, wenn bei Schritt g1 bei allen Maschenpunkten JA ist, wird angenommen, dass stationäre Lösungen erhalten wurden. Diese Lösungen sind von einer solchen Art, dass die Zeitdifferentialglieder auf den linken Seiten der Gleichungen (20-1), (20-2) und (20-3) Null werden, d. h. fp(M,N) = 0, fn(M,N) = 0, fψ(M,N) = 0; 1 ≤ M ≤ K, 1 ≤ N ≤ L. Dies bedeutet, dass die Lösungen der Gleichungen (8-1), (8-2) und (9) erhalten wurden.
Bis jetzt wurde das Grundkonzept des erfindungsgemäßen Verfahrens der sechsten Ausführungsform beschrieben. Im Vergleich mit dem herkömmlichen Verfahren ist das Verfahren der sechsten Ausführungsform durch die folgenden Aspekte gekennzeichnet:
1) Das Verfahren umfasst keine Prozesse zum Linearisieren von nichtlinearen Gleichungen;
2) Keine massive Matrix-Vektorgleichung, wie beispielsweise Gleichung (16) oder die aus Gleichung (16) hergeleitete Gleichung, die schematisch in Fig. 2 gezeigt ist, muss gelöst werden, und somit müssen keine Berechnungen, die so komplex sind, dass große Speicher erforderlich sind, um diese durchzuführen, ausgeführt werden. Es ist ausreichend, die Werte der rechten Seiten der Gleichungen (36-1), (36-2) und (36-3) zu finden und die variablen Grundgrößen zu korrigieren.
3) Da die Koeffizienten λp(M,N), λn(M,N) und λψ(M,N), die als Einstellfaktoren (d. h. die Empfindlichkeitskoeffizienten) verwendet werden, unterschiedliche Werte für die Maschenpunkte des Vorrichtungsraums annehmen können, können die Lösungen der Gleichungen wirksam konvergieren, in dem die Empfindlichkeitskoeffizienten auf Werte eingestellt werden, die für die physikalischen Eigenschaften der Vorrichtung geeignet sind.
4) Der Prozess des Lösens der Gleichungen (20-1), (20-2) und (20-3) ist nichts anderes als eine numerische Integration auf der Zeitachse. Die Integration kann natürlich mittels eines großen digitalen Computers durchgeführt werden, wobei sie jedoch ebenfalls durch ein Gerät mit elektronischen Schaltungen bzw. fester Verdrahtung gelöst werden kann, das elektronische Schaltungen umfasst, deren Zeitantwort (d. h. Einschwingphänomen) einen stationären Zustand erreichen, womit spezifische Gleichungen gelöst werden, wie es in Verbindung mit der Beschreibung der erfindungsgemäßen Geräte erläutert wird.
Wie es aus dem obigen ersichtlich ist, ist das Verfahren der sechsten Ausführungsform dieser Erfindung neuartig und unterscheidet sich stark von dem herkömmlichen Verfahren, das aus Lösen von Matrix-Vektoraufgabenstellungen besteht, wobei es möglicherweise ein neues Gebiet in der Vorrichtungs- Modellierungstechnologie öffnet.
Ein Gerät zum Durchführen des Verfahrens der sechsten Ausführungsform, wie es oben erläutert ist, wird nun mit Bezug auf Fig. 17 und 18 beschrieben. Dieses Gerät ist ausgestaltet, um Differentialgleichungen (19-1), (19-2) und (19-3) zu lösen, die Zeitdifferentialglieder und Empfindlichkeits-Koeffizienten enthält, die aus simultanen Gleichungen (d. h. Elektronen- und Loch-Transportgleichungen und Poisson-Gleichung) umgebildet werden, um die Modellierung einer Halbleitervorrichtung zu erreichen.
Fig. 17 ist eine schematische Darstellung dieses Geräts zum Analysieren eines Halbleitervorrichtungsmodells 200. Wie es in dieser Figur gezeigt ist, umfasst das Gerät eine Adresseneinheit 202, eine Maschenpunkteinheit 204, eine Speichereinheit 206 und eine Beurteilungseinheit 208.
Die Maschenpunkteinheit 204 umfasst drei Eingangsanschlüsse und drei Ausgangsanschlüsse. Variable Grundgrößen p(M,N), n(M,N) und ψ(M,N) der Maschenpunkte werden an die Eingangsanschlüsse der Einheit 204 geliefert, und die drei aktualisierten variablen Größen werden von den Ausgangsanschlüssen der Einheit 204 geliefert. Der interne Aufbau der Maschenpunkteinheit 204 wird später beschrieben. Zusätzliche Maschenpunkte können verwendet werden. Mit anderen Worten kann das Gerät eine bis m Maschenpunkteinheit(en) aufweisen. Falls es insgesamt m Maschenpunkteinheiten aufweist, kann es m parallele Berechnungen simultan durchführen. Bei der folgenden Beschreibung sei angenommen, dass das Gerät eine Maschenpunkteinheit aufweist.
Die Adresseneinheit 202 ist ausgestaltet, um eine beliebige Speicherstelle der Speichereinheit 206 zu kennzeichnen, um dadurch die Ausgabe von der Maschenpunkteinheit 204 bei dieser Speicherstelle zu speichern.
Die Beurteilungseinheit 208 wird verwendet, um zu bestimmen, ob sich die Differenz zwischen dem vorherigen Wert jeder variablen Größe und deren aktualisiertem Wert für alle Maschenpunkte des Vorrichtungsmodells 200 ausreichend verringert hat. Wenn sich diese Differenz auch nur für einen der Maschenpunkte nicht ausreichend verringert hat, wird der aktualisierte Werte als einer betrachtet, der vor dem nächsten Zeitpunkt liegt, und wird somit an den Eingang der Maschenpunkteinheit geliefert. Andererseits wird angenommen, wenn sich die Differenz für alle Maschenpunkte des Vorrichtungsmodells 200 ausreichend verringert hat, dass stationäre Lösungen erhalten wurden, und die Berechnung wird beendet.
Die interne Struktur und der Betrieb der Maschenpunkteinheit 204, die die Hauptkomponente des in Fig. 17 gezeigten Geräts ist, wird mit Bezug auf Fig. 18 beschrieben. Die alten Werte p&sup0;(M,N), n&sup0;(M,N) und ψ&sup0;(M,N) der variablen Grundgrößen p(M,N), n(M,N) und ψ(M,N) für einen spezifizierten Maschenpunkt (M,N) werden an die drei Eingangsanschlüsse 210, 212 und 213 der Maschenpunkteinheit 204 geliefert. Diese Werte p&sup0;(M,N), n&sup0;(M,N) und ψ&sup0;(M,N) werden in Blöcke 216, 218 und 220 zur Berechnung der Empfindlichkeitskoeffizienten -λp, λn und λψ, und ferner in Blöcke 222, 224 und 226 zur Berechnung der Funktionen fp, fn und fψ eingegeben. Insbesondere wird der Wert p&sup0;(M,N), der sich auf die Kontinuitätsgleichung von Löchern bezieht, an die Eingangsanschlüsse der Blöcke 216 und 222 zum Berechnen von - λp und fp(M,N) geliefert. Zur gleichen Zeit werden die Werte p&sup0;(K,L), die die Lochdichte bei den benachbarten Maschenpunkten darstellen, ebenfalls an die Blöcke 216 und 222 mittels der Adresseneinheit 202 geliefert.
Die Kontinuitätsgleichung der Löcher wird angewendet, um die Lochdichte zu bestimmen. Die beiden anderen variablen Größen n(M,N) und ψ(M,N) und die variablen Größen n(K,L) und ψ(K,L) für die benachbarten Maschenpunkte werden an Blöcke 216 und 222 mittels der Adresseneinheit 202 geliefert, um fp(M,N) und λp(M,N) zu berechnen.
Nachdem die Funktionen und Empfindlichkeitskoeffizienten für die drei variablen Größen berechnet wurden, werden sie in Multiplizierer 228, 230 und 232 eingegeben. Die Ausgaben dieser Multiplizierer 228, 230 und 232 werden in Addierer 234, 236 und 238 eingegeben. Die Addierer 234, 236 und 238 addieren die Ausgaben der Multiplizierer 228, 230 und 232 zu den alten Werten der entsprechenden variablen Größen, wobei aktualisierte oder neue Werte der variablen Größen p(M,N), n(M,N) und ψ(M,N) erzeugt werden, die an die Ausgangsanschlüsse 240, 242 und 244 der Maschenpunkteinheit 204 geliefert werden.
Die Maschenpunkteinheit 204 kann entweder eine digitale Datenverarbeitungseinheit oder eine analoge Datenverarbeitungseinheit oder eine Hybrideinheit sein, die sowohl digitale Daten als auch analoge Daten verarbeiten kann.
Noch ein weiteres Verfahren wird beschrieben, das das Verfahren der erfindungsgemäßen siebenten Ausführungsform ist und ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung zu analysieren.
Die Lochdichte p und die Elektronendichte n, die entweder mit der Boltzmann-Statistik, einem quasi-Fermi- Lochpotentials φp oder einem quasi-Fermi-Elektronenpotentials φn beschrieben werden, weisen die folgenden Beziehungen auf:
p(M,N) = niexp [θ(φp(M,N) - ψ(M,N))] ... (69-1)
n(M,N) = niexp [θ(ψp(M,N) - φn(M,N)] ... (69-2)
Das Einsetzen dieser Gleichungen (69-1) und (69-2) in die Gleichungen (19-1) und (19-2) ergibt die folgenden zwei Gleichungen:
dφp(M,N)/dt = -λp(M,N)fp(M,N)/[θp(M,N)] + λψ(M,N)fψ(M,N) ... (70-1)
dφn(M,N)/dt = -λn(M,N)fn(M,N)/[θn(M,N)] + λψ(M,N)fψ(M,N) ... (70-2)
Diese beiden Gleichungen (70-1) und (70-2) werden anstelle der Gleichungen (19-1) und (19-2) angewendet. Mit anderen Worten werden drei Variablen φp(M,N), φn(M,N) und ψ(M,N) als variable Grundgrößen anstelle von p(M,N), n(M,N) und ψ(M,N) verwendet. Dies kennzeichnet das Verfahren der siebenten Ausführungsform der Erfindung. Das Verfahren der siebenten Ausführungsform der Erfindung ist mit dem Verfahren der sechsten Ausführungsform identisch, insoweit dass die angewendete dritte Gleichung ebenfalls die Gleichung (19-3) ist, nämlich:
dψ(M,N)/dt = λψ(M,N)fψ(M,N) ... (70-3)
Im Gegensatz zu dem Verfahren der sechsten Ausführungsform können drei Potentiale φp, φn und ψ, die die gleiche physikalische Dimension aufweisen, als Grundvariablen angesehen werden. Folglich haben ihre Werte eine derart ähnliche Größe, dass sie ohne weiteres numerisch verarbeitet werden können.
Das Verfahren der siebenten Ausführungsform wird nun ausführlich mit Bezug auf das Ablaufdiagramm von Fig. 19 beschrieben. Zuerst werden bei Schritt a7 die Daten-Elemente, die die Struktur und Verunreinigungskonzentrationen der Halbleitervorrichtung darstellen, in den Speicher eines Computers eingegeben. Als nächstes werden bei Schritt b7 die Positionen der Maschenpunkte der Vorrichtung bestimmt. Bei Schritt c7 werden die Daten, die eine Vorspannung V aufweisen, in den Speicher des Computers eingegeben. Der Ablauf geht durch Schritt d7, bei dem anfängliche Probewerte φp&sup0;(M,N), φn&sup0;(M,N) und ψ&sup0;(M,N) für die variablen Grundgrößen φp(M,N), φn(M,N) und ψ(M,N) angewendet werden. Ferner wird bei Schritt e7 eine Integrationsstartzeit t¹ eingestellt.
Als nächstes werden bei Schritt f7 die Werte der rechten Seiten der folgenden Gleichungen (71-1), (71-2) und (71-3) für alle Maschenpunkte 1 ≤ M ≤ K, 1 ≤ N ≤ L erhalten. Die somit erhaltenen Werte werden als gleich φp¹(M,N), φn¹(M,N) und ψ¹(M,N) angesehen; 1 ≤ M ≤ K, 1 ≤ N ≤ L. Folglich rückt die Integration einen Schritt auf der Zeitachse vor.
φp¹(M,N) = φp&sup0;(M,N) - δtλp(M,N)fp&sup0;(M,N)/[θp&sup0;(M,N)] + δtλψ(M,N)fψ&sup0;(M,N) ... (71-1)
φn¹ = φn&sup0;(M,N) - δtλp(M,N)fn&sup0;(M,N)/[θn&sup0;(M,N)] + δtλψ(M,N)fψ&sup0;(M,N) ... (71-2)
ψ¹(M,N) = ψ&sup0;(M,N) + δtλψ(M,N)fψ&sup0;(M,N) ... (71-3)
Bei Schritt g7 werden die Absolutwerte der Änderungen in den variablen Größen φp¹(M,N) - φp&sup0;(M,N) , φn¹(M,N) - φn&sup0;(M,N) , ψ¹(M,N) - φn&sup0;(M,N) für alle Maschenpunkte erhalten. Ferner wird bei Schritt g7 bestimmt, ob die Absolutwerte dieser Änderungen für alle Maschenpunkte in vorbestimmte Fehlergrenzen fallen oder nicht.
Wenn es auch nur ein NEIN bei einem Maschenpunkt bei Schritt g7 gibt, geht der Ablauf zu Schritt h7, bei dem die korrigierten Werte φp¹(M,N), φn¹(M,N) und ψ¹(M,N) als neue Probewerte verwendet und folglich als gleich φp&sup0;(M,N), φn&sup0;(M,N) und ψ&sup0;(M,N) angesehen werden. Dann wird die Zeit t¹ auf t&sup0; eingestellt, und die Integration wird auf der Zeitachse auf die gleiche Art und Weise wie bei dem Verfahren der sechsten Ausführungsform durchgeführt.
Umgekehrt wird angenommen, wenn bei Schritt g1 bei allen Maschenpunkten JA ist, dass stationäre Lösungen erhalten wurden. Diese Lösungen sind von einer solchen Art, dass die Zeit-Differentialglieder auf den linken Seiten der Gleichungen (70-1), (70-2) und (70-3) Null werden, d. h. fp(M,N) = 0, fn(M,N) = 0, fψ(M,N) = 0; 1 ≤ M ≤ K, 1 ≤ N ≤ L. Dies bedeutet, dass die Lösungen der Gleichungen (8-1), (8-2) und 9 erhalten wurden.
Ein Gerät, das ausgestaltet ist, um das Verfahren der siebenten Ausführungsform durchzuführen, wird nun mit Bezug auf Fig. 20 beschrieben. Dieses Gerät ist mit dem in Fig. 17 dargestellten mit Ausnahme der Struktur der Maschenpunkteinheit 204 identisch. Daher wird nur die Maschenpunkteinheit 204 mit Bezug auf Fig. 20 erläutert.
Die alten Werte φp&sup0; (M,N), φn&sup0;(M,N) und ψ&sup0;(M,N) der variablen Grundgrößen φp(M,N), φn(M,N) und ψ(M,N) für einen spezifizierten Maschenpunkt (M,N) werden an die drei Eingangsanschlüsse 310, 312 und 314 der Maschenpunkteinheit 204 geliefert. Diese Werte φp&sup0;(M,N), φn&sup0;(M,N) und ψ&sup0;(M,N) werden in Blöcke 322, 324, 326, 316, 318 und 320 zum Berechnen entsprechender Funktionen fp/[θp(M,N)], fn/[θn(M,N)] und fψ und ferner der Empfindlichkeitskoeffizienten - λp, λn und λψ, eingegeben. Genauer gesagt wird der Wert φp&sup0;(M,N), der sich auf die Kontinuitätsgleichung von Löchern bezieht, an die Eingangsanschlüsse der Blöcke 316 und 322 zum Berechnen von -λp und fp/[θp(M,N)] geliefert. Zur gleichen Zeit werden ebenfalls die Werte p&sup0;(K,L), die die Lochdichte an den benachbarten Maschenpunkten darstellen, ebenfalls an die Blöcke 316 und 322 mittels der Adresseneinheit 202 geliefert.
Die Kontinuitätsgleichung der Löcher wird angewendet, um die Lochdichte zu bestimmen. Die beiden anderen variablen Größen φn(M,N) und ψ(M,N) und die variablen Größen n(K,L) und ψ(K,L) für die benachbarten Maschenpunkte werden an Blöcke 316 und 322 mittels der Adresseneinheit 202 geliefert, um fp(M,N) und λp(M,N) zu berechnen.
Nachdem die Funktionen und Empfindlichkeitskoeffizienten für die drei variablen Größen berechnet wurden, werden sie in Multiplizierer 328, 330 und 332 eingegeben. Die Ausgaben der Multiplizierer 328 und 332 werden durch einen Addierer 334 addiert, und die Ausgaben der Multiplizierer 330 und 332 werden durch einen Addierer 336 addiert. Ferner wird die Ausgabe des Addierers 334 und der entsprechende alte Wert in einen Addierer 338 eingegeben; die Ausgabe des Addierers 336 und des entsprechenden alten Werts werden in einen Addierer 340 eingegeben; und die Ausgabe des Multiplizierers 332 und der entsprechende alte Wert werden in einen Addierer 342 eingegeben. Die Addierer 338, 340 und 342 geben aktualisierte oder neue Werte aus, die an die Ausgangsanschlüsse 344, 346 und 348 der Maschenpunkteinheit 204 geliefert werden.
Die Maschenpunkteinheit 204 kann entweder eine digitale Datenverarbeitungseinheit oder eine analoge Datenverarbeitungseinheit oder eine Hybrideinheit sein, die sowohl digitale Daten als auch analoge Daten verarbeiten kann.
Ein erfindungsgemäßes Verfahren der achten Ausführungsform wird beschrieben, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung zu analysieren.
Wenn simultane Gleichungen (70-1), (70-2) und (70-3) mittels Zeitintegration gelöst werden, wodurch stationäre Lösungen erhalten werden, werden die Zeitdiffetentialglieder der Gleichung eliminiert. Somit ist fp(M,N) = 0, fn(M,N) = 0, fψ(M,N) = 0. Diese drei Gleichungen sind nichts anderes als die stationären Gleichungen der Halbleitervorrichtung. Umgekehrt können, wenn diese stationären Gleichungen erfüllt sind, die folgenden Gleichungen (72-1) und (72-2), die durch Löschen der zweiten Glieder von den rechten Seiten der Gleichungen (70-1) und (70-2) gebildet werden, und die Gleichung (70-3) zeitintegriert werden, um stationäre Lösungen zu erhalten.
dφp(M,N)/dt = -λp(M,N)fp(M,N)/[θp(M,N)] ... (72-1)
dφn(M,N)/dt = -λn(M,N)fn(M,N)/[θn(M,N)] ... (72-2)
Ferner ergibt sich, wenn stationäre Lösungen aus den Gleichungen (72-1), (72-2) und (70-3) erhalten werden: fp(M,N) = 0, fn(M,N) = 0, fψ(M,N) = 0, genau so als wenn wir die Gleichungen (70-1), (70-2) und (70-3) Zeitintegrieren. Das Verfahren der achten Ausführungsform der Erfindung startet mit dem Zeitintegrieren der Funktionen (72-1), (72-2) und (70-3).
Das Verfahren der achten Ausführungsform wird nun ausführlich mit Bezug auf das Ablaufdiagramm von Fig. 21 erläutert. Zuerst werden bei Schritt aß die Daten-Elemente, die die Strukturen und Verunreinigungskonzentrationen der Halbleitervorrichtung darstellen, in den Speicher eines Computers eingegeben. Als nächstes werden bei Schritt b8 die Positionen der Maschenpunkte der Vorrichtung bestimmt. Bei Schritt c8 werden die Daten, die eine Vorspannung V zeigen, in den Speicher des Computers eingegeben. Der Ablauf geht zu Schritt d8, bei dem die anfänglichen Probewerte φp&sup0;(M,N), φn&sup0;(M,N), ψ&sup0;(M,N) für die variablen Grundgrößen φp(M,N), φn(M,N), ψ(M,N) angewendet werden. Ferner wird bei Schritt e8 eine Integrationsstartzeit t¹ eingestellt.
Als nächstes werden bei Schritt f8 die Werte der rechten Seiten der folgenden Gleichungen (73-1), (73-2) und (73-3) für alle Maschenpunkte 1 ≤ M ≤ K, 1 ≤ N ≤ L erhalten. Die somit erhaltenen Werte werden als gleich φp¹(M,N), φn¹(M,N), ψ¹(M,N); 1 ≤ M ≤ K, 1 ≤ N ≤ L angesehen. Somit rückt die Integration einen Schritt auf der Zeitachse vor.
φp¹(M,N) = φp&sup0;(M,N) - δtλp(M,N)fp&sup0;(M,N)/[θp&sup0;(M,N)] ... (73-1)
φn¹(M,N) = φn&sup0;(M,N) - δtλp(M,N)fn&sup0;(M,N)/[θn&sup0;(M,N)] ... (73-2)
ψ¹(M,N) = ψ&sup0;(M,N) + δtλψ(M,N)fψ&sup0;(M,N) ... (73-3)
Bei Schritt g8 werden die Absolutwerte der Änderungen in den variablen Größen φp¹(M,N) - φp&sup0;(M,N) , φn¹(M,N) - φn&sup0;(M,N) , ψ¹(M,N) - ψ&sup0;(M,N) für alle Maschenpunkte erhalten. Ferner wird bei Schritt g8 bestimmt, ob die Absolutwerte dieser Änderungen für alle Maschenpunkte in vorbestimmte Fehlergrenzen fallen oder nicht.
Bei NEIN bei Schritt g7, geht der Ablauf zu Schritt h7, bei dem die korrigierten Werte φp¹(M,N), φn¹(M,N), ψ¹(M,N) als neue Probewerte verwendet werden und folglich als gleich φp&sup0;(M,N), φn&sup0;(M,N) und ψ&sup0;(M,N) angesehen werden. Dann wird die Zeit t¹ auf t&sup0; eingestellt, und die Integration auf der Zeitachse wird auf die gleiche Art und Weise wie bei dem Verfahren der siebenten Ausführungsform durchgeführt.
Umgekehrt wird angenommen, wenn bei Schritt g8 bei allen Maschenpunkten JA ist, dass stationäre Lösungen erhalten wurden. Diese Lösungen sind von einer solchen Art, dass die Zeitdifferentialglieder auf den linken Seiten der Gleichungen (72-1), (72-2) und (72-3) Null werden, d. h. fP(M,N) = 0, fn(M,N) = 0, fψ(M,N); 1 ≤ M ≤ K, 1 ≤ N ≤ L. Dies bedeutet, dass die Lösungen der Gleichungen (8-1), (8-2) und (9) erhalten wurden.
Ein Gerät, das ausgestaltet ist, um das Verfahren der achten Ausführungsform durchzuführen, wird nun mit Bezug auf Fig. 22 beschrieben. Dieses Gerät ist mit dem in Fig. 17 dargestellten mit Ausnahme der Struktur der Maschenpunkteinheit 204 identisch, die die Hauptkomponente des Geräts ist. Folglich wird nur die Maschenpunkteinheit 204 mit Bezug auf Fig. 22 erläutert.
Die alten Werte φp&sup0;(M,N), φn&sup0;(M,N) und ψ&sup0;(M,N) der variablen Grundgrößen φp(M,N), φn(M,N), ψ(M,N) für einen spezifizierten Maschenpunkt (M,N) werden an die drei Eingangsanschlüsse 410, 412 und 414 der Maschenpunkteinheit 204 geliefert. Diese Werte φp&sup0;(M,N), φn&sup0;(M,N) und ψ&sup0;(M,N) werden in Blöcke 422, 424, 426, 416, 418 und 420 zum Berechnen entsprechender Funktionen fp/[θp (M,N)], fn/[θn(M,N)] und fψ und ebenfalls der Empfindlichkeitskoeffizienten - λp, λn und λψ, eingegeben. Genauer gesagt wird der Wert φp(M,N), der sich auf die Kontinuitätsgleichung von Löchern bezieht, an die Eingangsanschlüsse der Blöcke 416 und 422 zum Berechnen von - λp und fp/[θp(M,N)] geliefert. Zur gleichen Zeit werden die Werte p&sup0;(K,L), die die Lochdichte an den benachbarten Maschenpunkten darstellen, ebenfalls an die Blöcke 416 und 422 mittels der in Fig. 17 gezeigten Adresseneinheit 202 geliefert.
Die Kontinuitätsgleichung von Löchern wird angewendet, um die Lochdichte zu bestimmen. Die beiden anderen variablen Größen φn(M,N) und ψ(M,N) und die variablen Größen n(K,L) und ψ(K,L) für die benachbarten Maschenpunkte werden an Blöcke 416 und 422 mittels der Adresseneinheit 202 geliefert, um fp(M,N) und λp(M,N) zu berechnen.
Nachdem die Funktionen und Empfindlichkeitskoeffizienten für die drei variablen Größen berechnet wurden, werden sie in Multiplizierer 428, 430 und 432 eingegeben. Die Ausgaben der Multiplizierer 428, 430 und 432 werden in Addierer 434, 436 und 438 eingegeben und zu den entsprechenden alten Werten der variablen Größen addiert, womit aktualisierte oder neue Werte derselben geliefert werden. Die Ausgabe der Addierer 434, 436 und 438, d. h. die neuen Werte, werden an die Ausgangsanschlüsse 440, 442 und 444 der Maschenpunkteinheit 204 geliefert.
Die Maschenpunkteinheit 202 kann entweder eine digitale Datenverarbeitungseinheit oder eine analoge Datenverarbeitungseinheit oder eine Hybrideinheit sein, die sowohl digitale Daten als auch analoge Daten verarbeiten kann.
Ein erfindungsgemäßes Verfahren der neunten Ausführungsform wird beschrieben, das ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung zu analysieren.
Wie es dargelegt wurde, spielt es keine Rolle, ob die simultanen Gleichungen (70-1), (70-2) und (70-3) oder die simultanen Gleichungen (72-1), (72-2) und (70-3) zeitintegriert werden, die Ergebnisse sind stationär und zueinander gleich. Somit werden bei dem Verfahren der neunten Ausführungsform beide Sätze von simultanen Gleichungen gelöst, um dadurch eine Konvergenz der stationären Lösung mit höherer Geschwindigkeit als bei den Ausführungsformen des siebenten und achten Verfahrens zu erreichen.
Das Verfahren der neunten Ausführungsform wird nun ausführlich mit Bezug auf das Ablaufdiagramm von Fig. 23 erläutert. Zuerst wird bei Schritt a9 die Daten-Elemente, die die Struktur und Verunreinigungskonzentrationen der Halbleitervorrichtung darstellen, in den Speicher eines Computers eingegeben. Als nächstes werden bei Schritt b9 die Positionen der Maschenpunkte der Vorrichtung bestimmt. Bei Schritt c9 werden die Daten, die eine Vorspannung V zeigen, in den Speicher des Computers eingegeben. Der Ablauf geht zu Schritt d9, bei dem anfängliche Probewerte φp&sup0;(M,N), φn&sup0;(M,N) und ψ&sup0;(M,N) für die variablen Grundgrößen φp(M,N), φn(M,N) und ψ(M,N) angewendet werden. Ferner wird bei Schritt e9 eine Integrationsstartzeit t¹ eingestellt.
Als nächstes wird bei Schritt f9 bestimmt, ob IPOIS gleich 1 oder 0 ist. Wenn IPOIS gleich 1 ist, geht der Ablauf zu Schritt g9, was dem Schritt f7 des Verfahren der siebenten Ausführungsform entspricht. Wenn IPOIS gleich 0 ist, geht der Ablauf zu Schritt h9, was dem Schritt f8 des Verfahrens der achten Ausführungsform entspricht.
Bei Schritt g9 oder h9 werden die Werte der rechten Seite der Gleichungen (71-1), (71-2) und (71-3) oder der Gleichungen (73-1), (73-2) und (73-3) für alle Maschenpunkte 1 ≤ M ≤ K, 1 ≤ N ≤ L erhalten. Die somit erhaltenen Werte werden als φp¹(M,N), φn¹(M,N) und ψn¹(M,N); 1 ≤ M ≤ K, 1 ≤ N ≤ L angesehen. Somit rückt die Integration einen Schritt auf der Zeitachse vor.
Bei Schritt 19 werden die Absolutwerte der Änderungen in den variablen Größen φp¹(M,N) - φp&sup0;(M,N) , φn¹(M,N) - φn&sup0;(M,N) und ψ¹(M,N) - ψ&sup0;(M,N) für alle Maschenpunkte erhalten. Ferner wird bei Schritt i9 bestimmt, ob die Absolutwerte dieser Änderungen für alle Maschenpunkte in vorbestimmte Fehlergrenzen fallen oder nicht.
Bei NEIN bei Schritt 19, geht der Ablauf zu Schritt j9, bei dem die korrigierten Werte φp¹(M,N), φn¹(M,N) und ψ¹(M,N) als neue Probewerte verwendet werden, und somit als gleich φp&sup0;(M,N), φn&sup0;(M,N) und ψ&sup0;(M,N) angesehen. Dann wird die Zeit t¹ auf tº eingestellt, und die Integration wird auf der Zeitachse auf die gleiche Art und Weise wie bei den Verfahren der siebenten und achten Ausführungsformen durchgeführt.
Umgekehrt wird angenommen, wenn JA bei Schritt 19 ist, dass stationäre Lösungen erhalten wurden. Diese Lösungen sind von einer solchen Art, dass die Zeitdifferentialglieder auf den linken Seiten der Gleichungen (70-1), (70-2) und (70-3) oder der Gleichungen (72-1), (72-2) und (72-3) Null werden, d. h. fp(M,N) = 0, fn(M,N) = 0, fψ(M,N) = 0; 1 ≤ M ≤ K, 1 ≤ N ≤ L. Dies bedeutet, dass die Lösungen der Gleichungen (8-1), (8-2) und (9) erhalten wurden.
Ein Gerät, das ausgestaltet ist, um das Verfahren der neunten Ausführungsform durchzuführen, wird nun mit Bezug auf Fig. 24 beschrieben. Dieses Gerät ist mit dem in Fig. 19 dargestellten mit Ausnahme der Struktur der Maschenpunkteinheit 204 identisch, die die Hauptkomponente ist. Die Einheit 204 ist die gleiche wie diejenige in den Fig. 20 und 22 mit Ausnahme der Bereitstellung von Schaltern SW1, SW2, SW3, SW4, SW5 und SW6. Wenn die beweglichen Kontakte dieser Schalter mit den stationären Kontakten A verbunden sind, arbeitet die Maschenpunkteinheit 204 auf die gleiche Art und Weise wie die in Fig. 20 gezeigte. Wenn die beweglichen Kontakte mit den stationären Kontakten B verbunden sind, arbeitet die Einheit 204 auf die gleiche Art und Weise wie die in Fig. 22 dargestellte.
Ein Verfahren der zehnten Ausführungsform der Erfindung wird beschrieben, das ebenfalls ausgestaltet ist, um die Funktionsweise einer Halbleitervorrichtung zu analysieren. Dieses Verfahren ist dadurch gekennzeichnet, dass die Empfindlichkeitskoeffizienten λp, λn und λψ auf eine spezifische Art und Weise durch Umschreiben der variablen Grundgrößen φp, φa und ψ bestimmt werden. Diese variablen Grundgrößen können schließlich umgeschrieben werden zu:
Wie bei dem Verfahren der sechsten Ausführungsform sind ωp, ωn und ωψ Konstanten, die zwischen 0 und 1 liegen. Das Verfahren der zehnten Ausführungsform ist mit den Verfahren der siebenten und zehnten Ausführungsform identisch, mit der Ausnahme, dass die durch Gleichungen (35-1), (35-2) und (35- 3) definierten Empfindlichkeitskoeffizienten λp, λn und λψ nach Schritt e7 oder Schritt e9 gekennzeichnet werden.
Ein spezifisches Beispiel einer numerischen Analyse, die die Verfahren der siebenten bis neunten Ausführungsformen durchführen können, wird nun beschrieben. Die in diesem Beispiel angewendeten variablen Grundgrößen sind: φp, φn und ψ. Die Gleichungen (35-1), (35-2) und (35-3) werden angewendet, um die Empfindlichkeitskoeffizienten λp, λn und λψ zu bestimmen.
Zuerst werden Daten A, die die Querschnittsstruktur eines bipolaren Siliziumtransistors darstellen, der in Fig. 25A gezeigt ist, und Daten B, die das Verunreinigungsprofil des Transistors darstellen, das in Fig. 25B gezeigt ist, in den Computer eingegeben, der benutzt wird, um die Analyse durchzuführen. Basierend auf diesen Datenelementen A und B wird eine zweidimensionale Raummodellierung ausgeführt. Somit wird die Querschnittsansicht (Fig. 25A) des Transistors in 806 Maschenpunkte aufgeteilt, die eine 26 (x-Achse) · 31 (y- Achse) Maschenmatrix bildet. Ferner werden in den Computer Verunreinigungskonzentrationen WnE = 10²&sup0; cm&supmin;³, NpB = 10¹&sup8; cm&supmin;³, NnB = 10¹&sup7; cm&supmin;³, Nnc = 10²&sup0; cm&supmin;³ und NPX,ext = 10¹&sup9; cm&supmin;³ eingegeben. Npx,ext ist die Oberflächenverunreinigungskonzentration der externen Basisdiffusionsschicht.
Tabelle 5 zeigt die Ergebnisse, die die numerische Analyse ergibt, wenn die Basis-Emitter-Spannung VBE von 0 V auf 1 V erhöht wird, während die Basis-Kollektor-Spannung VBC auf einen konstanten Wert von 0 V gehalten wird. Die Spannung VBE ist in der ersten Spalte von Tabelle 5 gezeigt. In der zweiten Spalte ist die Konstante ω gezeigt, wobei ωp = ωn = ωψ, erfüllt ist und die von 0 bis 1 reicht, wie es beschrieben wurde. In der dritten Spalte ist der maximale Absolutwert δ max der Änderungen in φp, φn und ψ spezifiziert. Wenn dieser Wert δ max unter den vorbestimmten Wert absinkt, beispielsweise 10&supmin;&sup9; (oder "1.E-9" in Tabelle 5), wird angenommen, dass die stationäre Lösung der Gleichungen (8-1), (8-2) und (9) erhalten wurde.
In der vierten Spalte von Tabelle 5 ist IPOIS , gezeigt, das entweder 0 oder 1 ist. Wenn IPOIS = 0 ist, werden die Gleichungen (72-1) und (72-2), die das Glied λψfψ nicht enthalten, angewendet. Wenn IPOIS = 1 ist, werden die Gleichungen (70-1) und (70-2), die das Glied λψfψ enthalten, angewendet. In der fünften Spalte von Tabelle 5 ist Nit spezifiziert, die die Iterationzahl, d. h. die Anzahl von Schritten auf der Zeitachse ist.
Der Kollektorstrom ist in der sechsten Spalte gezeigt, wohingegen die Differenz zwischen dem Emitterstrom und dem Basisstrom in der siebenten Spalte gezeigt ist. Die Tatsache, dass der Kollektorstrom und die Differenz zwischen den Emittier- und Basis-Strömen gleich ist, zeigt, dass das Gesetz der Stromkonservierung vollständig gültig ist. Somit gilt, dass je größer die Ähnlichkeiten des Kollektorstroms und der Stromdifferenz ist, desto genauer ist die durch die Analyse erhaltene numerische Lösung. Der in der achten Spalte von Tabelle 5 spezifizierte Fehler stellt Ic/(IE - IB) - 1 in Prozent dar. In der letzten Spalte, d. h. der neunten Spalte, ist der Anfangs-VBE gezeigt, was der Anfangswert ist, der zum Berechnen der entsprechenden Reihe von Tabelle 5 eingestellt ist. In dem Fall der sechsten Reihe, bei der VBE = 0,7 V ist, wird ein endgültiges Berechnungsergebnis von VBE auf 0,65 V auf die Berechnung des Falls angewendet, wobei VBE = 0,7 V ein Anfangswert ist.
Die Ergebnisse der numerischen Analyse werden genauer untersucht. Wenn die Basis-Emitter-Spannung VBE relativ niedrig ist, d. h. von 0 V bis 0,6 V reicht, während IPOIS = 0 ist, kann die Lösung mit einer praktisch vernünftigen Anzahl von Iterationen erhalten werden. Wenn die numerische Analyse von einem Super-Computer durchgeführt wurde, wurden etwa 8 Millisekunden benötigt, um einen Arbeitsschritt zu beenden. Wie es in Tabelle 5 gezeigt ist, wird die Lösung von VBE = 0,4 V durch Durchführen von 3557 Arbeitsschritten erhalten. Daher ist eine Zeitspanne von etwa 29 Sekunden erforderlich, was in der Praxis ausreichend kurz ist. Im Hinblick darauf ist der in Tabelle 5 gezeigte dritte Fall, bei dem Ipois = 0 ist, ein wirksames Verfahren bei relativ niedrigem VBE. Tabelle 5
Der Fall, bei dem VBE = 0,7 V ist, wird nicht nur mit Bezug auf Tabelle 5 sondern auch auf Fig. 26 erläutert. In Fig. 20b gibt die durchgezogene Kurve die Ergebnisse des Analyseverfahrens der achten Ausführungsform an, wobei IPOIS = 0 ist, und die gestrichelte Kurve gibt die Ergebnisse des Analyseverfahrens der siebenten Ausführungsform, wobei IPOIS = 1 ist. Bei dem Verfahren der achten Ausführungsform, bei dem IPOIS = 0 ist, verringert sich δ max anfangs schnell auf 3 · 10&supmin;&sup0;&sup7; cm&supmin;³, wobei es jedoch danach langsam konvergiert: Andererseits verringert sich bei dem Verfahren der siebenten Ausführungsform, wobei IPOIS = 1 ist, δ max anfangs nicht schnell, sondern verringert sich kontinuierlich in dem unteren Bereich von δ max. Somit wird bei dem Verfahren der neunten Ausführungsform die Analyse gestartet, wobei IPOIS auf 0 eingestellt ist, und IPOIS wird auf 1 eingestellt, wenn sich δ max auf 10&supmin;&sup6; verringert, wodurch das δ max sehr schnell konvergiert, wie es durch die in Fig. 26 gezeigte gepunkt- gestrichelte Kurve angegeben ist.
Diese Erfindung ist nicht auf die obigen Ausführungsformen begrenzt. Der bei einer Ausführungsform offenbarte technische Schutzumfang kann mit den in den anderen Ausführungsformen offenbarten Schutzumfang kombiniert werden.
Zusätzliche Vorteile und Modifikationen werden Fachleuten ohne weiteres in den Sinn kommen. Daher ist die Erfindung in ihren breiteren Aspekten nicht auf die spezifischen Einzelheiten, dargestellten Gerät und veranschaulichten Beispielen begrenzt, die hier gezeigt und beschrieben wurden.

Claims

1. Verfahren zur Modellierung einer Halbleitervorrichtung mittels eines Computers, wobei das Verfahren die Schritte umfaßt:

(a) Bestimmen eines Arrays von Gitterpunkten (M,N), wobei jeder Gitterpunkt einer unterschiedlichen Raumposition innerhalb der Halbleitervorrichtung entspricht;

(b) Ausdrücken des Verhaltens der Vorrichtung an jedem Gitterpunkt mit einer Gleichung oder einer Mehrzahl simultaner Gleichungen der Form f(x) = 0, wobei x ein physikalischer Parameter der Vorrichtung ist, der durch die Modellierung zu bestimmen ist;

(c) Umschreiben der Gleichung oder der simultanen Gleichungen in die folgende Gleichung (a), die einen Zeitdifferentialausdruck dx/dt und einen Empfindlichkeits- Koeffizienten λ(M,N) für jeden Gitterpunkt hat, so daß eine Gleichung der allgemeinen Formel hergestellt wird:

dx(M,N)/dt = λ(M,N)f(x(M,N))

(d) Zuordnen eines Wertes, der für die physikalischen Zustände geeignet ist, zu dem Empfindlichkeits-Koeffizienten λ(M,N) bei jedem Gitterpunkt,

(e) Linearisieren der Gleichung (i) und Diskretisieren in Bezug auf Zeit und physikalische Größen, wobei Gleichung (i) in die folgende Vektorgleichung (ii) für jeden Gitterpunkt umgewandelt wird:

δ x(M,N)(t + df) = A δ x(M,N)t (ii)

(f) Integrieren der Gleichung (ii) über der Zeit, nachdem der Empfindlichkeits-Koeffizient so erhalten wurde, daß der Eigenwert der Fehlerfortpflanzungsmatrix A, die in Gleichung (ii) enthalten ist, kleiner als 1 ist, um so den physikalischen Parameter x(M,N) für jeden Gitterpunkt (M,N) zu erhalten.

2. Verfahren zur Modellierung einer Halbleitervorrichtung nach Anspruch 1, bei dem die zu bestimmenden physikalischen Parameter Konstanten einer Äquivalentschaltung einer Halbleitervorrichtung sind, und wobei das Verfahren ausgestaltet ist, um solche Parameter mit dem Verfahren des kleinsten Fehlerquadrats zu bestimmen, wobei der Schritt (c) den Schritt aufweist:

Umschreiben des Problems der Bestimmung der Konstanten der Äquivalentschaltung auf die folgende Gleichung (i), die einen Zeitdifferentialausdruck dq/dt und einen Empfindlichkeits-Koeffizienten λ für jeden Gitterpunkt enthält:

(i)

dq/dt = -λ(Aq - b)

wobei q die Änderung bei einem Gitterpunkt (M,N) ist und wobei A und b Konstanten der Äquivalentschaltung sind.

3. Verfahren nach Anspruch 1, das ausgestaltet ist, um eine Halbleitervorrichtung mittels der Poisson-Gleichung zu modellieren, die die Beziehung zwischen dem elektrischen Potential und der Raumladung darstellt, wobei der Schritt (b) den Schritt aufweist:

Umschreiben der Poisson-Gleichung in die folgende Gleichung (a), die Zeitdifferentialausdrücke ∂ψ/∂t und einen Empfindlichkeits-Koeffizienten λ enthält:

4. Verfahren nach Anspruch 3, bei dem der Schritt (d) den Schritt aufweist:

Bestimmen der Empfindlichkeits-Koeffizienten λ(M,N) durch die folgende Gleichung:

λ(M,N)δt ≤ 2/2H + α)

wobei δt ein diskretes Zeitintervall auf der Zeitachse ist, und wobei (M,N) die Nummer irgendeines diskreten Raumgitterpunkts ist.

5. Verfahren nach Anspruch 1, das ausgestaltet ist, eine Halbleitervorrichtung mit Elektronen- und Loch- Transportgleichungen und der Poisson-Gleichung zu modellieren, wobei die simultanen Gleichungen, die den Loch- Transport, den Elektronen-Transport und das Potential ausdrücken, in der Form für jeden Gitterpunkt (M,N) umgeschrieben werden:

wobei der Schritt (e) die Schritte umfaßt:

Linearisieren der Gleichungen (ia), (ib) und (ic), Diskretisieren bezüglich der Zeit und Umwandeln in die folgende Vektorgleichung entsprechend der Raumpositionsabhängigkeit:

6. Verfahren nach Anspruch 1, bei dem die folgenden Beziehungen angenommen sind:

fp = 1/q div Jp - (G - U)

fn = 1/q div Jn + (G - U)

fψ = div grad ψ + q/e(Γ + p - n)

n = ni exp[θ(ψ - φn)]

p = pi exp[θ(φp - ψ)]

θ = q/KT

und das weiter den Schritt des Umschreibens der simultanen Gleichungen umfaßt, die für jeden Gitterpunkt (M,N) zu lösen sind, in:

wobei der Schritt (f) den Schritt umfaßt:

Linearisieren der Gleichungen (id), (ie) und (if), Diskretisieren bezüglich der Zeit und Umwandeln in die folgende Vektorgleichung entsprechend der Raumpositionsabhängigkeit:

7. Verfahren nach Anspruch 6, wobei das Verfahren weiterhin den Schritt des Eliminierens des Ausdrucks λψfψ aus Gleichungen (id) und (ie) umfaßt.

8. Verfahren nach Anspruch 6, bei dem der Schritt (d) den Schritt des Bestimmens der Empfindlichkeits-Koeffizienten λp, λn, λψ durch die folgenden Gleichungen umfaßt:

wobei δt die diskreten Zeitintervalle an der Zeitachse ist, (M,N) die Nummer eines jeden diskreten Raumgitterpunktes ist, und wobei (K,L) die Nummer jedes benachbarten Punktes ist, der um den Raumgitterpunkt (M,N) herum existiert.

9. Verfahren nach einem der Ansprüche 6 oder 7, wobei im Schritt (d) die Empfindlichkeits-Koeffizienten λp, λn, λψ, durch die folgenden Gleichungen festgelegt sind:

10. Gerät zum Modellieren einer Halbleitervorrichtung, wobei die Vorrichtung durch eine Gleichung oder simultane Gleichungen der Form f(x) = 0 modelliert werden kann, wobei x ein physikalischer Parameter der zu bestimmenden Vorrichtung ist, wobei die Gleichung oder simultanen Gleichungen in der Form geschrieben sind:

dx/dt = λf(x)

wobei die Vorrichtung umfaßt:

erste Mittel zum Bestimmen der Adresse der Gitterpunkte (M,N) der Halbleitervorrichtung, so daß die Vorrichtung durch eine Gleichung oder simultane Gleichungen der Form dargestellt ist:

zweite Mittel zur Durchführung einer Integration in Übereinstimmung mit x und λ an den Gitterpunkten (M,N), deren Adresse durch das erste Mittel bestimmt wurde, wobei die zweiten Mittel Versuchswerte für x und λ zu einem Zeitpunkt t verwenden, wobei die Integration x zu einem Zeitpunkt t+δt bestimmt wird,

dritte Mittel zum Speichern von x zur Zeit t+δt, das von dem zweiten Mittel erhalten wurde; und

vierte Mittel zum Bestimmen, ob das von dem zweiten Mitteln erhaltene x richtig ist oder nicht, und um durch Vergleich von x zur Zeit t mit x zur Zeit t+δt zu sehen, ob die Änderung in x bei jedem Gitterpunkt über die Zeit δt innerhalb der erlaubten Fehlergrenze ist, um x zur Zeit t+δt zurück zu dem zweiten Mittel zu koppeln, so daß das zweite Mittel die Integration erneut durchführt, wenn x außerhalb der erlaubten Fehlergrenze fällt, und zum Ausgeben von x, wenn x innerhalb einer erlaubten Fehlergrenze ist.

11. Gerät nach Anspruch 10, bei dem das zweite Mittel umfaßt:

m logische Zellen;

m Verstärker mit einstellbarem Verstärkungs- Koeffizienten λi oder einem Multiplizierer, der den m Verstärkern äquivalent ist, verbunden mit den logischen Zellen, wobei m gleich der Anzahl der Gitterpunkte ist, die von dem ersten Mittel bestimmt wurden, und wobei

das vierte Mittel Detektornetzwerkmittel umfaßt, das mit den Logikzellen verbunden ist, um zu erfassen, daß die logischen Zellen eine Einschwingstabilität erreicht haben;

wobei die Vorrichtung des weiteren Mittel zum Bestimmen der Verstärkungsfaktoren λ(M,N) enthält, so daß eine Fehlerfortpflanzungsmatrix den Eigenwert kleiner als 1 zu einer bestimmten Zeit während einer Einschwingzeitspanne hat.

12. Gerät nach Anspruch 10 oder 11, das ausgestaltet ist, einen Halbleiter mit der Poisson-Gleichung zu modellieren.

13. Gerät nach Anspruch 10, das ausgestaltet ist, um für jeden Gitterpunkt Konstanten der Äquivalentschaltung einer Halbleitervorrichtung durch Umschreiben des Problems in eine Differentialgleichung dq/dt = -λ(Aq - b) zu bestimmen und dann die Differentialgleichung durch das Verfahren des kleinsten Fehlerquadrats zu lösen:

wobei das zweite Mittel umfaßt:

m Logikzellen;

m Verstärker mit einstellbaren Verstärkungsfaktoren λi oder einem Multiplizierer, der den Verstärkern äquivalent ist, verbunden mit den Logikzellen, wobei m gleich der Anzahl von Gitterpunkten (M,N) ist, die von dem ersten Mittel bestimmt sind; und

wobei das vierte Mittel Detektornetzwerkmittel umfaßt, das mit den Logikzellen verbunden ist, um zu erfassen, daß die Logikzellen eine Einschwingstabilität erreicht haben.

14. Gerät nach Anspruch 11, das ausgestaltet ist, um die Vorrichtung mit Simultangleichungen zu modellieren, die aus Elektronen- und Loch-Transportgleichungen und der Poisson- Gleichung bestehen, indem die simultanen Gleichungen in Differentialgleichungen umgeschrieben werden: dp/dt = -λpfp, dn/dt = λnfn, und dψ/dt = λψfψ.

15. Gerät nach Anspruch 10, das ausgestaltet ist, um die Vorrichtung mit simultanen Gleichungen, die aus Elektronen- und Loch-Transportgleichungen und der Poisson-Gleichung bestehen, zu modellieren, indem die simultanen Gleichungen in Differentialgleichungen umgeschrieben werden:

dφp/dt = -λpfp/Θp + λψfψ

dfn/dt = -λnfn/Θn + λψfψ

dψp/dt = λψfψ

die Zeit-Differentialausdrücke und Empfindlichkeits- Koeffizienten enthalten, und dann Lösen der Differentialgleichungen, um so die Modellierung einer Halbleitervorrichtung zu bewerkstelligen, wobei:

die zweiten Mittel ausgestaltet sind, um die Integration entsprechend φp, φn, ψ, λp, λn und λψ entsprechend jedem Gitterpunkt durchzuführen, dessen Adresse durch das erste Mittel bestimmt wurde, wodurch φp, φn und ψ beim Gitterpunkt erhalten werden;

wobei das dritte Mittel ausgestaltet ist, um φp, φn und ψ zu speichern, die so durch das zweite Mittel erhalten wurden, und

wobei das vierte Mittel ausgestaltet ist, um zu bestimmen, ob φp, φn und ψ, die von dem zweiten Mittel erhalten wurden, in die erlaubten Fehlergrenzen fallen oder nicht, und um φp, φn und ψ zurück zu dem zweiten Mittel zu koppeln, so daß das zweite Mittel die Integration erneut durchführt, wenn φp, φn und ψ außerhalb der erlaubten Fehlergrenzen fallen, und um φp, φn und ψ auszugeben, wenn φp, φn und ω innerhalb der erlaubten Fehlergrenzen fallen.

16. Gerät nach Anspruch 10, das ausgestaltet ist, um die Vorrichtung mit simultanen Gleichungen zu modellieren, die aus Elektronen- und Loch-Transportgleichungen und der Poisson-Gleichung bestehen, mit dem in Anspruch 3 definierten Verfahren, indem die simultanen Gleichungen in Differentialgleichungen umgeschrieben werden:

dfp/dt = -λpfp/Θp

dφn/dt = -λnfn/Θn

dψp/dt = λψfψ

die Zeitdifferentialausdrücke und Empfindlichkeits- Koeffizienten enthalten, und dann durch Lösen der simultanen Gleichungen, um so das Modellieren einer Halbleitervorrichtung zu bewerkstelligen;

wobei das zweite Mittel ausgestaltet ist, um die Integration in Übereinstimmung mit φp, φn, ψ, λp, λn und λψ entsprechend jedem Gitterpunkt durchzuführen, dessen Adresse durch das erste Mittel bestimmt wurde, wodurch φp, φn und ψ an jedem Gitterpunkt erhalten werden;

wobei das vierte Mittel ausgestaltet ist, um zu bestimmen, ob φp, φn und ψ, die von dem zweiten Mittel erhalten wurden, in die erlaubten Fehlergrenzen fallen, und φp, φn und ψ zu dem zweiten Mittel zurückzukoppeln, so daß das zweite Mittel die Integration erneut durchführt, wenn φp, φn und ψ nicht innerhalb der erlaubten Fehlergrenzen fallen.

17. Gerät nach Anspruch 10, das ausgestaltet ist, um die Vorrichtung mit simultanen Gleichungen zu modellieren, die aus Elektronen- und Loch-Transportgleichungen und der Poisson-Gleichung bestehen, indem die simultanen Gleichungen in Differentialgleichungen umgeschrieben werden:

dp/dt = -λpfp

dn/dt = λnfn

dψ/dt = λψfψ

die Zeitdifferentialausdrücke und Empfindlichkeits- Koeffizienten enthalten, und dann durch Lösen der Differentialgleichungen, um so die Modellierung einer Halbleitervorrichtung zu bewerkstelligen,

wobei das zweite Mittel ausgestaltet ist, um die Integration in Übereinstimmung mit p, n, ψ, λp, λn und λψ entsprechend jedem Gitterpunkt durchzuführen, dessen Adresse von dem ersten Mittel bestimmt wurde, wodurch p, n und ψ an dem Gitterpunkt erhalten werden; und

wobei das vierte Mittel ausgestaltet ist, um zu bestimmen, ob p, n und ψ, die von dem zweiten Mittel erhalten wurden, in erlaubte Fehlergrenzen fallen, und um p, n und ψ zurück zu dem zweiten Mittel zu koppeln, so daß das zweite Mittel die Integration erneut ausführt, wenn p, n und ψ außerhalb der erlaubten Fehlergrenzen fallen, und um p, n und ψ auszugeben, wenn p, n und y in die erlaubten Fehlergrenzen fallen.