(57)【特許請求の範囲】
【請求項１】
半導体装置の内部及び周囲に有限個の格子点を設置し、複数の格子を生成するステップと、
該格子の内部のそれぞれに、線形に変化するキャリアの擬フェルミポテンシャルを設定するステップと、
該擬フェルミポテンシャルから、それぞれの格子内の任意の点におけるキャリア濃度を求めるステップと、
該キャリア濃度を用いて、ポアッソンの方程式及びキャリアの連続の方程式を、所定の境界条件の下で解き、半導体装置の電気的特性を得るステップ
とを少なくとも有するシミュレーション方法。
【請求項２】
Ｑをキャリアの電荷、ｋをボルツマン定数、Ｔを温度、τpをキャリアの運動量緩和時間、ｖthをキャリアの熱速度、φを擬フェルミポテンシャル、Ｅrを位置ベクトルｒにおける電界ベクトル、φＥrを擬フェルミポテンシャルφを考慮した位置ベクトルｒにおける電界ベクトル、UＥrを位置ベクトルｒにおける単位電界ベクトル、UφＥrを擬フェルミポテンシャルφrを考慮した位置ベクトルｒにおける単位電界ベクトル、Ｃ（ｒ）を位置ベクトルｒにおけるキャリア濃度、Λ（ｒ）をキャリア等濃度面に垂直方向のドリフト速度ベクトルを与える係数、ξ（ｒ）をキャリア等濃度面に平行方向のドリフト速度を与える係数としたとき、位置ベクトルｒにおけるキャリア電流密度Ｊ（ｒ）が、
【数１】
0003660137.tif 000005
で与えられることを用いて、前記キャリアの連続の方程式を解くことを特徴とする請求項１記載のシミュレーション方法。
【請求項３】
Ｑをキャリアの電荷、τpをキャリアの運動量緩和時間、Ｃ（ｒ）を位置ベクトルｒにおけるキャリア濃度、Λ（ｒ）をキャリア等濃度面に垂直方向のドリフト速度ベクトルを与える係数としたとき、位置ベクトルｒにおけるキャリア電流密度ベクトルＪ（ｒ）の発散▽・Ｊ（ｒ）が（Ｑ／τp）Λ（ｒ）Ｃ（ｒ）で与えられることを用いて前記キャリアの連続の方程式を解くことを特徴とする請求項１記載のシミュレーション方法。
【請求項４】
半導体装置の内部及び周囲に有限個の格子点を設置し、複数の格子を生成する格子生成手段と、
該格子の内部のそれぞれに、線形に変化するキャリアの擬フェルミポテンシャルを設定する線形擬フェルミポテンシャル設定手段と、
半導体装置の所定の電極領域に印加するバイアスを設定するバイアス設定手段と、
前記擬フェルミポテンシャルから、それぞれの格子内の任意の点におけるキャリア濃度を求め、ポアッソンの方程式及びキャリアの連続の方程式に対する係数行列と残差ベクトルを設定する係数行列と残差ベクトル設定手段と、
前記ポアッソンの方程式及びキャリアの連続の方程式の解に、前記係数行列を解いて得られた修正値を加えて、解を修正し、該解が収束するまで計算を繰り返す行列計算手段
とを少なくとも有するシミュレータ。
【請求項５】
半導体装置の内部及び周囲に有限個の格子点を設置し、複数の格子を生成するステップと、
該格子の内部のそれぞれに、線形に変化するキャリアの擬フェルミポテンシャルを設定するステップと、
該擬フェルミポテンシャルから、それぞれの格子内の任意の点におけるキャリア濃度を求めるステップと、
該キャリア濃度を用いて、ポアッソンの方程式及びキャリアの連続の方程式を、所定の境界条件の下で解き、半導体装置の電気的特性を得るステップ
とを少なくとも有するシミュレーションプログラムを記録した記録媒体。
【請求項６】
製造しようとする半導体装置の内部及び周囲に有限個の格子点を設置し、複数の格子を生成し、該格子の内部のそれぞれに、線形に変化するキャリアの擬フェルミポテンシャルを設定し、該擬フェルミポテンシャルから、それぞれの格子内の任意の点におけるキャリア濃度を求め、該キャリア濃度を用いて、ポアッソンの方程式及びキャリアの連続の方程式を、所定の境界条件の下で解き、半導体装置の電気的特性を得るステップと、
該電気的特性が所望の半導体装置の電気的特性であるか否か、評価・検討するステップと、
該電気的特性が所望の半導体装置の電気的特性であると判断された場合には、半導体材料に対し、所定の不純物ドーピング工程及び形状加工工程を含む一連の半導体製造プロセスを実行し、半導体装置を得るステップと
を少なくとも有する半導体装置の製造方法。
【発明の詳細な説明】
【０００１】
【発明の属する技術分野】
本発明はＬＳＩ等の半導体集積回路、高周波トランジスタ、電力用半導体装置、発光素子、受光素子等の種々の半導体装置のシミュレーション技術及びこれを用いた半導体装置の製造方法に関し、特に半導体装置のデバイスシミュレーション方法、デバイスシミュレータ、デバイスシミュレーションプログラムを記録した機械読みとり可能な記録媒体およびこれらを用いた半導体装置の製造方法に係る。
【０００２】
【従来の技術】
昨今の複雑高度化した半導体装置においては、半導体装置の特性評価を、実物を作成して、その実物の特性を測定し評価を行うのでは、膨大な時間の無駄と、膨大な費用の無駄を発生させる恐れがあり、好ましくない。このため、微細化され高集積密度化されたＬＳＩ等の半導体装置の製造工程においては、製造しようとする半導体装置の物理的、電気的な挙動を、これとほぼ同じ法則に支配されるコンピュータシステムの挙動によって模擬（シミュレーション）し、模擬的な半導体装置に対して特性評価を行い、特性を確認してから現実の製造工程を開始するのが必須となりつつある。
【０００３】
しかし、コンピュータシステム等でのシミュレーションによる特性評価の結果に相当量の計算誤差が含まれていれば、当然、これを用いた半導体装置の製造工程に不具合を発生させる。例えば、数週間程度以上を必要とする長く、しかも複雑な製造工程を経て試作品が完成された後に、その特性等に不具合が発生すれば、再び設計をし直し、再び、長時間かつ複雑な製造工程を繰り返さなければならず、時間とランニングコストが無駄に消費されてしまう。このように、シミュレーションの精度の向上は、直接半導体装置の製造の期間の短縮化にもつながるため非常に重要である。もとより、設計の失敗によって失われた時間は回復することはできないのである。特に半導体産業においては、より高性能の半導体装置の開発が競われているが、最も重要なのはこの半導体装置の開発の速度であるため、設計、開発期間の短縮は極めて重要である。
【０００４】
この半導体装置の特性評価段階で用いられる装置（デバイス・シミュレータ）は半導体装置の形状、半導体中の不純物分布、端子に印加する電圧などの入力データを受け取る。そして、このデバイス・シミュレータは半導体の電気特性の基本方程式であるポアソン方程式と電子電流連続式と正孔電流連続式との解を数値計算により求め、さらに、半導体装置の電気的特性を出力する。
【０００５】
半導体装置の電気的特性を得るために必要な基本方程式は、以下の（１）式乃至（３）式に示すような電位分布、キャリア濃度（具体的には電子濃度ｎ、正孔濃度ｐ）の連立微分方程式である。
【０００６】
▽・Ｄ＝ｑ（ｐ−ｎ＋ＮD−ＮA）                         …(１)
∂（−ｑｎ）／∂ｔ＋▽・Ｊe＝−ｑＵ                    …(２）
∂（ｑｐ）／∂ｔ＋▽・Ｊn＝ｑＵ                        …(３)
ここで、Ｄは電束密度ベクトル、ｑは素電荷、ＮDはドナー濃度、ＮAはアクセプター濃度、ｔは時間、Ｊeは電子電流密度ベクトル、Ｊnは正孔電流密度ベクトル、Ｕはキャリア生成再結合率である。この基本方程式に境界条件として、半導体装置の各点における電位Ψ、電子濃度ｎ、正孔濃度ｐを任意に与え、それらの値が上記基本方程式を満たすまで修正を加えて、コンピュータシステム等で解を求める。しかし、（１）式乃至（３）式に示す基本方程式をコンピュータシステムで解くためには、キャリア電流密度ベクトルの離散化を行う必要がある。
【０００７】
そこで、図７に示すように、デバイスシミュレーションでは、まず、シミュレーションの対象となる半導体装置の内部及び周囲（解析領域と呼ぶ。）に有限個の格子点を設置する。そして、２次元空間の場合は、各格子点を頂点とする多角形を設置する（３次元空間の場合には、各格子点を頂点とする多面体を設置する）。この多角形は格子と呼ばれる。これらの多角形または多面体は隙間なく、かつ重なりなく解析領域を覆う。格子点には必要に応じて物理量を対応させる。また、それぞれの格子点に対して、この格子点を内部に含む、コントロールボリュームと呼ばれる多角形または多面体を対応させる。コントロールボリュームも隙間なく、かつ重なりなく解析領域を覆う。キャリア電流密度ベクトルの離散化とは、格子点の物理量を基にしてコントロールボリュームから流出する（或いは流入する）キャリア電流密度ベクトルを算出する技術のことである。
【０００８】
このキャリア電流密度ベクトルの離散化に関して、２つの従来技術が知られている。
【０００９】
１つは、シークフリート・セルベルハー（Siegfried Selberherr）による「半導体装置の解析とシミュレーション」（“Analysis and Simulation of Semiconductor Devices”）、スプリンガー・フェアラーク（Springer-Verlag Wien New York）社、1984年発行.の6.3節：有限要素（ Finite Elements）,第181頁乃至第191頁に示されている技術である。以下、この技術を「第１の従来技術」と呼ぶ。
【００１０】
第１の従来技術は、１次元の線分上を流れるキャリア電流密度を計算する方法として知られているＳＧ(Scharfetter-Gummel)法を、２次元の平面に設定した三角形または四角形の格子系に対して応用した方法である。すなわち、第１の従来技術では、格子系のそれぞれの線分（三角形または四角形の辺）を流れるキャリア電流密度を、キャリア電流は格子系のそれぞれの線分上のみを流れると仮定してＳＧ法で計算する。この２次元の格子系は、同じ形で、導線をつなぎあわせた回路のようなものである。第１の従来技術においては、この回路の結び目と結び目の間をつなぐ導線を流れる電流を計算する。この計算は、以下のようである。
【００１１】
まず、格子点Ｉに対応するコントロールボリュームから流出する電流SelＩO,Iを次式で与える。ここで、図８と図９に示すような記号の定義を用いる。
【００１２】
SelＩO,I＝Σ［Ｊ＝Ｏ，ＮI］SGＩIJ                      …(４)
（４）式において、Σ［Ｊ＝Ｏ，ＮI］はSGＩIJ  のＪ＝０からＪ＝ＮIまでの和を意味する（本発明では、図１２に定義するような、一般的な数学的表現とは異なる表記を便宜上採用していることに留意されたい）。
【００１３】
SGＩIJは格子点Ｉと格子点Ｊを端点とする辺に対して算出される、ＩからＪへ向かう電流で、次式で定義される。
【００１４】
【数２】
SGＩIJ＝（ｄIJ・ＤIJ／ｌIJ）（ＢIJＣI−Ｂ*IJＣJ）      …(５)
ＢIJ＝Ｂ（ＱΨIJ／ｋＴ）                               …(６)
Ｂ*IJ＝Ｂ（−ＱΨIJ／ｋＴ）                            …(７)
ここで、ｋはボルツマン定数、Ｔは温度、ＣI，ＣJはそれぞれ格子点Ｉ，Ｊにおけるキャリア濃度、ΨIJはそれぞれ格子点Ｉ，Ｊにおける電位差（ΨIJ＝ΨI−ΨJ）、ｄIJは三角形の格子の外接円の中心と辺ＩＪとの距離、ＤIJはＩからＪへ向かう拡散係数、ｌIJは格子点ＩＪ間の距離である。また、Ｑはキャリアの電荷で、電子に対しては−ｑ、正孔（ホール）に対しては＋ｑである。
【００１５】
第１の従来技術では電流は格子点を結ぶ辺を流れるとしてＳＧ法で計算している。しかし、それらの辺はシミュレーションを行うために仮想的に設置したものである。現実の半導体装置では、電流は半導体装置を構成している半導体領域の任意の点を流れる。このシミュレーションと現実の電流の流れ方の差異から、第１の従来技術の問題点が生じる。それは、格子点の設置の仕方に依存して計算誤差が発生するという問題である。この計算誤差は「擬拡散」と呼ばれる。擬拡散が発生するとシミュレーション結果が正しくなくなるので、擬拡散は重大な問題である。
【００１６】
擬拡散の問題は多くの研究者や技術者によって指摘されている。例えば、アイ・イー・イー・イー・トランズアクション・オン・コンピューター・エイデッド・デザイン・オブ・インテグレイテッドサーキッツ・アンド・システムズ(IEEE Transactions on Computer-Aided Design of Integrated Circuits and Systems)第１２巻（ 1993年）、第１０号、の第1535頁乃至第1541頁のヒー・イーとテン・ティーメン（He Yie and Teng Zhimeng）による「ドリフト-拡散方程式に対する流れに逆らった非振動の定式化」（“Nonosillatory Streamline Upwind Formulatioins for Drift-Diffusion Equation”）に示されている（以下、この文献に記載された技術を「第２の従来技術」と呼ぶ）。第２の従来技術の目的は第１の従来技術において発生する擬拡散を緩和することである。この第２の従来技術においては、キャリア電流は格子の線上のみではなく、任意の方向に流れると仮定する。キャリア電流は、一般に拡散（ディフュージョン）電流と移流（ドリフト）電流の和で与えられる。そして、第２の従来技術では、移流電流の方向に格子とは別の仮想的な線分を考え、この仮想的な線分に付いてＳＧ法を使って電流密度を計算する。しかし、ＳＧ法を使うためには線分の両端のキャリア濃度が分かっていることが必要であり、仮想的な線分の両端のキャリア濃度は実在する格子点のキャリア濃度から内挿して求める。このようにしてＳＧ法を用いて算出したキャリア電流密度は、２次元平面におけるキャリア電流密度ベクトルの仮想線分方向の成分であると仮定している。これとは別に、第２の従来技術では、仮想的な線分とは垂直な方向に、濃度勾配に比例した拡散電流が流れると仮定する。そして、この二つの電流のベクトル和が第２の従来技術によって算出されるキャリア電流密度ベクトルとなる。
【００１７】
具体的には、第２の従来技術においては、擬拡散を緩和するために図１０に示すように、キャリア濃度の補間を行っている。格子点（ｉ，ｊ）における電流密度ベクトルは、第１の従来技術のように周囲の格子点のキャリア密度のみから算出するのではなく、格子点（ｉ，ｊ）の上流側のキャリア濃度ｎUを基にしてキャリア電流密度ベクトルを計算する。ｎUは周囲の格子点（ｉ−１，ｊ）と（ｉ−１，ｊ−１）のキャリア濃度から算出する。即ち、３つの格子点のキャリア濃度からキャリア電流密度ベクトルを算出する。
【００１８】
しかし、第２の従来技術でも擬拡散をなくすことはできない。図１１は第２の従来技術から引用したもので、第１の従来技術と第２の従来技術による計算誤差を示す図である。第２の従来技術でも擬拡散が生じていることが分かる。
【００１９】
【発明が解決しようとする課題】
１次元の線分上を流れるキャリア電流密度を計算する方法であるＳＧ法を、２次元乃至３次元のデバイスシミュレーションに適用する際の問題点は、算出されるデバイス特性やデバイスシミュレーションの収束性が格子の設定方法に依存して変化してしまうことである。デバイスシミュレーションにおいて生成する格子は、実際の半導体装置にはない、計算の便宜上のものである。理想的には、格子の設定方法に依存せずに計算結果が得られることが、デバイスシミュレーションとしては望ましいということになる。
【００２０】
この格子の設定方法に対するキャリア電流密度ベクトルの依存性は、ＳＧ法が１次元の辺上の電流密度ベクトルしか算出できないことから生じている。電流が２次元や３次元空間に分布するデバイスシミュレーションをＳＧ法で実行するためには、格子の辺上しか電流が流れないと仮定（第１の従来技術）、若しくは、実在する格子点のキャリア濃度から内挿した、格子とは別の仮想的な線分を考え、この仮想的な線分についてＳＧ法を使って電流密度を計算する。しかし、この仮定のために、算出されるデバイス特性やデバイスシミュレーションの収束性が格子の設定方法に強く依存するという問題が生じる。
【００２１】
このように、第１の従来技術や第２の従来技術においては、格子内の任意の点におけるキャリア電流密度ベクトルを与えることができなかった。このため、格子の辺上を流れるキャリア電流密度ベクトルを基に、近似的にキャリア電流連続式を解いていた。このため、擬拡散の問題が発生し、、基本方程式の解を正しく算出できず、半導体装置の電気特性を不正確にしかシミュレーションできなかった。
【００２２】
上記問題点を鑑み、本発明の目的は、格子内の任意の点におけるキャリア電流密度ベクトルを正確に算出し、(１)式乃至（３）式で示された基本方程式の解を正確に算出することが出来るシミュレーション方法を提供することである。
【００２３】
本発明の他の目的は、格子内の任意の点におけるキャリア電流密度の発散を正確に算出し、基本方程式の解を正確に算出することが出来るシミュレーション方法を提供することである。
【００２４】
本発明のさらに他の目的は、格子内の任意の点におけるキャリア電流密度ベクトルを正確に算出し、基本方程式の解を正確に算出することが出来るシミュレータを提供することである。
【００２５】
本発明のさらに他の目的は、格子内の任意の点におけるキャリア電流密度の発散を正確に算出し、基本方程式の解を正確に算出することが出来るシミュレータを提供することである。
【００２６】
本発明のさらに他の目的は、格子内の任意の点におけるキャリア電流密度ベクトルを正確に算出し、基本方程式の解を正確に算出することが出来るシミュレーションプログラムを記録した記録媒体を提供することである。
【００２７】
本発明のさらに他の目的は、格子内の任意の点におけるキャリア電流密度の発散を正確に算出し、基本方程式の解を正確に算出することが出来るシミュレーションプログラムを記録した記録媒体を提供することである。
【００２８】
本発明のさらに他の目的は、試作に先立つ半導体装置の特性予測を正確に行うことにより、所望の特性の半導体装置を短期間で開発できる半導体装置の製造方法を提供することである。
【００２９】
【課題を解決するための手段】
上記目的を達成するため、本発明の第１の特徴は、半導体装置の内部及び周囲に有限個の格子点を設置し、複数の格子を生成するステップと、この格子の内部のそれぞれに、線形に変化するキャリアの擬フェルミポテンシャルを設定するステップと、この擬フェルミポテンシャルから、それぞれの格子内の任意の点におけるキャリア濃度を求めるステップと、このキャリア濃度を用いて、(１)式乃至（３）式で示されたポアッソンの方程式及びキャリアの連続の方程式を、所定の境界条件の下で解き、半導体装置の電気的特性を得るステップとを少なくとも有するシミュレーション方法であることである。ここで、「キャリア」とは、電子及び正孔（ホール）を意味し、「キャリアの擬フェルミポテンシャル」とは、電子の擬フェルミポテンシャルφn及び正孔の擬フェルミポテンシャルφpを意味する。
【００３０】
本発明の第１の特徴に係るシミュレーション方法によれば、生成した格子内の任意の点におけるキャリア電流密度ベクトルを、正確に算出できる。このため、(１)式乃至（３）式で示されたデバイスシミュレーションの基本方程式の解を正確に算出できる。つまり、本発明の第１の特徴に係るシミュレーション方法によれば、求める半導体装置の電気特性を正確にシミュレーションできる。
【００３１】
本発明の第１の特徴においては、半導体装置を構成する半導体領域中のキャリアは熱速度ｖthでブラウン運動をしていると仮定し、衝突のタイムスケールを運動量緩和時間τpで与える。そして、半導体領域内部の微小領域（格子）を想定して、微小領域（格子）の内部では、電界ベクトルＥと擬フェルミポテンシャルの勾配▽φn、▽φpは定数であることを仮定している。また、キャリア濃度はボルツマン統計に従うと仮定する。
【００３２】
従って、正孔の擬フェルミポテンシャルφp及び電子の擬フェルミポテンシャルφnを用いると、正孔濃度ｐ及び電子濃度ｎは次の(８)式及び（９）式で与えられる。
【００３３】
ｐ＝ｎiexp（−ｑ（Ψ−φp）／ｋＴ）                   …(８)
ｎ＝ｎiexp（ｑ（Ψ−φn）／ｋＴ）                     …(９)
(１)式乃至（３）式で示されたポアッソンの方程式及びキャリアの連続の方程式を、所定の境界条件の下で解くためには、ドナー濃度ＮD、アクセプター濃度ＮA、キャリア生成再結合率Ｕについても、生成された格子の内部のそれぞれにおいて、線形若しくは指数関数的に変化する等の所定の仮定をすることが好ましい。  (１)式乃至（３）式で示された基本方程式は連立偏微分方程式であるが、電子計算機によって解くためにはそれらを連立一次方程式に変換しなければならない。
【００３４】
詳細な計算手順は後述するが、本発明の第１の特徴においては、積分方程式を代数方程式に変換して、位置ベクトルｒにおけるキャリア電流密度ベクトルＪ（ｒ）を、キャリアの運動量緩和時間τp、キャリアの熱速度ｖth、擬フェルミポテンシャルφrを用いて
【数３】
Ｊ（ｒ）＝（ｖdv（ｒ）＋ｖdp（ｒ））ＱＣ（ｒ）         …(10)
ｖdv（ｒ）＝（ｋＴ／Ｑτp）Λ（ｒ）UφＥr／｜φＥr｜     …(11)
ｖdp（ｒ）＝ｖthξ（ｒ）（UＥr−（UＥr・UφＥr）UφＥr）  …(12)
で計算することが好ましい。(11)式のΛ（ｒ）はキャリア等濃度面に垂直方向のドリフト速度ベクトルを与える係数で、
【数４】
Λ（ｒ）＝１−｜Ｅr／２｜sinｈ（ζrφ）／（ηφ・sinｈ（ζr））…(13)
と表わされる。ここで２ζrは熱電位（ｋＴ／ｑ）で規格化した位置ベクトルｒにおけるポテンシャルで、
２ζr＝Ｑ｜Ｅr｜ ｜ｌ｜／ｋＴ                          …(14)
と表わされる。(13)式におけるηφは次式で与えられる。
【００３５】
ηφ＝｜Ｅr／２＋▽φr｜                                …(15)
(13)式におけるζrφは熱電位で規格化した位置ベクトルｒにおけるηφに対応したポテンシャルで、
ζrφ＝Ｑηφ｜ｌ｜／ｋＴ                               …(16)
と表わされる。一方、(12)式のξ（ｒ）はキャリア等濃度面に平行方向のドリフト速度を与える係数で、以下の式で示される。
【００３６】
ξ（ｒ）＝coth（ζr）−１／ζr                        …(17)
さらに(11)式におけるφＥrは擬フェルミポテンシャルφrを考慮した位置ベクトルｒにおける電界ベクトルで
φＥr＝Ｅ＋▽φr                                       …(18)
と表わされる。ここで、Ｅは位置ベクトルｒにおける電界ベクトルである。またUＥrは位置ベクトルｒにおける単位電界ベクトル、UφＥrは位置ベクトルｒにおける擬フェルミポテンシャルφrを考慮した単位電界ベクトルで、それぞれ以下のように表わされる。
【００３７】
UＥr＝Ｅr／｜Ｅr｜                                    …(19)
UφＥr＝φＥr／｜φＥr｜                                 …(20)
また、本発明の第１の特徴において、位置ベクトルｒにおけるキャリア電流密度ベクトルＪ（ｒ）の発散▽・Ｊ（ｒ）が
▽・Ｊ（ｒ）＝（Ｑ／τp）Λ（ｒ）Ｃ（ｒ）             …(21)
で与えられることを用いて半導体装置の電気的特性を得るようにしてもよい。  上述のように、本発明の第１の特徴においては、半導体装置を構成する半導体領域中のキャリアは熱速度ｖthでブラウン運動をしていると仮定している。このような熱速度ｖthでブラウン運動の下では、ある時刻ｔ0にある点ｘにあったキャリアはτp時間後の時刻ｔ1＝ｔ0＋τpには平均自由行程ｌ＝ｖth・τpだけ離れた点にあることになる。３次元空間で考えると、もしも電界が０であれば、キャリアは球面上に均一に分布するはずである。しかし、一般には電界が０でないので、キャリアは電界の方向に偏って分布する。本発明の第１の特徴においては、ボルツマン統計を仮定しているので、この偏った分布は電界Ｅとキャリア温度Ｔと球面の半径である平均自由行程ｌの関数として求まる。τpの間のキャリアの変位を平均自由行程ベクトルｌで表すと、キャリアが平均自由行程ベクトルｌの方向に変位する単位時間単位立体角当たりの確率Ｐ（ｌ）は次式で与えられる。
【００３８】
【数５】
Ｐ（ｌ）＝１／（４πτp）・ζr・exp（ζr）／sinｈζr      …(22)
平均自由行程ベクトルｌよりもずっと小さな微小領域を点ｘの周りに取ると、ｔ0に同領域にあったキャリアはｔ1にはすべて外に出ている。逆にｔ1に同領域にいるキャリアは元のｔ0にはｘの周りの半径ｌの球面上にあったはずである。球面上のキャリア濃度は次式で与えられる。
【００３９】
【数６】
Ｃ（ｌ）＝Ｃ（０）exp（Ｑ（Ｅ＋▽φ）・ｌ／ｋＴ）          …(23)
それらのキャリアがτp後に点ｘにある確率はＰ（−ｌ）に比例するので、それらの積を積分すると電流密度の発散が求まる。
【００４０】
▽・Ｊ＝（Ｑ／τp）Ｃ（０）Λ（ｒ）                      …(24)
(10)式の電流密度ベクトルＪ（ｒ）は(24)式を利用して求められる。
【００４１】
本発明の第２の特徴は、格子生成手段、線形擬フェルミポテンシャル設定手段、バイアス設定手段，係数行列と残差ベクトル設定手段及び行列計算手段とを少なくとも有するシミュレータであることである。ここで、格子生成手段は、半導体装置の内部及び周囲に有限個の格子点を設置し、複数の格子を生成する。線形擬フェルミポテンシャル設定手段は、生成された格子の内部のそれぞれに、線形に変化するキャリアの擬フェルミポテンシャルを設定する。バイアス設定手段は、ソース電極領域，ドレイン電極領域、ゲート電極領域等の半導体装置の所定の電極領域に印加するバイアスを、シミュレーションする電気特性に応じて入力する。係数行列と残差ベクトル設定手段は、擬フェルミポテンシャルから、それぞれの格子内の任意の点におけるキャリア濃度を求め、基本方程式に対する係数行列と残差ベクトルを設定する。行列計算手段は、周知のニュートン法によって解を求める。すなわち、基本方程式の解に、係数行列を解いて得られた修正値を加えて、解を修正し、基本方程式の解が収束するまで計算を繰り返す。そして、解が収束していれば、所定のバイアスの解としてその解を保存する。このようにして、各バイアスについて、バイアス設定手段，係数行列と残差ベクトル設定手段及び行列計算手段とを用いて、ポアッソンの方程式及びキャリアの連続の方程式から変換された連立一次方程式を、所定の境界条件の下で解き、計算誤差の少ない半導体装置の電気的特性を短時間に得ることができる。
【００４２】
本発明の第２の特徴に係るシミュレータによれば、生成した格子内の任意の点におけるキャリア電流密度ベクトルを、正確に算出できる。このため、(１)式乃至（３）式で示された半導体デバイスシミュレーションの基本方程式の解を正確に算出できる。つまり、本発明の第２の特徴に係るシミュレータによれば、求める半導体装置の電気特性を正確にシミュレーションできる。
【００４３】
本発明の第１の特徴で述べたシミュレーション方法を実現するためのプログラムは、コンピュータ読取り可能な記録媒体に保存し、この記録媒体をコンピュータシステムによって読み込ませることにより、本発明のシミュレーション方法を実行することができる。すなわち、本発明の第３の特徴は、半導体装置の内部及び周囲に有限個の格子点を設置し、複数の格子を生成するステップと、この格子の内部のそれぞれに、線形に変化するキャリアの擬フェルミポテンシャルを設定するステップと、この擬フェルミポテンシャルから、それぞれの格子内の任意の点におけるキャリア濃度を求めるステップと、このキャリア濃度を用いて、ポアッソンの方程式及びキャリアの連続の方程式を、所定の境界条件の下で解き、半導体装置の電気的特性を得るステップとを少なくとも有するシミュレーションプログラムを記録した記録媒体であることである。本発明の第３の特徴によれば、半導体装置の電気特性を正確にシミュレーションできるプログラムを記録した記録媒体を提供できる。なぜなら、格子内の任意の点におけるキャリア電流密度ベクトルを正確に算出できるプログラムであるため、キャリア電流連続式を含む半導体デバイスシミュレーションの基本方程式の解を正確に算出できるからである。ここで、「記録媒体」とは、例えばコンピュータの外部メモリ装置、半導体メモリ、磁気ディスク、光ディスク、光磁気ディスク、磁気テープなどのプログラムを記録することができるような機械読みとり可能な記録媒体を意味する。具体的には、フロッピーディスク、ＣＤ−ＲＯＭ，ＭＯディスク、カセットテープ、オープンリールテープなどが「記録媒体」に含まれる。
【００４４】
本発明の第４の特徴は、製造しようとする半導体装置の内部及び周囲に有限個の格子点を設置し、複数の格子を生成し、格子の内部のそれぞれに、線形に変化するキャリアの擬フェルミポテンシャルを設定し、擬フェルミポテンシャルから、それぞれの格子内の任意の点におけるキャリア濃度を求め、キャリア濃度を用いて、ポアッソンの方程式及びキャリアの連続の方程式を、所定の境界条件の下で解き、半導体装置の電気的特性を得るステップと、この得られた電気的特性が所望の半導体装置の電気的特性であるか否か、評価・検討するステップと、得られた電気的特性が所望の半導体装置の電気的特性であると判断された場合には、半導体材料に対し、所定の不純物ドーピング工程及び形状加工工程を含む一連の半導体製造プロセスを実行し、半導体装置を得るステップとを少なくとも有する半導体装置の製造方法であることである。
【００４５】
本発明の第４の特徴によれば、所望の特性の半導体装置を従来よりも短期間で開発できる。なぜなら、本発明のデバイスシミュレーション方法により半導体デバイスの電気的特性を正確にシミュレーションし、より精度の高い半導体装置の評価を行うことができるからである。すなわち、本発明の第４の特徴によれば、試作に先立つ半導体装置の特性予測を従来よりも正確に行えるので、より少ない試作回数で半導体装置を開発でき、研究・開発の効率が大幅に改善される。
【００４６】
【発明の実施の形態】
以下、図面を参照して本発明の実施の形態を説明する。
【００４７】
（第１の実施の形態）
本発明の第１の実施の形態に係るデバイスシミュレータ（デバイスシミュレーション装置）は、プロセスシミュレーションで得られた素子構造と不純物分布等の結果を入力し、入力した構造に対して印加電圧、電流などの電気的な境界条件を設定して素子の電気的特性をシミュレーションする機能を有する。さらに、本発明の第１の実施の形態に係るデバイスシミュレータは、その結果得られる半導体装置内の電位・電界・電流分布や電子・正孔などのキャリア分布、あるいは電流−電圧特性などを出力装置上に表示する機能を有する。
【００４８】
すなわち図１（ａ）に示すように、本発明のデバイスシミュレータ７０は、操作者からのデータや命令などの入力を受け付ける入力部６１と、半導体装置の電気的特性をシミュレーションする処理制御部７１と、シミュレーション結果を出力する出力部６４と、半導体装置の特性の解析に必要な所定のデータなどを格納したデータ記憶部６２と、デバイスシミュレーションプログラムなどを格納したプログラム記憶部６３とから少なくとも構成されている。処理制御部７１は端子電圧または電流条件を設定する電圧／電流設定手段７２および素子特性計算手段７３とを少なくとも有している。
【００４９】
図１に示す素子特性計算手段７３は、(１)式乃至（３）式で示されたポアソン方程式、電子電流連続式と正孔電流連続式を連立一次方程式に変換して解き、電界効果トランジスタ（ＦＥＴ）、バイポーラトランジスタ（ＢＪＴ）、静電誘導トランジスタ（ＳＩＴ）、サイリスタなどの半導体装置中の電子や正孔等のふるまいを計算することができる。このため、図１（ｂ）に示すように、本発明の素子特性計算手段７３は、格子生成手段５１、線形擬フェルミポテンシャル設定手段５２、バイアス設定手段５３，係数行列と残差ベクトル設定手段及び行列計算手段を少なくとも具備している。ここで、格子生成手段５１は、対象とした半導体装置の内部及び周囲に有限個の格子点を設置し、複数の格子を生成する。線形擬フェルミポテンシャル設定手段５２は、格子生成手段５１により生成された格子の内部のそれぞれに、線形に変化するキャリアの擬フェルミポテンシャルを設定する。バイアス設定手段５３は、ソース電極領域，ドレイン電極領域、ゲート電極領域等の半導体装置の所定の電極領域に印加するバイアスを、シミュレーションする電気特性に応じて入力する。係数行列と残差ベクトル設定手段５４は、擬フェルミポテンシャルから、それぞれの格子内の任意の点におけるキャリア濃度を求め、(１)式乃至（３）式で示された基本方程式に対する係数行列と残差ベクトルを設定する。そして、行列計算手段５５は、周知の行列解法によって解を求める。すなわち、基本方程式の解に、行列を解いて得られた修正値を加えて、解を修正し、基本方程式の解が収束するまで計算を繰り返す。そして、解が収束していれば、所定のバイアスの解としてその解を保存する。
【００５０】
このようにして、各バイアスについて、バイアス設定手段５３，係数行列と残差ベクトル設定手段５４及び行列計算手段５５とを用いて、ポアッソンの方程式及びキャリアの連続の方程式から変換された連立一次方程式を、所定の境界条件の下で解き、求める半導体装置の電気的特性を得ることができる。
【００５１】
図１（ａ）において、入力部６１はキーボード、マウス、ライトペンまたはフロッピーディスク装置などで構成される。処理制御部７１、データ記憶部６２およびプログラム記憶部６３はＣＰＵ、及びこのＣＰＵに接続されたＲＯＭ、ＲＡＭ、磁気ディスクなどの記憶装置を含む通常のコンピュータシステムで構成される。また出力部６４はディスプレイ装置やプリンタ装置などにより構成されている。
【００５２】
なお、本発明のデバイスシミュレーション方法を実現するためのプログラムはコンピュータ読取り可能な記録媒体に保存しておいてもよい。この記録媒体をコンピュータシステムによって読み込ませ、プログラム記憶部６３に格納し、このプログラムを処理制御部７１で実行して本発明のデバイスシミュレーション方法を実現することもできる。ここで、記録媒体とは、例えばコンピュータの外部メモリ装置、磁気ディスク、光ディスク、光磁気ディスク、磁気テープなどのプログラムを記録することができるような記録媒体などが含まれる。
【００５３】
図２は本発明のデバイスシミュレータの外観を示す鳥瞰図である。このデバイスシミュレータ８０の本体は、フロッピーディスク装置（フロッピーディスクドライブ）８１および光ディスク装置（光ディスクドライブ）８２を具備している。フロッピーディスクドライブ８１に対してはフロッピーディスク８３を、また光ディスクドライブ８２に対してはＣＤ−ＲＯＭ８４をその挿入口から挿入し、所定の読み出し操作を行うことにより、これらの記録媒体に格納されたプログラムをシステム内にインストールすることができる。また、所定のドライブ装置を接続することにより、例えばゲームパック等に利用されている半導体メモリとしてのＲＯＭ８５や、磁気テープとしてのカセットテープ８６等を用いることもできる。
【００５４】
次に、図１に示す素子特性計算手段７３における計算手順について述べる。デバイスシミュレーションの基本方程式は、(１)式乃至（３）式で示されたポアソン方程式、電子電流連続式と正孔電流連続式からなる連立偏微分方程式である。これらの基本方程式は連立偏微分方程式であるけれども、電子計算機によって解くためにはそれらを連立一次方程式に変換しなければならない。この手順は以下の通りである。
【００５５】
(１)式乃至（３）式をそれぞれコントロールボリュームで積分する。なお以下の説明において“INTG”は図１２に定義するように積分記号を表わす。すなわちaINTGbｆ（ｘ）ｄｘは関数ｆ（ｘ）の区間（ａ，ｂ）における線積分（定積分）を意味し、INTGCvｆ（ν）ｄ（ν）は関数ｆ（ν）のコントロールボリュームＣＶにおける体積積分を意味するものとする。またINTGCv,iｆ（ν）ｄ（ν）はｉ番目のコントロールボリュームＣv,iにおける体積積分を意味するものとする。このような表記のもとで
【数７】
INTGCv,i▽・Ｄｄν＝ｑ（Ｐi−Ｎi＋ＮDi−ＮAi）         …(25)
∂（−ｑＮi）／∂ｔ＋INTGCv,i▽・Ｊeｄν＝−ｑＵi      …(26)
∂（ｑＰi）／∂ｔ＋INTGCv,i▽・Ｊnｄν＝ｑＵi          …(27)
ここで
Ｐi＝INTGCv,iｐｄν                                    …(28)
Ｎi＝INTGCv,iｎｄν                                    …(29)
は、正孔濃度ｐ及び電子濃度ｎをｉ番目のコントロールボリュームで積分した量である。ここで、正孔濃度ｐ及び電子濃度ｎは、前述の(８)式及び（９）式で与えられる。
【００５６】
本発明の第１の実施の形態では電位Ψと正孔の擬フェルミポテンシャルφp及び電子の擬フェルミポテンシャルφnは各格子内では座標に関して線形に変化することを仮定している。また、電位Ψは半導体領域及び絶縁体領域に属する各格子点において値を与え、正孔の擬フェルミポテンシャルφp及び電子の擬フェルミポテンシャルφnは半導体領域に属する各格子点において値を与える。各格子は三角形なので、Ψ、φp、φnとも、３つの頂点での値から各格子内の任意の点における値を算出できる。更に、(28)式及び(29)式の積分は解析的に実行できる。
【００５７】
また、ドナー濃度ＮD、アクセプター濃度ＮA及びキャリア生成再結合率Ｕについても、格子の内部のそれぞれにおいて、線形若しくは指数関数的に変化すると仮定すると、
ＮD,i＝INTGCv,iＮDｄν                                …(30)
ＮA,i＝INTGCv,iＮAｄν                                …(31)
Ｕi＝INTGCv,iＵｄν                                   …(32)
である。ガウスの定理を用いると(25)式〜(27)式は以下のように変換できる。
【００５８】
【数８】
INTGΩCv,iＤ・ｄＳ−ｑ（Ｐi−Ｎi＋ＮD,i−ＮA,i）＝０   …(33)
∂（−ｑＮi）／∂ｔ＋INTGΩCv,iＪe・ｄＳ＝−ｑＵi      …(34)
∂（ｑＰi）／∂ｔ＋INTGΩCv,iＪn・ｄＳ＝ｑＵi          …(35)
図６に示したコントロールボリュームに対して(33)式の左辺第１項は次式で与えられる。
【００５９】
【数９】
0003660137.tif 000002
ここでＤi（ｉ＝２〜６）は格子ｉ（ｉ＝２〜６）における電束密度ベクトルである。すなわち格子ｉにおける誘電率をεｉ、格子ｉにおける電界ベクトルをＥiとすれば
Ｄ2＝ε2・Ｅ2                                         …(37)
Ｄ3＝ε3・Ｅ3                                         …(38)
Ｄ4＝ε4・Ｅ4                                         …(39)
Ｄ5＝ε5・Ｅ5                                         …(40)
Ｄ6＝ε6・Ｅ6                                         …(41)
である。各格子内で電位Ψは座標に対して線形に変化することを仮定しているので、各格子の内部で電界は一定である。
【００６０】
Ｅi＝−▽Ψi                                          …(42)
たとえば、格子２の電界Ｅ2は格子点２での電位Ψ2、格子点１での電位Ψ1、格子点６での電位Ψ6から求めることができる。また(36)式におけるｄ6,2,2、ｄ2,3,3、ｄ3,4,4、ｄ4,5,5、ｄ5,6,6は以下の(43)式〜(47)式で与えられる。すなわち
【数１０】
ｄ6,2,2＝ｄM6,C2ｎ1,6＋ｄC2,M2ｎ1,2                   …(43)
ｄ2,3,3＝ｄM2,C3ｎ1,2＋ｄC3,M3ｎ1,3                   …(44)
ｄ3,4,4＝ｄM3,C4ｎ1,3＋ｄC4,M4ｎ1,4                   …(45)
ｄ4,5,5＝ｄM4,C5ｎ1,4＋ｄC5,M5ｎ1,5                   …(46)
ｄ5,6,6＝ｄM5,C6ｎ1,5＋ｄC6,M6ｎ1,6                   …(47)
ここでｄMk,Cl（ｋ＝２〜６；ｌ＝２〜６）は点Ｍｋから点Ｃｌへ向かう辺Ｍｋ，Ｃｌの長さで、ｎ1,kは点１から格子点ｋへ向かう単位ベクトルである。
【００６１】
たとえば、ｄM6,C2は点Ｍ６から点Ｃ２へ向かう辺Ｍ６，Ｃ２の長さ、ｎ1,6は格子点１から格子点６へ向かう方向の単位ベクトルである。(36)式を代入すれば、(33)式は格子点におけるφp，φnの代数方程式で記述できる。
【００６２】
次に(34)式の積分方程式を代数方程式に変換する方法を説明する。図６に示したコントロールボリュームに対して(34)式の左辺第２項は次式で与えられる。
【００６３】
【数１１】
0003660137.tif 000003
ここで、ｖdem（ｍ＝２〜６）は格子ｍにおける電子のドリフト速度ベクトルであり、φn,mを電子の擬フェルミポテンシャルとして、次式で与えられる。
【００６４】
【数１２】
ｖdem＝（−ｋＴ／ｑτpe）Λem（UφＥm／｜φＥm｜）＋ξem・ｖth  …(49)
Λem＝１−｜Ｅm／２｜sinｈ（ζm,φ）／（ηφ・sinｈ（ζm））   …(50)
ξem＝（coth（ζm）−１／ζm）（UＥm−（UＥm・UφＥm）UφＥm） …(51)
ここで、τpeは電子の運動量緩和時間である。また、ζm及びζm,pは、それぞれ熱電位（ｋＴ／ｑ）で規格化した格子ｍにおけるポテンシャル、及び電子の擬フェルミポテンシャルφn,mを考慮した、熱電位で規格化した格子ｍにおけるポテンシャルであり、次式で与えられる。
【００６５】
ζm＝−ｑ｜Ｅm｜ ｜ｌ｜／２ｋＴ                       …(52)
ηφ＝｜Ｅm／２＋▽φn,m｜                            …(53)
ζm,φ＝−ｑηφ｜ｌ｜／ｋＴ                           …(54)
さらに、φＥmは擬フェルミポテンシャルを考慮した格子ｍにおける電界ベクトルで、
φＥm＝Ｅm＋▽φn,m                                    …(55)
で与えられる。また、UφＥm及びUＥmは、それぞれ擬フェルミポテンシャルを考慮した格子ｍにおける単位電界ベクトル、及び格子ｍにおける単位電界ベクトルで、
【数１３】
UφＥm＝（Ｅm＋▽φn,m）／｜Ｅm＋▽φn,m｜              …(56)
UＥm＝Ｅm／｜Ｅm｜                                      …(57)
で与えられる。
【００６６】
先に述べたように、Ｅmは格子ｍ（ｍ＝２〜６）内で一定である。また、電子の擬フェルミポテンシャルφn,mは各格子内で座標の線形関数であるから▽φn,mも格子ｍ内で一定のベクトルである。即ち、ｖdemは定数ベクトルである。
【００６７】
一方(48)式におけるＮ6,2,2、Ｎ2,3,3、Ｎ3,4,4、Ｎ4,5,5、Ｎ5,6,6は以下の(58)式〜(62)式で与えられる。
【００６８】
【数１４】
Ｎ6,2,2＝ＮM6,C2ｎ1,6＋ＮC2,M2ｎ1,2                   …(58)
Ｎ2,3,3＝ＮM2,C3ｎ1,2＋ＮC3,M3ｎ1,3                   …(59)
Ｎ3,4,4＝ＮM3,C4ｎ1,3＋ＮC4,M4ｎ1,4                   …(60)
Ｎ4,5,5＝ＮM4,C5ｎ1,5＋ＮC5,M5ｎ1,5                   …(61)
Ｎ5,6,6＝ＮM5,C6ｎ1,5＋ＮC6,M6ｎ1,6                   …(62)
ここで、ＮMi,Cjは電子濃度ｎを辺Ｍｉ，Ｃｊ（ｉ＝２から６，ｊ＝２〜６）で線積分した関数で、次の(63)式で定義される。
【００６９】
ＮMi,Cj＝MiINTGCjｎ（ｘ）ｄｘ                         …(63)
辺Ｍｉ，Ｃｊ上でｎは指数関数的に変化するので(63)式は解析的に積分できる。
【００７０】
したがって、(48)式は格子点における電位Ψ及び電子の擬フェルミポテンシャルφnの代数方程式で記述できる。従って、(34)式も格子点における電位Ψと正孔の擬フェルミポテンシャルφp及び電子の擬フェルミポテンシャルφnの代数方程式で記述できることが説明できた。
【００７１】
(35)式も(34)式と同様に格子点における電位Ψと正孔の擬フェルミポテンシャルφp及び電子の擬フェルミポテンシャルφnの代数方程式で記述できる。
【００７２】
本発明の第１の実施の形態において解くべきデバイスシミュレーションの基本方程式は、全てのコントロールボリュームに対する(33)式〜(35)式の連立代数方程式である。この(33)式〜(35)式の連立代数方程式を解くに当たって、所定のバイアスにおける解の近似値を何らかの方法で推定し、解の初期値として設定する必要がある。本発明の第１の実施の形態では次の３つの方法を用いた。
【００７３】
（ａ）１つは所定のバイアスの直前のバイアスにおける解と、更にもう１回前のバイアスの解から線形外挿によって求める方法である。この方法は３回目以降のバイアスの初期値を推定する際に用いた。
【００７４】
（ｂ）もう１つの方法は直前のバイアスの解をそのまま用いる方法である。この方法は２回目のバイアスの初期値を求める際に用いた。
【００７５】
（ｃ）最後の方法は１回目の解、即ち熱平衡状態の解を求める際に用いる方法である。
【００７６】
図３は本発明の第１の実施の形態に係るデバイスシミュレーション方法のフローチャートである。手順を順に説明する。
【００７７】
（イ）「開始」の手順ではデバイスシミュレーションを開始する。本発明の第１の実施の形態では、ハードディスク等のプログラム記憶部６３に格納されたソフトウェアプログラムの起動を行う。
【００７８】
（ロ）ステップＳ３１においては、シミュレーション対象とする半導体装置の形状や不純物分布、各電極領域の名称、比誘電率、仕事関数などをデバイスシミュレータに入力する。本発明の第１の実施の形態では、周知のプロセスシミュレーターにより得られたる半導体装置の形状や不純物分布等の値を入力する。図４は本発明の第１の実施の形態において入力したデバイス構造を示す図である。この半導体装置はＭＯＳＦＥＴである。基板領域１４の材質はシリコンである。基板領域１４の形状は、基板領域１４の境界上の格子点の座標を反時計周りに順番に指定してデバイスシミュレータに入力した。なお、以下に示す他の領域についても、それぞれの領域の形状の入力方法は基板領域１４と同じである。シリコンの比誘電率は１１．７であるので、この値も入力した。不純物としては燐（Ｐ）、砒素（Ａｓ）、ホウ素（Ｂ）がシリコン中に分布している。シリコン領域の内部及び境界上の格子点において、燐、砒素、ホウ素のそれぞれの濃度を入力することにより、不純物分布をデバイスシミュレータに入力した。第１の絶縁体領域１６，第２の絶縁体領域１７及び第３の絶縁体領域１８の材質はシリコン酸化物（ＳｉＯ2）である。第１の絶縁体領域１６，第２の絶縁体領域１７及び第３の絶縁体領域１８の形状も基板領域１４と同じ方法で入力した。また、シリコン酸化物の比誘電率は３．９であるので、この値も入力した。本発明の第１の実施の形態における半導体装置には、図４示すようにソース電極領域１１，ドレイン電極領域１２、ゲート電極領域１３、及び基板電極領域１５の４つの電極領域がある。それらの電極領域１１，１２，１３，１５の形状も基板領域１４と同じ方法でそれぞれ入力した。このうち、ソース電極領域１１，ドレイン電極領域１２、基板電極領域１５の材質はアルミニウムであるので、仕事関数の値として４．２ｅＶを入力した。ゲート電極領域１３の材質はポリシリコンである。そこで、比誘電率の値として１１．７を入力した。また、ゲート電極領域１３における不純物分布を、基板領域１４と同じ方法で入力した。
【００７９】
（ハ）ステップＳ３２においては、周知の格子生成技術を用いて三角形格子を生成する。図５は生成した格子を示す図である。基板領域１４の内部にある三角形格子は基板領域１４に属する。その他の領域についても同様である。また、基板領域１４の内部及び境界上にある格子点は基板領域１４に属する。基板領域１４、第１の絶縁体領域１６，第２の絶縁体領域１７、第３の絶縁体領域１８、ゲート電極領域１３に属する三角形格子には、それぞれの領域の比誘電率を設定する。ソース電極領域１１、ドレイン電極領域１２、基板電極領域１５に属する三角形格子には、それぞれの領域１１，１２，１５の仕事関数を設定する。また、基板領域１４、ゲート電極領域１３に属する格子には不純物濃度を設定する。
【００８０】
（ニ）ステップＳ３３においては、所定の「物理モデル」を設定する。本発明の第１の実施の形態では、例えば、基板領域１４に属する格子点で移動度を算出するための移動モデルとそのパラメーターを設定する。例えば、周知の移動度モデルとそのパラメーターを設定することが出来る。更に、基板領域１４に属する格子点でキャリア生成再結合率Ｕを算出するためのキャリア生成再結合率モデルとそのパラメーターを設定する。本発明の第１の実施の形態では、生成された格子の内部のそれぞれにおいて、線形若しくは指数関数的に変化するキャリア生成再結合率を設定する。
【００８１】
（ホ）ステップＳ３４においては、ステップＳ３２において生成した三角形格子の内部に、線形に変化するキャリアの擬フェルミポテンシャルを設定する。即ち、それぞれの格子の内部に、電子の擬フェルミポテンシャルφe及び正孔の擬フェルミポテンシャルφpを設定する。さらに、ドナー濃度ＮD、アクセプター濃度ＮAについても、生成された格子の内部のそれぞれにおいて、線形若しくは指数関数的に変化する値を設定する。
【００８２】
（へ）ステップＳ３５においては、ソース電極領域１１，ドレイン電極領域１２、ゲート電極領域１３、及び基板電極領域１５の４つの電極領域に印加するバイアスを、シミュレーションする電気特性に応じて入力する。例えば、ＭＯＳＦＥＴの三極管特性をシミュレーションする場合には、ドレイン電極領域１２に例えば０．１Ｖを印加した後、ゲート電極領域１３を０Ｖから例えば３Ｖまで０．１Ｖずつ変化させるコマンドをデバイスシミュレータに入力する。すると、デバイスシミュレータはまず、全ての電極領域１１，１２，１３，１５のバイアスが０Ｖである状態の解、即ち熱平衡状態の解を以下で述べる方法により算出し、次にドレイン電極領域１２のバイアス即ちドレインバイアスだけを０．１Ｖに変化させた状態の解を算出する。次に、ゲートバイアスを０．１Ｖに変化させた解、０．２Ｖに変化させた解を順に算出し、ゲートバイアスが３Ｖの解まで算出する。
【００８３】
（ト）あるバイアスにおける解はステップＳ３6の「初期値の推定」からステップＳ４1の「解の保存」までの手順で算出する。これらの手順で行うことは、擬フェルミポテンシャルから、それぞれの格子内の任意の点におけるキャリア濃度を求め、これを用いて、(１)式乃至（３）式で示されたポアッソンの方程式及びキャリアの連続の方程式を、所定の境界条件の下で解くことである。この際、(１)式乃至（３）式で示された基本方程式を,(33)式〜(35)式で示された連立一次方程式に変換し、全てのコントロールボリュームに対してこの(33)式〜(35)式の連立一次方程式を解く。具体的には、周知のニュートン法によって、ステップ３６からステップＳ４１までの手順によって解く。ステップＳ３６においては、所定のバイアスにおける解の近似値を何らかの方法で推定し、解の初期値として設定する必要がある。本発明の第１の実施の形態では前述した（ａ）乃至（ｃ）の３つの方法を用いた。まず、絶縁体領域に属する格子点の電位を０Ｖとする。次に、基板領域及びゲート電極領域に属する格子点では、擬フェルミポテンシャルを０Ｖとし、ｐ型半導体領域であれば正孔濃度ｐが、ｎ型半導体領域であれば電子濃度ｎが不純物濃度に一致するように、電位Ψを決めて初期値とする。ステップＳ３７では、ステップＳ３４で設定した擬フェルミポテンシャルから、それぞれの格子内の任意の点におけるキャリア濃度を求め、基本方程式に対する係数行列と残差ベクトルを設定する。この係数行列と残差ベクトルを設定に際しては、ステップＳ３３，３４で、線形若しくは指数関数的に変化すると仮定して設定されたドナー濃度ＮD、アクセプター濃度ＮA、キャリア生成再結合率Ｕ等も用いられる。そして、ステップＳ３８では、ステップＳ３７で設定した係数行列と残差ベクトルを基に、周知の技術により解の修正値を算出する。さらに、ステップＳ３９では、周知の技術により、基本方程式の解にステップＳ３８で算出した解の修正値を加えて、解を修正する。その後、ステップＳ４０では、周知の技術により、基本方程式の解が収束したか否かを判定する。収束していなければステップＳ３７に戻り「係数行列と残差ベクトルの設定」の手順を次に行う。収束していればステップＳ４１の「解の保存」の手順を次に行う。ステップＳ４１では、所定のバイアスの解としてその解を保存する。
【００８４】
（チ）ステップＳ４２では、全バイアスの計算が終了したか否かを判定する。終了していなければステップＳ３５に戻り、「バイアスの設定」の手順を行い、ステップＳ３６以降の手順を繰り返す。ステップＳ４２で、全バイアスの計算が終了していれば、ステップＳ４３の「結果の出力」の手順を次に行う。ステップＳ４３では、シミュレーション結果の出力を行う。例えば、各バイアスにおける各電極領域の端子電流や、各格子点における電位やキャリア濃度などを出力する。ステップＳ４３の「結果の出力」の手順が終われば、求める半導体装置の電気的特性が得られ、本発明の第１の実施の形態に係るシミュレーションプログラムを終了する。
【００８５】
（第２の実施の形態）
本発明の第２の実施の形態は先に示した第１の実施の形態と共通する内容を含むので、異なる点のみ説明する。
【００８６】
第２の実施の形態においては、前述した(34)式の左辺第２項を次式によって算出する。
【００８７】
【数１５】
0003660137.tif 000004
ここで、(50）式のΛemを用いて、
Ｆe,m＝（−ｑ／τpe）Λem                             …(65)
である。またＮMj,Ck,Mk1は頂点がＭｊ，Ｃｋ，Ｍｋ１である四角形領域での電子濃度ｎの積分であり、次式で定義される。
【００８８】
ＮMj,Ck,Mk1＝INTGMj,Ck,Mk1ｎｄν                      …(66)
電位Ψと電子の擬フェルミポテンシャルφnは各格子において位置の線形関数であり、しかも、積分領域が多角形なので、(66)式は解析的に積分できる。
【００８９】
したがって(48)式の右辺は格子点における電位Ψと正孔の擬フェルミポテンシャルφp及び電子の擬フェルミポテンシャルφnの代数方程式で記述できる。即ち、第１の実施形態と同じ手順でデバイスシミュレーションを実行できる。
【００９０】
（半導体装置の製造方法）
本発明の半導体装置の製造方法の各工程の流れは以下のようになる。まず、デバイスシミュレーションに先立ち、プロセスシミュレーションが実行される。プロセスシミュレーションによって得られた半導体装置中の不純物や欠陥の分布等のデータは、デバイスシミュレータに入力される。このデバイスシミュレーションを行う際には、プロセスシミュレーションで得られた素子構造と不純物分布等の結果と同時に、印加電圧、電流などの電気的な境界条件を与えるための入力データが加えられる。
【００９１】
そして、デバイスシミュレータにおいては、第１及び第２の実施の形態で説明したように、製造しようとする半導体装置の内部及び周囲に有限個の格子点を設置し、複数の格子を生成し、格子の内部のそれぞれに、線形に変化するキャリアの擬フェルミポテンシャルを設定し、擬フェルミポテンシャルから、それぞれの格子内の任意の点におけるキャリア濃度を求め、キャリア濃度を用いて、ポアッソンの方程式及びキャリアの連続の方程式を、所定の境界条件の下で解き、半導体装置の電気的特性を得る。「半導体装置の電気的特性」とは、例えば、電流−電圧特性、インピーダンス特性、Ｓパラメータや高周波特性等のデバイス特性である。
【００９２】
さらに、必要に応じてデバイスシミュレーションの結果としてのデバイス特性は、回路シミュレーションの入力データとなり回路特性を得るようにしても良い。
【００９３】
デバイスシミュレーションまたは回路シミュレーションにより、得られた特性が作ろうとしている所望の特性になるかどうかを調べ、所望の特性であれば実際の半導体装置の製造工程に取りかかる。もし、所望の特性にならないときには、考えた製造工程では作りたい素子は作れないので、製造工程の条件を変更したり、工程の順番など手順そのものを変更し再度プロセスシミュレーションを行い、このプロセスシミュレーションの結果を入力データとしてデバイスシミュレーションを行う。
【００９４】
さらに、実際の半導体装置の製造工程の結果得られた現実の半導体装置の特性を測定し、当初の要求仕様を満足するか否か評価する。この評価により実際に製造された半導体装置の特性が要求仕様を満足しなければ、設計変更がなされ、再度プロセスシミュレーションを行う。そして、このプロセスシミュレーションの結果を入力データとしてデバイスシミュレーションを行うという一連の手順からなるループが繰り返される。
【００９５】
ＬＳＩ等の半導体装置の分野では、研究（設計）から開発までの期間の短さを競っている。このような半導体産業における競争の現実を考慮すれば、シミュレーション期間はなるべく短期、且つ正確に行わなければならない。本発明によれば高精度なデバイスシミュレーションが必要な複雑な半導体装置の研究（設計）から開発までのループの周期が飛躍的に短縮される。
【００９６】
【発明の効果】
本発明によれば複雑且つ微細な半導体装置の電気特性を正確にシミュレーションできるシミュレーション方法を提供できる。
【００９７】
本発明によれば複雑且つ微細な半導体装置の電気特性を正確にシミュレーションできるシミュレータを提供できる。
【００９８】
本発明によれば、複雑且つ微細な半導体装置の電気特性を正確にシミュレーションできるプログラムを記録した記録媒体を提供できる。
【００９９】
本発明によれば所望の特性の半導体装置を従来よりも短期間で開発できる。なぜなら、試作に先立つ半導体装置の特性予測を従来よりも正確に行えるので、より少ない試作回数で複雑且つ微細な半導体装置を開発できるからである。
【０１００】
従って、本発明の工業的利益およびその需要性は極めて高い。
【図面の簡単な説明】
【図１】本発明の第１の実施の形態に係るシミュレータの模式的なブロック図である。
【図２】本発明の第１の実施の形態に係るシミュレータの鳥瞰図である。
【図３】本発明の第１の実施の形態に係るシミュレーション方法を示すフローチャートである。
【図４】本発明の第１の実施の形態に係るシミュレーション方法において用いたデバイス構造を示す図である。
【図５】本発明の第１の実施の形態に係るシミュレーション方法において生成した三角形格子を示す図である。
【図６】本発明の第１の実施の形態に係るシミュレーション方法において用いた格子と格子点とコントロールボリュームの関係を示す図である。
【図７】格子点と格子とコントロールボリュームの関係を示す図である。
【図８】第１の従来技術を説明するための図で、格子点の配置を示す図である。
【図９】第１の従来技術を説明するための図で、記号の定義を示す図である。
【図１０】第２の従来技術を示す図である。斜めの直線は移流電流ベクトルの方向を示す。
【図１１】第１の従来技術と第２の従来技術の計算誤差を示す図である。横軸は図１１におけるθである。
【図１２】本発明における数学的表現の定義を示す図である。
【符号の説明】
１１  ソース電極領域
１２  ドレイン電極領域
１３  ゲート電極領域
１４  基板領域
１５  基板電極領域
１６  第１の絶縁体領域
１７  第２の絶縁体領域
１８  第３の絶縁体領域
【図１】
0003660137.tif 000006
【図２】
0003660137.tif 000007
【図３】
0003660137.tif 000008
【図４】
0003660137.tif 000009
【図５】
0003660137.tif 000010
【図６】
0003660137.tif 000011
【図７】
0003660137.tif 000012
【図８】
0003660137.tif 000013
【図９】
0003660137.tif 000014
【図１０】
0003660137.tif 000015
【図１１】
0003660137.tif 000016
【図１２】
0003660137.tif 000017