DE4207435A1

DE4207435A1 - Topographie-simulationsverfahren

Info

Publication number: DE4207435A1
Application number: DE4207435A
Authority: DE
Inventors: Masato Fujinaga
Original assignee: Mitsubishi Electric Corp
Current assignee: Mitsubishi Electric Corp
Priority date: 1991-03-08
Filing date: 1992-03-09
Publication date: 1992-10-01
Also published as: KR960013371B1; JP2851447B2; JPH04280450A; KR920018924A; US5471403A

Description

Die Erfindung betrifft ein Topographie- Simulationsverfahren. Insbesondere betrifft die Erfindung ein Topographie-Simulationsverfahren, das den Rückfluß-Prozeß in der Halbleiter-Industrie dreidimensional simuliert.

Die hohe Integration von ULSI-Schaltungen ist fortgeschritten und als Ergebnis davon müssen die Elementstrukturen zu komplexen dreidimensionalen Strukturen mit Graben- oder Schicht-Technik übergehen. Um mit einem derartigen Trend fertig zu werden, wurde der Rückfluß von BPSG (Borphosphorsilikat-Glas) zur Einebnung von großen Stufen vor Verdrahtung und Zusammenschaltung verwendet. Nachdem die Oberflächen-Topographie des Rückflusses von der des abgelagerten BPSG abhängt, ist es notwendig, die Topographie von Ablagerung und Rückfluß in aufeinanderfolgenden Prozessen auszuwerten. Vor kurzem sind die Temperaturen von thermischen Prozessen auf weniger als 900° gefallen und die zu einebnenden Löcher sind kleiner als wenige Mikrometer. Daher wurden 3D-Modelle des Rückflusses für die niedrige Temperatur erforderlich.

Der Rückfluß von Dünnschicht-BPSG ist von F.A. Leon (1988) unter Verwendung eines Oberflächen-Diffusionsmodells simuliert worden.

Das bekannte Rückflußmodell ist in F.A. Leon, IEEE Trans. Electron CAD, Vol. 7, No. 2, Seiten 168 bis 173, 1988, beschrieben.

Unter Bezugnahme auf die Fig. 6 wird eine Topographie eines Borphosphorsilikat-Glases (BPSG) 1 durch Drahtsegment 2 beschrieben. Dieses BPSG 1 wird eine flüssigkeitsähnliche Flüssigkeit, wenn es heißer als 600° c wird. Die Fig. 7 ist eine Vergrößerung von A in der Fig. 6. Die Fig. 7 zeigt, wie der Krümmungsradius für diese Darstellung des Profils definiert wird. "R" beschreibt den Krümmungsradius in der Ecke des BPSG 1. Dieses Rückflußmodell zeigt, daß ein Materialteilchen von einem großen Krümmungsradius durch Oberflächendiffusion zu einem kleinen fließt. Und falls die Oberfläche durchgedrückt ist, wird der negative Krümmungsradius verwendet. Die nächste Gleichung zeigt die freie Energie:

F = γA,

wobei F die freie Energie der Oberfläche ist, γ die Oberflächenenergie pro Einheit und A die Fläche ist. Ein lokales chemisches Potential kann wie folgt definiert werden:

wobei N die Teilchenanzahl und Ω das molare Volumen ist. Im allgemeinen ändert sich der Krümmungsradius über der Flüssigkeitsoberfläche, wodurch Gradienten des chemischen Potentials erzeugt werden, die als Kräfte wirken. Diese Kräfte erzeugen einen Materialfluß, der wie folgt ausgedrückt werden kann:

wobei J der Oberflächen-Volumenfluß (Volumen des transportierten Materials durch eine Einheitslänge pro Einheitszeit) ist, γ die molare Oberflächenkonzentration, D der Oberflächen-Diffusionskoeffizient (mit Einheiten von Länge über Zeit), k die Boltzmann-Konstante, T die Temperatur und ∇s der Gradient bezüglich einer Strecke entlang der Oberfläche. R1(s) und R2(s) beschreiben den Krümmungsradius R1 und R2 in orthogonalen Richtungen wie in Fig. 9 gezeigt. Das Diagramm in Fig. 8 zeigt, wie die Bewegung des Punktes P₁₀ berechnet wird. Wenn R_i-1, R₁₀ und R_i-1 den Krümmungsradius an den Punkten P_i-1, bzw. P₁₀ bzw. P_i+1 bezeichnen, dann ist der Fluß vom Punkt P_i-1 zum Punkt P₁₀ gegeben durch

In der Fig. 8 wird gezeigt, wie sich der Punkt P_i0 zum Punkt P_i während eines Zeitinkrements Δt weiterbewegt. J_i-1,1 ist der Fluß vom Punkt P_i-1 zum Punkt P₁ und J_{i, i+1} ist der Fluß vom Punkt P_i zum Punkt P_i+1. Die Differenz der beiden Flüsse bestimmt die Strecke, die der Punkt P₁₀ zurücklegt. Bekannte Rückflußmodelle müssen den Krümmungsradius nach Unterteilung eines Abschnitts einer gekrümmten Oberfläche berechnen, die einen lokalen Krümmungsradius in orthogonalen Richtungen aufweist. Jedoch ist es sehr schwierig, in dieser Art und Weise den Krümmungsradius auf einer durch polygonartige Dreiecke beschriebenen dreidimensionalen Oberfläche zu berechnen. Das bekannte Rückflußmodell ist nämlich im Grunde genommen ein Modell für zwei, aber nicht für drei Dimensionen.

Die Erfindung betrifft ein Verfahren zur Topographie-Simulation, das den Schritt umfaßt das in einige Zellen unterteilte Feld zu analysieren, den Schritt des Zeichnens der Iso-Konzentrationsfläche, die durch den Satz von einigen Polygonen aus der Material-Konzentration von jeder Zelle beschrieben wird, den Schritt des Zeichnens des chemischen Potentials der Zelle aus dem Verhältnis des Material-Volumen-Zuwaches zum Material-Oberflächen-Zuwachs, wenn sich die Polygone in der Zelle, die die Iso-Konzentrationsfläche aufweist, in normaler Richtung von jedem Polygon bewegen, den Schritt des Zeichnens der Material-Stromdichte aus der Differenz von jedem chemischen Potential der Zelle, den Schritt des Zeichnens der Material-Konzentration der Zelle aus der Material-Stromdichte, den Schritt des Zeichnens einer neuen Iso-Konzentrationsfläche aus der Material-Konzentration. In der vorliegenden Erfindung kann daher das chemische Potential ohne den Krümmungsradius gezeichnet werden. In den Zeichnungen zeigt:

Fig. 1 ein Flußdiagramm, das das Topographie- Simulations-Verfahren des bevorzugten Ausführungsbeispiels der Erfindung beschreibt;

Fig. 2 eine Figur, die das Unterteilen des Analyse-Feldes in einige Zellen beschreibt;

Fig. 3 eine Figur, in der die würfelförmige Zelle in 24 Tetraheder unterteilt ist;

Fig. 4 eine Figur, die das Verfahren zum Zeichnen des chemischen Potentials erklärt;

Fig. 5 eine Figur, die beschreibt, wie die BPSG Moleküle von einer Position mit einem großen chemischen Potential zu demjenigen eines kleinen fließen;

Fig. 6 eine praktische Berechnung der Fig. 5;

Fig. 7 das Rückfluß-Modell im Vergleich mit einer SEM-Aufnahme;

Fig. 8 eine graphische Darstellung, die die Beziehung des maximalen Winkels der Oberfläche des BPSG-Films und der UND-Gesamtkonzentration von B203 und P205 zeigt;

Fig. 9 eine graphische Darstellung, die die Beziehung des Oberflächen- Diffusionskoeffzientens und der maximalen Winkel zeigt;

Fig. 10, 11 die Ergebnisse von Simulation und Experiment;

Fig. 12 eine Figur, die die bekannte Topographie-Simulation beschreibt;

Fig. 13 eine Vergrößerung von A in der Fig. 12;

Fig. 14 eine Figur von Drahtbewegungen; und

Fig. 15 eine Figur, die den Krümmungsradius des bekannten Verfahrens beschreibt.

Die Fig. 1 zeigt ein Flußdiagramm, das das Topographie-Simulationsverfahren des bevorzugten Ausführungsbeispiels der Erfindung beschreibt. Zunächst wird, wie in Fig. 2 gezeigt, das Analyse-Feld in einige rechteckige feste Zellen unterteilt, und durch den Schritt S1 eine Kontur-Oberfläche berechnet.

Die Material-Oberfläche wird als Kontur-Oberfläche aus gedrückt. Wie in Fig. 2 gezeigt, wird ein Analysebereich in würfelförmige Zellen unterteilt und die Material-Existenzrate (oder die Dichte) in dem Mittel jeder würfelförmigen Zelle angegeben. Als erstes werden die Existenz raten an den Ecken durch eine lineare Interpolation von den acht Werten in den Zellenmitten wie folgt berechnet:

C (x,y,z) = δ1 · C1 + δ2 · C2 + δ3 · C3 + δ4 · C4 + δ5 · C5 + δ6 · C6 + Δ7 · C7 + δ 8 · C8 = Vi/V0,

wobei angenommen wird, daß V0 das Würfelzellenvolumen ist, und Vi das Volumen des Würfels, der durch den Punkt x, y, z und durch die gegenüberliegende Ecke der Zellenmittel zur Ecke (i) ist, Ci die Existenzrate in dem Mittel der würfelförmigen, die benachbart zur Ecke x, y, z angeordnet ist. Als nächstes wird die würfelförmige Zelle in 24 Tetraheder unterteilt, wie in Fig. 3 dargestellt. Es gibt ein anderes Unterteilungsverfahren in 6 Tetraheder mit einer schnellen Berechnungsgeschwindigkeit, jedoch wird es hier nicht verwendet, weil die Topographie asymetrisch wird. Als drittes wird die Existenzrate in den Ecken jedes Tetraheders abgeleitet. Als viertes werden die Punkte von CO auf den Seiten in jedem Tetraheder durch lineare Interpolation berechnet, wenn CO eine Interface-Existenzrate (0,5) bezeichnet. Wenn diese Punkte untereinander durch Geraden in einen Tetraheder verbunden werden, wird ein Dreieck oder ein Quadrat abgeleitet. Die Kontur-Oberfläche wird durch diese Dreiecke und Quadrate zusammengebaut und eine 3D-Topographie wird als Kontur-Oberfläche definiert.

Das chemische Potential in jeder Oberflächenzelle wird im Schritt S3 berechnet.

Es wird angenommen, daß F die freie Helmholtz-Energie ist, Γ die Oberflächenspannung, A die Fläche der Material-Oberfläche, N die Anzahl der Moleküle und µ das chemische Potential einer kleinen Oberfläche. Temperatur und Druck werden als konstant angenommen.

dF = µ · dN = Γ · dA,

wobei dF, dN und dA die Inkremente von F bzw. N bzw. A sind, und außerdem die Existenz von thermischem Gleichgewicht in dem kleinen Oberflächenbereich. Hier wird das Volumeninkrement nicht beinhaltet, nachdem er der Bulk-Energie (Elastizitätsenergie) entspricht. Das chemische Oberflächenpotential kann wie folgt ausgedrückt werden:

µ = Γ · dA/dN = Γ · (dA/dV)(dV/dN)
= (Γ · Ω) · (dA/dV),

wobei Γ die Oberflächenspannung ist, Ω das Volumen eines Moleküls und dV das Volumeninkrement. Krümmung wird als das Verhältnis des Flächeninkrements zum Volumeninkrement definiert. Wenn diese Definition verwendet wird, ist es möglich, den Krümmungsradius für eine allgemeine 3D-Topographie zu berechnen, die aus Dreicken und Polygonen zusammengesetzt ist.

Wie in Fig. 4 gezeigt, unter besonderem Augenmerk auf die ΔABC(1) und ΔDBC(2), falls diese Dreiecke in Richtung der Normalen um einen Abstand verschoben werden, werden die vorhergehenden gestreunt und ein Zwischenraum oder eine Duplikation wird zwischen den zwei Dreiecken erzeugt. Der Zwischenraum oder Duplikation ist der Zuwachs der Oberflächen. Falls kein Zwischenraum oder Duplikation auftritt, ist der Reziprokwert des Krümmungsradius gleich Null. Für den Fall eines Zwischenraums ist das Flächeninkrement positiv. Für den Fall einer Duplikation ist jenes negativ.

Entsprechend dem obigen Flächeninkrement dA ist das Volumeninkrement dV als die Fläche der Dreiecke (ΔG1BC und D G2BC) multipliziert mit ε gegeben, wobei der Punkt G1 der Massenmittelpunkt (oder Innenmittelpunkt) des Dreiecks 1 ist und der Punkt G2 jener des Dreiecks 2.

Somit wird die Krümmung einer Oberflächenzelle wie folgt ausgedrückt:

dA = Σi < j bÿ · |ni - nj|
dV = Σi < j Uÿ · ε
dA/dV = (Σi < j bÿ · |ni - nj|]/(Σi < j Uÿ),

wobei angenommen wird, daß ni und nj die normalen Einheitsvektoren der Dreiecke (i und j) sind, b die Länge der Seite (BC) und U die Fläche der Dreiecke (ΔG1BC + ΔG2BC).

Unter der Annahme, daß C (x, y, z) eine Material-Existenz bezeichnet, daß die Materialoberfläche stetig ist und daß C(x,y,z) zweimal bezüglich x,y,z differenzierbar ist, wird C(x,y,z) wie folgt entwickelt:

C(r) = C(r₀) + ΣiC/ri|(r = r0)(ri - ri0)
+ 1/2 [ΣiC/ri|(r = r0)(ri - ri0)2]
+ 2 · Σi < j [C/(ri rj)|(r = r0)(ri - ri0)(rj - rj0)]
+ . . .
r = (x,y,z), r0 = (y0, y0, 70).

Der erste Term auf der rechten Seite dieses Ausdruck bedeutet die Position der Existenz von Material, und der zweite Term stellt den Normalen-Vektor auf den Oberflächen von Material-Topographie und den Gradient der Material-Existenzrate dar. Der dritte und vierte Term bezeichnen den Gradient (Raum r) von Grad (C) und beziehen sich auf die Krümmung der Material-Oberfläche.

In dem herkömmlichen Verfahren teilt sich die 3D-Krümmung in zwei Hauptkrümmungsradien (R1, R2), die senkrecht zueinander sind, wie folgt:

dA/dV = 1/R1 + 1/R2.

Jedoch kann das herkömmliche Verfahren nicht selber die allgemeine mit Dreiecken verbundene Topographiebeschreibung berechnen, weil die 3D-Krümmung in zwei annähernd senkrecht zueinander stehenden Radien unterteilt werden muß. Somit ist die hier vorgeschlagene Definition von Krümmung nützlich, wenn es schwierig ist, in den Senkrechten zu unterteilen.

Als nächstes wird im Schritt 4 der Oberflächen-Volumenfluß entsprechend der Differenz des chemischen Potentials zwischen einer Oberflächenzelle und einer nächsten Oberflächenzelle berechnet.

Durch die Nernst & Einstein-Beziehung ist die Durchschnittsgeschwindigkeit von Oberflächen-Teilchen gegeben durch

V = - {D/(k · T)} · grad (µ).

Die Volumen-Flußdichte kann wie folgt ausgedrückt werden:

J = Ω · R · v
= [D′/(k · T)] · grad (dA/dV)
D′/(k · T) = D · Γ · Ω² · R/(k · T) |cm⁴/Sek.|

O: Oberflächendichte von Teilchen auf der Oberfläche.

Wie in Fig. 5 dargestellt, fließen die BPSG-Moleküle von Positionen eines großen chemischen Potentials zu jenen von kleinen. Für die praktische Berechnung, wie in Fig. 4 gezeigt, arbeitet der Massentransport-Algorithmus derart, daß sich das Volumen von Oberflächenteilchen von einer Zelle zu der benachbarten Zelle über die Kontaktlinien (fett-gedruckte Linien von Fig. 6) von Zellrand und Kontur-Oberfläche bewegt. Der Gradient von U wird durch die Differenz zwischen dem chemischen Potential von Massenzentrum von der Kontur-Oberfläche in einer Zelle und demjenigen in der benachbarten Zelle geleitet.

Danach wird der Material-Volumenfluß während einer winzigen Zeit Δt im Schritt S5 berechnet, um die Materialdichte zu zeichnen. Und die Kontur-Oberfläche wird durch Linear-Interpolation der Materialdichte berechnet. Die Kontur-Oberfläche der festen Konzentration wird die neue Oberfläche. Derart wird ein nächstes Δt gesetzt, wenn die Oberfläche nach Δt Minuten simuliert ist. Danach kehrt der Schritt von S7 nach S2 zurück, und es wird dieselbe Simulation, wie oben beschrieben, ausgeführt. Wenn die Dispositionszeit des Rückfluß-Prozesses t0 vom Schritt S7 vorüber ist, ist diese Simulation beendet.

In diesem erfindungsgemäßen Verfahren wird das chemische Potential durch Oberflächenenergie aus dem sich ändernden Material-Volumen und Oberflächen ohne den Krümmungsradius berechnet. Daher ist es möglich, eine vollständige dreidimensionale Simulation durchzuführen. Diese Simulation ist zweckdienlich für die Entwicklung von ULSI.

Das erste Experiment wird zur Anpassung des Oberflächen-Diffusions-Koeffizienten durchgeführt: D′/kT (t bezeichnet Zeit/Sek.), das zweite wird zur Auswertung des 3D-Rückflußmodells durchgeführt.

Erstes Experiment: Eine Stufe von 1 µm wird auf einem Siliziumsubstrat gebildet, und ein BPSG (Borphosphatsilikat-Glas) Film (1 µm) eines TEOS (Tetraethylorthosilikat-Glas) System wird auf dem Substrat bei 375°C durch APCVD (chemische Verdampfungsablagerung bei Atmosphärendruck) abgelagert. Gleichzeitig werden BPSG-Filme mit einigen verschiedenen Konzentrationen von Borglas (B203) und Phosphorglas (P205) abgelagert. Im Vergleich mit BPSG von SiH4-O2-System wird ein TEOS BPSG-Film einer isotropischen Ablagerung wie in der SEM-Aufnahme aus Fig. 7 gezeigt, und dann wird die Oberfläche gleichmäßiger wie dasjenige eines SiH4-O2-Systems.

Die abgelagerten BPSG-Filme werden zurückgeflossen unter der Bedingung, daß die Atmosphäre N2 ist, die Temperatur 850°C und die Rückflußzeit 30 Minuten.

Fig. 8 zeigt die Beziehung des maximalen Winkels der Oberfläche des BPSG Films und die Gesamtkonzentration (Mol-%) von B203 und P205, entsprechend dem obigen Experiment. Die Referenzwerte werden als Beispiel in Fig. 8 gezeigt. Jedoch ist es unmöglich, die Werte der Erfinder direkt zu vergleichen, wegen des Temperaturunterschiedes und der anfänglichen Topographie. Jedoch ist die Abweichung in unserem Ergebnis hoch und die maximalen Winkel hängen nicht von der einzigen Gesamtkonzentration ab, wurde eine Anpassung des Oberflächen-Diffusionskoeffizienten (D′/kT) durchgeführt.

Als Ergebnis der Simulation zeigt Fig. 9 die Beziehung des Oberflächen-Diffusionskoeffizienten und die maximalen Winkel. Diese Beziehung kann als eine Exponentialfunktion beschrieben werden. Unter Bezugnahme auf die Fig. 8 und 9 kann der Oberflächen-Diffusionskoeffizient wie folgt abgeleitet werden:

D′ · t/(k · T) × 10 - 28 · exp (h - Cimp) [cm⁴]
h = 1,6 [l/(Mol.-%)]

Fig. 6 zeigt einige Vergleiche von Simulationen und Experimenten und die Topographien werden ausgewertet. Die Simulation von Ablagerung wird als vollkommen isotrop angenommen. Wie in Fig. 6(a) gezeigt stimmen die Simulation und das Experiment gut überein. Isotropische Ablagerungen weisen auf eine starke Oberflächen-Diffusion oder auf eine perfekte Steuerung durch chemische Reaktion hin. Für die Rückflußsimulation zeigen Fig. 6(b) und 6(c) zwei Fälle der Oberflächen-Diffusionskoeffizienten (D′t/kT). Die Simulationen von D′t/kT: 9,6 × 10^-19 (cm⁴/Sek.) und 1.8 × 10^-17 (cm⁴/Sek.) entsprechen dem Experiment von Cimp: 16 Mol.-% bzw. 17,6 Mol.-%. Diese Topographien von Simulationen und Experimenten sind in guter Übereinstimmung.

Als nächstes wird unter Verwendung dieser abgeleiteten Koeffizienten der Rückfluß in drei Dimensionen simuliert. Das Experiment ist folgendes, nämlich daß zunächst BPSG auf einem ebenen Siliziumsubstrat abgelagert wird, als zweites das Loch eines Quadrats von 2.5 um durch Ätzen von 0,5 um Dicke gebildet wird und als letztes der BPSG-Film in der N2-Atmosphäre bei 850°C für 10 Minuten zurückgeflossen wird. Die Ergebnisse von Simulation und Experiment sind in Fig. 10 bzw. Fig. 11 gezeigt. In diesen Querschnittsansichten ist die Oberflächen-Topographie planarisiert, die Kantenform des Quadrats wird rund und die Seitenmauer wird schräg.

Claims

Topographie-Simulationsverfahren, gekennzeichnet durch folgende Schritte:
Unterteilen des Analysefeldes in einige Zellen;
Zeichnen der Iso-Konzentrationsfläche, die durch einen Satz von einigen Polygonen von der Material-Konzentration von jeder Zelle beschrieben wird;
Zeichnen des chemischen Potentials der Zelle aus dem Verhältnis des Material-Volumenzuwachses zum Material-Oberflächen-Flächenzuwachs, wenn die Polygone in der Zelle, die die Iso-Konzentrationsfläche aufweist, sich in der Normalrichtung von jedem Polygon bewegen;
Zeichnen der Material-Stromdichte aus der Differenz von jedem chemischen Potential der Zelle;
Zeichnen der Material-Konzentration der Zelle aus der Material-Stromdichte; und
Zeichnen einer neuen Iso-Konzentrationsfläche aus der Material-Konzentration.