JP4106880B2 - 単結晶内欠陥の密度分布及びサイズ分布のシミュレーション方法 - Google Patents
単結晶内欠陥の密度分布及びサイズ分布のシミュレーション方法 Download PDFInfo
- Publication number
- JP4106880B2 JP4106880B2 JP2001227192A JP2001227192A JP4106880B2 JP 4106880 B2 JP4106880 B2 JP 4106880B2 JP 2001227192 A JP2001227192 A JP 2001227192A JP 2001227192 A JP2001227192 A JP 2001227192A JP 4106880 B2 JP4106880 B2 JP 4106880B2
- Authority
- JP
- Japan
- Prior art keywords
- single crystal
- temperature
- melt
- distribution
- mesh
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Expired - Lifetime
Links
Images
Landscapes
- Crystals, And After-Treatments Of Crystals (AREA)
Description
【０００１】
【発明の属する技術分野】
本発明は、チョクラルスキー（以下、ＣＺという。）法にて引上げられるシリコン等の単結晶内欠陥の密度分布及びサイズ分布をコンピュータシミュレーションする方法に関するものである。
【０００２】
【従来の技術】
従来、この種のシミュレーション方法として、図８に示すように、総合伝熱シミュレータを用いてＣＺ法によるシリコン単結晶４引上げ時の引上げ機１内のホットゾーン構造及びそのシリコン単結晶４の引上げ速度に基づいて、シリコン融液２の熱伝導率を操作することによりシリコン融液２の内部温度分布を予測し、この内部温度分布からシリコン単結晶４のメッシュの座標及び温度をそれぞれ求め、更にシリコン単結晶４内の格子間シリコン及び空孔の拡散係数及び境界条件に基づいて拡散方程式を解くことにより、上記格子間シリコン及び空孔の密度分布をコンピュータを用いて求める方法が知られている。このシミュレーション方法では、ホットゾーンの各部材がメッシュ分割されてモデル化される。特にシリコン融液２のメッシュは計算時間を短くするために１０ｍｍ程度と比較的粗く設定される。
【０００３】
【発明が解決しようとする課題】
しかし、上記従来の格子間シリコン及び空孔の密度分布のシミュレーション方法では、実際の引上げ機においては発生するシリコン融液の対流を考慮しておらず、またシリコン融液のメッシュが比較的粗いため、固液界面形状の再現性が悪く、精度の良い結晶内温度分布を提供できない。そのため、シリコン単結晶内の格子間シリコン及び空孔の密度分布が実測値と大幅に相違する問題点があった。また上記従来のシミュレーション方法では、結晶冷却過程に比較的高温で生成されるシリコン単結晶内の欠陥、例えばボイドや酸素析出物の種類や、密度分布及びサイズ分布は判らなかった。
本発明の目的は、融液の対流を考慮して成長中の単結晶内の温度分布を解析した後に、融液から切離された単結晶の冷却過程を考慮して解析することにより、単結晶内のボイド及び高温酸素析出物の密度分布及びサイズ分布を正確に予測できる、単結晶内欠陥の密度分布及びサイズ分布のシミュレーション方法を提供することにある。
【０００４】
【課題を解決するための手段】
請求項１に係る発明は、図１〜図３及び図７に示すように、引上げ機１１による単結晶１４の融液１２からの引上げ開始時から単結晶１４の冷却完了時までの引上げ機１１のホットゾーンをメッシュ構造でモデル化する第１ステップと、ホットゾーンの各部材毎にメッシュをまとめかつこのまとめられたメッシュに対する各部材の物性値とともに単結晶１４の引上げ長及びこの引上げ長に対応する単結晶１４の引上げ速度をそれぞれコンピュータに入力する第２ステップと、各部材の表面温度分布をヒータの発熱量及び各部材の輻射率に基づいて求める第３ステップと、各部材の表面温度分布及び熱伝導率に基づいて熱伝導方程式を解くことにより各部材の内部温度分布を求めた後に融液１２が乱流であると仮定して得られた乱流モデル式及びナビエ・ストークスの方程式を連結して解くことにより対流を考慮した融液１２の内部温度分布を更に求める第４ステップと、単結晶１４及び融液１２の固液界面形状を単結晶の三重点Ｓを含む等温線に合せて求める第５ステップと、第３ステップから第５ステップを三重点Ｓが単結晶１４の融点になるまで繰返し引上げ機１１内の温度分布を計算して単結晶１４のメッシュの座標及び温度を求めこれらのデータをそれぞれコンピュータに入力する第６ステップと、単結晶１４の引上げ長及び引上げ高さを段階的に変えて第１ステップから第６ステップまでを繰返し引上げ機１１内の温度分布を計算して単結晶１４のメッシュの座標及び温度を求めこれらのデータをそれぞれコンピュータに入力する第７ステップと、単結晶１４の融液１２からの引上げ開始時から単結晶１４の冷却完了時までの時間を所定の間隔毎に区切りこの区切られた時間間隔毎に第７ステップで求めた単結晶１４のメッシュの座標及び温度のデータから単結晶１４の引上げ長及び引上げ高さと単結晶１４内の温度分布とを求める第８ステップと、単結晶１４内の空孔及び格子間原子の拡散係数及び境界条件に基づいて拡散方程式を解くことにより所定の時間間隔の経過した後の空孔及び格子間原子の密度分布を求める第９ステップと、空孔の密度分布に基づいてボイドの形成開始温度を求める第１０ステップと、単結晶１４内のそれぞれのメッシュの格子点における温度が次第に低下してボイドの形成開始温度になったときのボイドの密度を求める第１１ステップと、単結晶１４内のそれぞれのメッシュの格子点における温度がボイドの形成開始温度より低いときのボイドの半径を求める第１２ステップと、第８ステップから第１２ステップを単結晶１４の冷却が完了するまで繰返す第１３ステップとを含むコンピュータを用いて単結晶内欠陥の密度分布及びサイズ分布のシミュレーションを行う方法であって、乱流モデル式が次の式（２）で表されるｋｌ−モデル式であり、このモデル式の乱流パラメータＣとして０．４〜０．６の範囲内の任意の値が用いられたことを特徴とする、コンピュータを用いて単結晶内欠陥の密度分布及びサイズ分布のシミュレーションを行う方法である。
【数２】
ここで、κ t は融液の乱流熱伝導率であり、ｃは融液の比熱であり、Ｐｒ t はプラントル数であり、ρは融液の密度であり、ｄは融液を貯留するるつぼ壁からの距離であり、ｋは融液の平均流速に対する変動成分の二乗和である。
【０００５】
この請求項１に記載されたシミュレーション方法では、融液１２の対流を考慮して融液１２から成長する単結晶１４内の温度分布を求めるだけでなく、更に冷却過程における単結晶１４内の温度分布までも求めることによって、即ち融液１２から切離された単結晶１４の冷却過程における単結晶１４の徐冷及び急冷の効果を考慮して解析することによって、単結晶１４内のボイドの密度分布及びサイズ分布を正確に予測できる。
【０００６】
請求項２に係る発明は、図４〜図７に示すように、請求項１記載の第１〜第９ステップと、空孔の密度分布に基づいて高温酸素析出物の形成開始温度を求める第１０’ステップと、単結晶１４内のそれぞれのメッシュの格子点における温度が次第に低下して高温酸素析出物の形成開始温度になったときの高温酸素析出物の密度を求める第１１’ステップと、単結晶１４内のそれぞれのメッシュの格子点における温度が高温酸素析出物の形成開始温度より低いときの高温酸素析出物の半径を求める第１２’ステップと、第８、第９、第１０’、第１１’及び第１２’ステップを単結晶１４の冷却が完了するまで繰返す第１３’ステップとを含むコンピュータを用いて単結晶内欠陥の密度分布及びサイズ分布のシミュレーションを行う方法であって、乱流モデル式が次の式（２）で表されるｋｌ−モデル式であり、このモデル式の乱流パラメータＣとして０．４〜０．６の範囲内の任意の値が用いられたことを特徴とする、コンピュータを用いて単結晶内欠陥の密度分布及びサイズ分布のシミュレーションを行う方法である。
【数２】
ここで、κ t は融液の乱流熱伝導率であり、ｃは融液の比熱であり、Ｐｒ t はプラントル数であり、ρは融液の密度であり、ｄは融液を貯留するるつぼ壁からの距離であり、ｋは融液の平均流速に対する変動成分の二乗和である。
【０００７】
この請求項２に記載されたシミュレーション方法では、融液１２の対流を考慮して融液１２から成長する単結晶１４内の温度分布を求めるだけでなく、更に冷却過程における単結晶１４内の温度分布までも求めることによって、即ち融液１２から切離された単結晶１４の冷却過程における単結晶１４の徐冷及び急冷の効果を考慮して解析することによって、単結晶１４内の高温酸素析出物の密度分布及びサイズ分布を正確に予測できる。
なお、本明細書において「高温酸素析出物」とは、ボイド欠陥の内側に析出する内壁酸化膜を球形に近似して扱ったものをいう。また酸素析出物を高温と限定したのは、ボイド欠陥より低温で析出するＢＭＤ（Bulk Micro Defect）等の析出核となるものと区別するためである。
【０００８】
【発明の実施の形態】
次に本発明の実施の形態を図面に基づいて説明する。
図７に示すように、シリコン単結晶引上げ機１１のチャンバ内には、シリコン融液１２を貯留する石英るつぼ１３が設けられる。この石英るつぼ１３は図示しないが黒鉛サセプタ及び支軸を介してるつぼ駆動手段に接続され、るつぼ駆動手段は石英るつぼ１３を回転させるとともに昇降させるように構成される。また石英るつぼ１３の外周面は石英るつぼ１３から所定の間隔をあけてヒータ（図示せず）により包囲され、このヒータは保温筒（図示せず）により包囲される。ヒータは石英るつぼ１３に投入された高純度のシリコン多結晶体を加熱・溶融してシリコン融液１２にする。またチャンバの上端には図示しないが円筒状のケーシングが接続され、このケーシングには引上げ手段が設けられる。引上げ手段はシリコン単結晶１４を回転させながら引上げるように構成される。
【０００９】
このように構成されたシリコン単結晶引上げ機１１におけるシリコン単結晶１４内のボイドの密度分布及びサイズ分布のシミュレーション方法を図１〜図３及び図７に基づいて説明する。
先ず第１ステップとして、シリコン単結晶１４を所定長さＬ1（例えば１００ｍｍ）まで引上げた状態におけるシリコン単結晶引上げ機１１のホットゾーンの各部材、即ちチャンバ，石英るつぼ１３，シリコン融液１２，シリコン単結晶１４，黒鉛サセプタ，保温筒等をメッシュ分割してモデル化する。具体的には上記ホットゾーンの各部材のメッシュ点の座標データをコンピュータに入力する。このときシリコン融液１２のメッシュのうちシリコン単結晶１４の径方向のメッシュであってかつシリコン融液１２のシリコン単結晶１４直下の一部又は全部のメッシュ（以下、径方向メッシュという。）を０．０１〜５．００ｍｍ、好ましくは０．２５〜１．００ｍｍに設定する。またシリコン融液１２のメッシュのうちシリコン単結晶１４の長手方向のメッシュであってかつシリコン融液１２の一部又は全部のメッシュ（以下、長手方向メッシュという。）を０．０１〜５．００ｍｍ、好ましくは０．１〜０．５ｍｍに設定する。
【００１０】
径方向メッシュを０．０１〜５．００ｍｍの範囲に限定したのは、０．０１ｍｍ未満では計算時間が極めて長くなり、５．００ｍｍを越えると計算が不安定になり、繰返し計算を行っても固液界面形状が一定に定まらなくなるからである。また長手方向メッシュを０．０１〜５．００ｍｍの範囲に限定したのは、０．０１ｍｍ未満では計算時間が極めて長くなり、５．００ｍｍを越えると固液界面形状の計算値が実測値と一致しなくなるからである。なお、径方向メッシュの一部を０．０１〜５．００の範囲に限定する場合には、シリコン単結晶１４直下のシリコン融液１２のうちシリコン単結晶１４外周縁近傍のシリコン融液１２を上記範囲に限定することが好ましく、長手方向メッシュの一部を０．０１〜５．００の範囲に限定する場合には、シリコン融液１２の液面近傍及び底近傍を上記範囲に限定することが好ましい。
【００１１】
第２ステップとして上記ホットゾーンの各部材毎にメッシュをまとめ、かつこのまとめられたメッシュに対して各部材の物性値をそれぞれコンピュータに入力する。例えば、チャンバがステンレス鋼にて形成されていれば、そのステンレス鋼の熱伝導率，輻射率，粘性率，体積膨張係数，密度及び比熱がコンピュータに入力される。またシリコン単結晶１４の引上げ長及びこの引上げ長に対応するシリコン単結晶１４の引上げ速度と、後述する乱流モデル式（１）の乱流パラメータＣとをコンピュータに入力する。
【００１２】
第３ステップとして、ホットゾーンの各部材の表面温度分布をヒータの発熱量及び各部材の輻射率に基づいてコンピュータを用いて求める。即ち、ヒータの発熱量を任意に設定してコンピュータに入力するとともに、各部材の輻射率から各部材の表面温度分布をコンピュータを用いて求める。次に第４ステップとしてホットゾーンの各部材の表面温度分布及び熱伝導率に基づいて熱伝導方程式（１）をコンピュータを用いて解くことにより各部材の内部温度分布を求める。ここでは、記述を簡単にするためｘｙｚ直交座標系を用いたが、実際の計算では円筒座標系を用いる。
【００１３】
【数１】
ここで、ρは各部材の密度であり、ｃは各部材の比熱であり、Ｔは各部材の各メッシュ点での絶対温度であり、ｔは時間であり、λx，λy及びλzは各部材の熱伝導率のｘ，ｙ及びｚ方向成分であり、ｑはヒータの発熱量である。
一方、シリコン融液１２に関しては、上記熱伝導方程式（１）でシリコン融液１２の内部温度分布を求めた後に、このシリコン融液１２の内部温度分布に基づき、シリコン融液１２が乱流であると仮定して得られた乱流モデル式（２）及びナビエ・ストークスの方程式（３）〜（５）を連結して、シリコン融液１２の内部流速分布をコンピュータを用いて求める。
【００１４】
【数２】
ここで、κtはシリコン融液１２の乱流熱伝導率であり、ｃはシリコン融液１２の比熱であり、Ｐｒtはプラントル数であり、ρはシリコン融液１２の密度であり、Ｃは乱流パラメータであり、ｄはシリコン融液１２を貯留する石英るつぼ１３壁からの距離であり、ｋはシリコン融液１２の平均流速に対する変動成分の二乗和である。
【００１５】
【数３】
ここで、ｕ，ｖ及びｗはシリコン融液１２の各メッシュ点での流速のｘ，ｙ及びｚ方向成分であり、νlはシリコン融液１２の分子動粘性係数（物性値）であり、νtはシリコン融液１２の乱流の効果による動粘性係数であり、Ｆx，Ｆy及びＦzはシリコン融液１２に作用する体積力のｘ，ｙ及びｚ方向成分である。
【００１６】
上記乱流モデル式（２）はｋｌ（ケイエル）−モデル式と呼ばれ、このモデル式の乱流パラメータＣは０．４〜０．６の範囲内の任意の値が用いられることが好ましい。乱流パラメータＣを０．４〜０．６の範囲に限定したのは、０．４未満又は０．６を越えると計算により求めた界面形状が実測値と一致しないという不具合があるからである。また上記ナビエ・ストークスの方程式（３）〜（５）はシリコン融液１２が非圧縮性であって粘度が一定である流体としたときの運動方程式である。
上記求められたシリコン融液１２の内部流速分布に基づいて熱エネルギ方程式（６）を解くことにより、シリコン融液１２の対流を考慮したシリコン融液１２の内部温度分布をコンピュータを用いて更に求める。
【００１７】
【数４】
ここで、ｕ，ｖ及びｗはシリコン融液１２の各メッシュ点での流速のｘ，ｙ及びｚ方向成分であり、Ｔはシリコン融液１２の各メッシュ点での絶対温度であり、ρはシリコン融液１２の密度であり、ｃはシリコン融液１２の比熱であり、κlは分子熱伝導率（物性値）であり、κtは式（１）を用いて計算される乱流熱伝導率である。
【００１８】
次いで第５ステップとして、シリコン単結晶１４及びシリコン融液１２の固液界面形状を図２の点Ｓで示すシリコンの三重点Ｓ（固体と液体と気体の三重点(tri-junction)）を含む等温線に合せてコンピュータを用いて求める。第６ステップとして、コンピュータに入力するヒータの発熱量を変更し（次第に増大し）、上記第３ステップから第５ステップを三重点がシリコン単結晶１４の融点になるまで繰返した後に、引上げ機１１内の温度分布を計算してシリコン単結晶のメッシュの座標及び温度を求め、これらのデータをコンピュータに記憶させる。
【００１９】
次に第７ステップとして、シリコン単結晶１４の引上げ長Ｌ1にδ（例えば５０ｍｍ）だけ加えて上記第１ステップから第６ステップまでを繰返した後に、引上げ機１１内の温度分布を計算してシリコン単結晶１４のメッシュの座標及び温度を求め、これらのデータをコンピュータに記憶させる。この第７ステップは、シリコン単結晶１４の引上げ長Ｌ1が長さＬ2（Ｌ2はシリコン融液１２から切離されたときのシリコン単結晶１４の長さ（成長完了時の結晶長）である。）に達してシリコン単結晶１４がシリコン融液１２から切離された後、更にシリコン単結晶１４が引上げられてその高さＨ1（Ｈ1はシリコン単結晶１４の直胴開始部からシリコン融液１２の液面までの距離である（図７）。）がＨ2（Ｈ2は冷却完了時のシリコン単結晶１４の直胴開始部からシリコン融液１２の液面までの距離である。）に達するまで、即ちシリコン単結晶１４の冷却が完了するまで行われる。なお、シリコン単結晶１４がシリコン融液１２から切離された後は、シリコン単結晶１４の引上げ高さＨ1にδ（例えば５０ｍｍ）だけ加え、上記と同様に上記第１ステップから第６ステップまでを繰返す。
【００２０】
シリコン単結晶１４の引上げ高さＨ1がＨ2に達すると、第８ステップに移行する。第８ステップでは、シリコン単結晶１４をシリコン融液１２から成長させて引上げ始めたときｔ0から、シリコン単結晶１４をシリコン融液１２から切離して更にシリコン単結晶１４を引上げ、その冷却が完了したときｔ1までの時間を、所定の間隔Δｔ秒（微小時間間隔）毎に区切る。このときシリコン単結晶１４内の格子間シリコン及び空孔の拡散係数及び境界条件のみならず、後述するボイド及び高温酸素析出物の密度分布及びサイズ分布を求めるための式に用いられる定数をそれぞれコンピュータに入力する。上記区切られた時間間隔Δｔ秒毎に、第７ステップで求めたシリコン単結晶１４のメッシュの座標及び温度のデータから、シリコン単結晶１４の引上げ長Ｌ1及び引上げ高さＨ1と、シリコン単結晶１４内の温度分布とを求める。
【００２１】
即ち、第１〜第７ステップでシリコン単結晶のメッシュの座標及び温度を引上げ長δ毎に求め、シリコン単結晶を例えば５０ｍｍ引上げるのに数十分要するため、この数十分間でのシリコン単結晶のメッシュの温度変化を時間の関数として微分することにより、時刻ｔ0からΔｔ秒後におけるシリコン単結晶１４の引上げ長Ｌ1及び引上げ高さＨ1とシリコン単結晶１４内の温度分布を算出する。次にシリコン単結晶１４内の空孔及び格子間シリコンの拡散係数及び境界条件に基づいて拡散方程式を解くことにより、Δｔ秒経過後の空孔及び格子間シリコンの密度分布を求める（第９ステップ）。
具体的には、空孔の密度Ｃvの計算式が次の式（７）で示され、格子間シリコンの密度Ｃiの計算式が次の式（８）で示される。式（７）及び式（８）において、密度Ｃv及び密度Ｃiの経時的進展を計算するために、空孔と格子間シリコンの熱平衡がシリコン単結晶の全表面で維持されると仮定する。
【００２２】
【数５】
ここで、Ｋ1及びＫ2は定数であり、Ｅi及びＥvはそれぞれ格子間シリコン及び空孔の形成エネルギーであり、Ｃv e及びＣi eの肩付き文字ｅは平衡量、ｋはボルツマン定数、Ｔは絶対温度を意味する。
上記平衡式は時間で微分され、空孔及び格子間シリコンに対してそれぞれ次の式（９）及び式（１０）になる。
【００２３】
【数６】
【００２４】
ここで、Θ(ｘ)はヘビサイド関数（Heaviside function）である。即ち、ｘ＜０のときΘ(ｘ)＝０であり、かつｘ＞０のときΘ(ｘ)＝１である。またＴnは高温酸素析出物の形成開始温度Ｔpとボイドの形成開始温度Ｔv0とを比較したときの高い方の温度である。更に式（９）及び（１０）のそれぞれ右側の第１項はフィックの拡散式であり、右側の第１項中のＤv及びＤiは、次の式（１１）及び（１２）に表される拡散係数である。
【００２５】
【数７】
【００２６】
ここで△Ｅv及び△Ｅiはそれぞれ空孔及び格子間シリコンの活性化エネルギーであり、ｄv及びｄiはそれぞれ定数である。また式（９）及び式（１０）のそれぞれ右側の第２項中の
【数８】
は熱拡散による空孔及び格子間シリコンの活性化エネルギーであり、ｋBはボルツマン定数である。式（９）及び式（１０）のそれぞれ右側の第３項のｋivは空孔及び格子間シリコンペアの再結合定数である。式（９）及び式（１０）のそれぞれ右側の第４項のＮv0はボイドの密度であり、ｒv0はボイドの半径であり、更に式（９）の右側の第５項のＮpは高温酸素析出物の密度であり、Ｒpは高温酸素析出物の半径であり、γはＳｉＯ2析出物が歪まずに析出するために必要な酸素１原子当りの空孔消費量である。
【００２７】
一方、上記式（９）が成立つのは、空孔が析出するための流速が十分大きく、ＳｉマトリックスとＳｉＯ2との単位質量当りの体積差を埋められる場合、即ちＤv（Ｃv−Ｃv e）≧γＤ0Ｃ0の場合である。上記以外の場合には、次式（１３）が成立つ。
【００２８】
【数９】
次に第１０ステップとして、上記拡散方程式を解くことにより求めた空孔の密度Ｃv分布に基づいて、ボイドの形成開始温度Ｔv0を次の式（１４）から求める。
【００２９】
【数１０】
ここで、Ｃvmはシリコン融液の融点Ｔmでの空孔平衡濃度であり、Ｅvは空孔形成エネルギーであり、Ｔmはシリコン単結晶１４の融点温度である。またσvはシリコン単結晶１４の結晶面（１１１）における界面エネルギーであり、ρはシリコン単結晶１４の密度であり、ｋBはボルツマン定数である。
【００３０】
第１１ステップとして、シリコン単結晶１４内のそれぞれのメッシュの格子点における温度が次第に低下してボイドの形成開始温度Ｔv0になったときに、次の近似式（１５）を用いてボイドの密度Ｎv0を求める。
【００３１】
【数１１】
更に第１２ステップとして、シリコン単結晶１４内のそれぞれのメッシュの格子点における温度がボイドの形成開始温度Ｔv0より低いときのボイドの半径ｒv0を、次の式（１６）から求める。
【００３２】
【数１２】
ここで、ｔ1はシリコン単結晶１４のメッシュの格子点における温度がボイドの形成開始温度Ｔv0まで低下したときの時刻であり、ｒcrはボイドの臨界径である。上記第８ステップから第１２ステップをシリコン単結晶１４の冷却が完了するまで繰返す（第１３ステップ）。なお、式（９）〜式（１６）は連成してコンピュータにより解く。
【００３３】
一方、上記シリコン単結晶引上げ機１１におけるシリコン単結晶１４内の高温酸素析出物の密度分布及びサイズ分布のシミュレーション方法を図４〜図７に基づいて説明する。
第１〜第９ステップまでは、上記ボイドの密度分布及びサイズ分布を求める場合と同一である。
次いで第１０’ステップとして、拡散方程式（９）及び（１０）を解くことにより求めた空孔の密度Ｃv分布に基づいて、高温酸素析出物の形成開始温度Ｔpを次の式（１７）から求める。
【００３４】
【数１３】
ここで、Ｃ0は酸素濃度であり、Ｃ0mはシリコン融液１２の融点Ｔmでの酸素平衡濃度であり、Ｅ0は酸素溶解エネルギーである。またＥvは空孔形成エネルギーであり、σpはシリコン単結晶１４内のＳｉとＳｉＯ2との界面エネルギーである。更にγはＳｉＯ2析出物が歪まずに析出するために必要な酸素１原子当りの空孔消費量であり、その値は０．６８である。
【００３５】
次に第１１’ステップとして、シリコン単結晶１４内のそれぞれのメッシュの格子点における温度が次第に低下して高温酸素析出物の形成開始温度Ｔpになったときに、次の近似式（１８）を用いて高温酸素析出物の密度Ｎpを求める。なお、式（１８）において、ａ2は定数である。
【００３６】
【数１４】
更に第１２’ステップとして、シリコン単結晶１４内のそれぞれのメッシュの格子点における温度が高温酸素析出物の形成開始温度Ｔpより低いときの高温酸素析出物の半径Ｒpを、次の式（１９）から求める。
【００３７】
【数１５】
ここで、ｔ2はシリコン単結晶１４のメッシュの格子点における温度が高温酸素析出物の形成開始温度Ｔpまで低下したときの時刻であり、Ｒcrは高温酸素析出物の臨界径である。
【００３８】
一方、上記式（１９）が成立つのは、空孔が析出するための流速が十分大きく、ＳｉマトリックスとＳｉＯ2との単位質量当りの体積差を埋められる場合、即ちＤv（Ｃv−Ｃv e）≧γＤ0Ｃ0の場合である。上記以外の場合には、次式（２０）が成立つ。
【００３９】
【数１６】
上記第８、第９、第１０’、第１１’及び第１２’ステップをシリコン単結晶１４の冷却が完了するまで繰返す（第１３’ステップ）。なお、上記式（９）〜式（１３）及び式（１７）〜式（２０）は連成してコンピュータにより解く。
【００４０】
上述のように、シリコン融液１２の対流を考慮してシリコン融液１２から成長するシリコン単結晶１４内の温度分布を求めた後に、シリコン融液１２から切離されたシリコン単結晶１４の冷却過程、即ちシリコン単結晶１４がシリコン融液１２から切離された後のシリコン単結晶１４の引上げ速度を考慮し、シリコン単結晶１４の徐冷及び急冷の効果を結果に反映して解析することにより、シリコン単結晶１４内のボイド及び高温酸素析出物の密度分布及びサイズ分布を正確に予測できる。この結果、シリコン融液１２から引上げられるシリコン単結晶１４内の欠陥の密度分布及びサイズ分布を所望の分布にするために、引上げ機１１の設計段階で引上げ機１１の構造を検討することができる。
なお、この実施の形態では、単結晶としてシリコン単結晶を挙げたが、ＧａＡｓ単結晶，ＩｎＰ単結晶，ＺｎＳ単結晶若しくはＺｎＳｅ単結晶でもよい。
【００４１】
【実施例】
次に本発明の実施例を比較例とともに詳しく説明する。
＜実施例１＞
図７に示すように、石英るつぼ１３に貯留されたシリコン融液１２から直径６インチのシリコン単結晶１４を引上げる場合の、シリコン単結晶１４内のボイドの密度分布及びサイズ陥分布を、図１〜図３のフローチャートに基づくシミュレーション方法により求め、更にシリコン単結晶１４内の高温酸素析出物の密度分布及びサイズ陥分布を、図４〜図６のフローチャートに基づくシミュレーション方法により求めた。
【００４２】
即ち、シリコン単結晶引上げ機１１のホットゾーンをメッシュ構造でモデル化した。ここで、シリコン融液１２のシリコン単結晶１４直下のシリコン単結晶１４の径方向のメッシュを０．７５ｍｍに設定し、シリコン融液１２のシリコン単結晶１４直下以外のシリコン単結晶１４の径方向のメッシュを１〜５ｍｍに設定した。またシリコン融液１２のシリコン単結晶１４の長手方向のメッシュを０．２５〜５ｍｍに設定し、乱流モデル式の乱流パラメータＣとして０．４５を用いた。更にシリコン単結晶１４の引上げ開始時ｔ0から冷却完了時ｔ1までの引上げ長及び引上げ高さの段階的な変更を５０ｍｍずつとした。このような条件下で、シリコン融液１２の対流を考慮してシリコン単結晶１４内の温度分布を求めた。
【００４３】
次にシリコン単結晶１４の引上げ開始時ｔ0から冷却完了時ｔ1までの時間を所定の間隔Δｔ秒（微小な時間の間隔）毎に区切り、この時間間隔Δｔ秒毎に上記シリコン単結晶１４のメッシュの座標及び温度のデータから、シリコン単結晶１４の引上げ長及び引上げ高さとシリコン単結晶内１４の温度分布を求め、更にボイドの形成開始温度を求めてボイドの密度分布及びサイズ分布を求めた。即ち、シリコン単結晶１４をシリコン融液１２から切離した後のシリコン単結晶１４の徐冷及び急冷を考慮して、シリコン単結晶１４内のボイドの密度分布及びサイズ陥分布をコンピュータを用いてそれぞれ求めた。その結果を図９（ａ）及び（ｂ）に示す。
【００４４】
一方、上記と同様に、シリコン単結晶１４の引上げ開始時ｔ0から冷却完了時ｔ1までの時間を所定の間隔Δｔ秒（微小な時間の間隔）毎に区切り、この時間間隔Δｔ秒毎に上記シリコン単結晶１４のメッシュの座標及び温度のデータから、シリコン単結晶１４の引上げ長及び引上げ高さとシリコン単結晶内１４の温度分布を求め、更に高温酸素析出物の形成開始温度を求めて高温酸素析出物の密度分布及びサイズ分布を求めた。即ち、シリコン単結晶１４をシリコン融液１２から切離した後のシリコン単結晶１４の徐冷及び急冷を考慮して、シリコン単結晶１４内の高温酸素析出物の密度分布及びサイズ陥分布をコンピュータを用いてそれぞれ求めた。その結果を図９（ｃ）及び（ｄ）に示す。
【００４５】
＜比較例１＞
実施例１と同一形状のシリコン単結晶を引上げたときにライフタイムを実験により求めた。その結果を図９（ｅ）に示す。
【００４６】
図９（ａ）及び図９（ｂ）から明らかなように、シリコン単結晶内のボイドの密度分布は下方に向うに従って多くなり、シリコン単結晶内のボイドのサイズ分布は上方に向うに従って大きくなった。
図９（ｃ）及び図９（ｄ）から明らかなように、シリコン単結晶内の高温酸素析出物の密度分布は略Ｖ字状に最も多くなる領域が現れ、シリコン単結晶内の高温酸素析出物の密度分布も略Ｖ字状に最も多くなる領域が現れた。
上記結果を図９（ｅ）の実際に測定した例と比較すると、実測例中のボイドなどが現れる領域と、図９（ａ）及び図９（ｂ）におけるボイドが存在する領域とが一致し、実測例中のＯＳＦリングが観察される領域と、図９（ｃ）及び図９（ｄ）でのＶ字状の領域とが一致することが判った。
【００４７】
【発明の効果】
以上述べたように、本発明によれば、融液の対流を考慮して融液から成長する単結晶内の温度分布をコンピュータを用いて求めるだけでなく、更に冷却過程における単結晶内の温度分布までも求めることによって、即ち融液から切離された単結晶の冷却過程における単結晶の徐冷及び急冷の効果を考慮することによって、単結晶内の欠陥の密度分布及びサイズ分布をコンピュータを用いて求める。この結果、単結晶内のボイド及び酸素析出物の密度分布及びサイズ分布を正確に予測できる。
【図面の簡単な説明】
【図１】本発明実施形態シリコン単結晶内のボイドの密度分布及びサイズ分布のシミュレーション方法の前半を示すフローチャート。
【図２】そのボイドの密度分布及びサイズ分布のシミュレーション方法の中盤を示すフローチャート。
【図３】そのボイドの密度分布及びサイズ分布のシミュレーション方法の後半を示すフローチャート。
【図４】シリコン単結晶内の高温酸素析出物の密度分布及びサイズ分布のシミュレーション方法の前半を示すフローチャート。
【図５】その高温酸素析出物の密度分布及びサイズ分布のシミュレーション方法の中盤を示すフローチャート。
【図６】その高温酸素析出物の密度分布及びサイズ分布のシミュレーション方法の後半を示すフローチャート。
【図７】本発明のシリコン融液をメッシュ構造としたシリコン単結晶の引上げ機の要部断面図。
【図８】従来例のシリコン融液をメッシュ構造としたシリコン単結晶の引上げ機の要部断面図。
【図９】実施例１のボイド及び高温酸素析出物の密度分布及びサイズ分布と、比較例１の実際に測定した例を示す図。
【符号の説明】
１１ シリコン単結晶引上げ機
１２ シリコン融液
１４ シリコン単結晶
Ｓ シリコンの三重点
【発明の属する技術分野】
本発明は、チョクラルスキー（以下、ＣＺという。）法にて引上げられるシリコン等の単結晶内欠陥の密度分布及びサイズ分布をコンピュータシミュレーションする方法に関するものである。
【０００２】
【従来の技術】
従来、この種のシミュレーション方法として、図８に示すように、総合伝熱シミュレータを用いてＣＺ法によるシリコン単結晶４引上げ時の引上げ機１内のホットゾーン構造及びそのシリコン単結晶４の引上げ速度に基づいて、シリコン融液２の熱伝導率を操作することによりシリコン融液２の内部温度分布を予測し、この内部温度分布からシリコン単結晶４のメッシュの座標及び温度をそれぞれ求め、更にシリコン単結晶４内の格子間シリコン及び空孔の拡散係数及び境界条件に基づいて拡散方程式を解くことにより、上記格子間シリコン及び空孔の密度分布をコンピュータを用いて求める方法が知られている。このシミュレーション方法では、ホットゾーンの各部材がメッシュ分割されてモデル化される。特にシリコン融液２のメッシュは計算時間を短くするために１０ｍｍ程度と比較的粗く設定される。
【０００３】
【発明が解決しようとする課題】
しかし、上記従来の格子間シリコン及び空孔の密度分布のシミュレーション方法では、実際の引上げ機においては発生するシリコン融液の対流を考慮しておらず、またシリコン融液のメッシュが比較的粗いため、固液界面形状の再現性が悪く、精度の良い結晶内温度分布を提供できない。そのため、シリコン単結晶内の格子間シリコン及び空孔の密度分布が実測値と大幅に相違する問題点があった。また上記従来のシミュレーション方法では、結晶冷却過程に比較的高温で生成されるシリコン単結晶内の欠陥、例えばボイドや酸素析出物の種類や、密度分布及びサイズ分布は判らなかった。
本発明の目的は、融液の対流を考慮して成長中の単結晶内の温度分布を解析した後に、融液から切離された単結晶の冷却過程を考慮して解析することにより、単結晶内のボイド及び高温酸素析出物の密度分布及びサイズ分布を正確に予測できる、単結晶内欠陥の密度分布及びサイズ分布のシミュレーション方法を提供することにある。
【０００４】
【課題を解決するための手段】
請求項１に係る発明は、図１〜図３及び図７に示すように、引上げ機１１による単結晶１４の融液１２からの引上げ開始時から単結晶１４の冷却完了時までの引上げ機１１のホットゾーンをメッシュ構造でモデル化する第１ステップと、ホットゾーンの各部材毎にメッシュをまとめかつこのまとめられたメッシュに対する各部材の物性値とともに単結晶１４の引上げ長及びこの引上げ長に対応する単結晶１４の引上げ速度をそれぞれコンピュータに入力する第２ステップと、各部材の表面温度分布をヒータの発熱量及び各部材の輻射率に基づいて求める第３ステップと、各部材の表面温度分布及び熱伝導率に基づいて熱伝導方程式を解くことにより各部材の内部温度分布を求めた後に融液１２が乱流であると仮定して得られた乱流モデル式及びナビエ・ストークスの方程式を連結して解くことにより対流を考慮した融液１２の内部温度分布を更に求める第４ステップと、単結晶１４及び融液１２の固液界面形状を単結晶の三重点Ｓを含む等温線に合せて求める第５ステップと、第３ステップから第５ステップを三重点Ｓが単結晶１４の融点になるまで繰返し引上げ機１１内の温度分布を計算して単結晶１４のメッシュの座標及び温度を求めこれらのデータをそれぞれコンピュータに入力する第６ステップと、単結晶１４の引上げ長及び引上げ高さを段階的に変えて第１ステップから第６ステップまでを繰返し引上げ機１１内の温度分布を計算して単結晶１４のメッシュの座標及び温度を求めこれらのデータをそれぞれコンピュータに入力する第７ステップと、単結晶１４の融液１２からの引上げ開始時から単結晶１４の冷却完了時までの時間を所定の間隔毎に区切りこの区切られた時間間隔毎に第７ステップで求めた単結晶１４のメッシュの座標及び温度のデータから単結晶１４の引上げ長及び引上げ高さと単結晶１４内の温度分布とを求める第８ステップと、単結晶１４内の空孔及び格子間原子の拡散係数及び境界条件に基づいて拡散方程式を解くことにより所定の時間間隔の経過した後の空孔及び格子間原子の密度分布を求める第９ステップと、空孔の密度分布に基づいてボイドの形成開始温度を求める第１０ステップと、単結晶１４内のそれぞれのメッシュの格子点における温度が次第に低下してボイドの形成開始温度になったときのボイドの密度を求める第１１ステップと、単結晶１４内のそれぞれのメッシュの格子点における温度がボイドの形成開始温度より低いときのボイドの半径を求める第１２ステップと、第８ステップから第１２ステップを単結晶１４の冷却が完了するまで繰返す第１３ステップとを含むコンピュータを用いて単結晶内欠陥の密度分布及びサイズ分布のシミュレーションを行う方法であって、乱流モデル式が次の式（２）で表されるｋｌ−モデル式であり、このモデル式の乱流パラメータＣとして０．４〜０．６の範囲内の任意の値が用いられたことを特徴とする、コンピュータを用いて単結晶内欠陥の密度分布及びサイズ分布のシミュレーションを行う方法である。
【数２】
【０００５】
この請求項１に記載されたシミュレーション方法では、融液１２の対流を考慮して融液１２から成長する単結晶１４内の温度分布を求めるだけでなく、更に冷却過程における単結晶１４内の温度分布までも求めることによって、即ち融液１２から切離された単結晶１４の冷却過程における単結晶１４の徐冷及び急冷の効果を考慮して解析することによって、単結晶１４内のボイドの密度分布及びサイズ分布を正確に予測できる。
【０００６】
請求項２に係る発明は、図４〜図７に示すように、請求項１記載の第１〜第９ステップと、空孔の密度分布に基づいて高温酸素析出物の形成開始温度を求める第１０’ステップと、単結晶１４内のそれぞれのメッシュの格子点における温度が次第に低下して高温酸素析出物の形成開始温度になったときの高温酸素析出物の密度を求める第１１’ステップと、単結晶１４内のそれぞれのメッシュの格子点における温度が高温酸素析出物の形成開始温度より低いときの高温酸素析出物の半径を求める第１２’ステップと、第８、第９、第１０’、第１１’及び第１２’ステップを単結晶１４の冷却が完了するまで繰返す第１３’ステップとを含むコンピュータを用いて単結晶内欠陥の密度分布及びサイズ分布のシミュレーションを行う方法であって、乱流モデル式が次の式（２）で表されるｋｌ−モデル式であり、このモデル式の乱流パラメータＣとして０．４〜０．６の範囲内の任意の値が用いられたことを特徴とする、コンピュータを用いて単結晶内欠陥の密度分布及びサイズ分布のシミュレーションを行う方法である。
【数２】
【０００７】
この請求項２に記載されたシミュレーション方法では、融液１２の対流を考慮して融液１２から成長する単結晶１４内の温度分布を求めるだけでなく、更に冷却過程における単結晶１４内の温度分布までも求めることによって、即ち融液１２から切離された単結晶１４の冷却過程における単結晶１４の徐冷及び急冷の効果を考慮して解析することによって、単結晶１４内の高温酸素析出物の密度分布及びサイズ分布を正確に予測できる。
なお、本明細書において「高温酸素析出物」とは、ボイド欠陥の内側に析出する内壁酸化膜を球形に近似して扱ったものをいう。また酸素析出物を高温と限定したのは、ボイド欠陥より低温で析出するＢＭＤ（Bulk Micro Defect）等の析出核となるものと区別するためである。
【０００８】
【発明の実施の形態】
次に本発明の実施の形態を図面に基づいて説明する。
図７に示すように、シリコン単結晶引上げ機１１のチャンバ内には、シリコン融液１２を貯留する石英るつぼ１３が設けられる。この石英るつぼ１３は図示しないが黒鉛サセプタ及び支軸を介してるつぼ駆動手段に接続され、るつぼ駆動手段は石英るつぼ１３を回転させるとともに昇降させるように構成される。また石英るつぼ１３の外周面は石英るつぼ１３から所定の間隔をあけてヒータ（図示せず）により包囲され、このヒータは保温筒（図示せず）により包囲される。ヒータは石英るつぼ１３に投入された高純度のシリコン多結晶体を加熱・溶融してシリコン融液１２にする。またチャンバの上端には図示しないが円筒状のケーシングが接続され、このケーシングには引上げ手段が設けられる。引上げ手段はシリコン単結晶１４を回転させながら引上げるように構成される。
【０００９】
このように構成されたシリコン単結晶引上げ機１１におけるシリコン単結晶１４内のボイドの密度分布及びサイズ分布のシミュレーション方法を図１〜図３及び図７に基づいて説明する。
先ず第１ステップとして、シリコン単結晶１４を所定長さＬ1（例えば１００ｍｍ）まで引上げた状態におけるシリコン単結晶引上げ機１１のホットゾーンの各部材、即ちチャンバ，石英るつぼ１３，シリコン融液１２，シリコン単結晶１４，黒鉛サセプタ，保温筒等をメッシュ分割してモデル化する。具体的には上記ホットゾーンの各部材のメッシュ点の座標データをコンピュータに入力する。このときシリコン融液１２のメッシュのうちシリコン単結晶１４の径方向のメッシュであってかつシリコン融液１２のシリコン単結晶１４直下の一部又は全部のメッシュ（以下、径方向メッシュという。）を０．０１〜５．００ｍｍ、好ましくは０．２５〜１．００ｍｍに設定する。またシリコン融液１２のメッシュのうちシリコン単結晶１４の長手方向のメッシュであってかつシリコン融液１２の一部又は全部のメッシュ（以下、長手方向メッシュという。）を０．０１〜５．００ｍｍ、好ましくは０．１〜０．５ｍｍに設定する。
【００１０】
径方向メッシュを０．０１〜５．００ｍｍの範囲に限定したのは、０．０１ｍｍ未満では計算時間が極めて長くなり、５．００ｍｍを越えると計算が不安定になり、繰返し計算を行っても固液界面形状が一定に定まらなくなるからである。また長手方向メッシュを０．０１〜５．００ｍｍの範囲に限定したのは、０．０１ｍｍ未満では計算時間が極めて長くなり、５．００ｍｍを越えると固液界面形状の計算値が実測値と一致しなくなるからである。なお、径方向メッシュの一部を０．０１〜５．００の範囲に限定する場合には、シリコン単結晶１４直下のシリコン融液１２のうちシリコン単結晶１４外周縁近傍のシリコン融液１２を上記範囲に限定することが好ましく、長手方向メッシュの一部を０．０１〜５．００の範囲に限定する場合には、シリコン融液１２の液面近傍及び底近傍を上記範囲に限定することが好ましい。
【００１１】
第２ステップとして上記ホットゾーンの各部材毎にメッシュをまとめ、かつこのまとめられたメッシュに対して各部材の物性値をそれぞれコンピュータに入力する。例えば、チャンバがステンレス鋼にて形成されていれば、そのステンレス鋼の熱伝導率，輻射率，粘性率，体積膨張係数，密度及び比熱がコンピュータに入力される。またシリコン単結晶１４の引上げ長及びこの引上げ長に対応するシリコン単結晶１４の引上げ速度と、後述する乱流モデル式（１）の乱流パラメータＣとをコンピュータに入力する。
【００１２】
第３ステップとして、ホットゾーンの各部材の表面温度分布をヒータの発熱量及び各部材の輻射率に基づいてコンピュータを用いて求める。即ち、ヒータの発熱量を任意に設定してコンピュータに入力するとともに、各部材の輻射率から各部材の表面温度分布をコンピュータを用いて求める。次に第４ステップとしてホットゾーンの各部材の表面温度分布及び熱伝導率に基づいて熱伝導方程式（１）をコンピュータを用いて解くことにより各部材の内部温度分布を求める。ここでは、記述を簡単にするためｘｙｚ直交座標系を用いたが、実際の計算では円筒座標系を用いる。
【００１３】
【数１】
一方、シリコン融液１２に関しては、上記熱伝導方程式（１）でシリコン融液１２の内部温度分布を求めた後に、このシリコン融液１２の内部温度分布に基づき、シリコン融液１２が乱流であると仮定して得られた乱流モデル式（２）及びナビエ・ストークスの方程式（３）〜（５）を連結して、シリコン融液１２の内部流速分布をコンピュータを用いて求める。
【００１４】
【数２】
【００１５】
【数３】
【００１６】
上記乱流モデル式（２）はｋｌ（ケイエル）−モデル式と呼ばれ、このモデル式の乱流パラメータＣは０．４〜０．６の範囲内の任意の値が用いられることが好ましい。乱流パラメータＣを０．４〜０．６の範囲に限定したのは、０．４未満又は０．６を越えると計算により求めた界面形状が実測値と一致しないという不具合があるからである。また上記ナビエ・ストークスの方程式（３）〜（５）はシリコン融液１２が非圧縮性であって粘度が一定である流体としたときの運動方程式である。
上記求められたシリコン融液１２の内部流速分布に基づいて熱エネルギ方程式（６）を解くことにより、シリコン融液１２の対流を考慮したシリコン融液１２の内部温度分布をコンピュータを用いて更に求める。
【００１７】
【数４】
【００１８】
次いで第５ステップとして、シリコン単結晶１４及びシリコン融液１２の固液界面形状を図２の点Ｓで示すシリコンの三重点Ｓ（固体と液体と気体の三重点(tri-junction)）を含む等温線に合せてコンピュータを用いて求める。第６ステップとして、コンピュータに入力するヒータの発熱量を変更し（次第に増大し）、上記第３ステップから第５ステップを三重点がシリコン単結晶１４の融点になるまで繰返した後に、引上げ機１１内の温度分布を計算してシリコン単結晶のメッシュの座標及び温度を求め、これらのデータをコンピュータに記憶させる。
【００１９】
次に第７ステップとして、シリコン単結晶１４の引上げ長Ｌ1にδ（例えば５０ｍｍ）だけ加えて上記第１ステップから第６ステップまでを繰返した後に、引上げ機１１内の温度分布を計算してシリコン単結晶１４のメッシュの座標及び温度を求め、これらのデータをコンピュータに記憶させる。この第７ステップは、シリコン単結晶１４の引上げ長Ｌ1が長さＬ2（Ｌ2はシリコン融液１２から切離されたときのシリコン単結晶１４の長さ（成長完了時の結晶長）である。）に達してシリコン単結晶１４がシリコン融液１２から切離された後、更にシリコン単結晶１４が引上げられてその高さＨ1（Ｈ1はシリコン単結晶１４の直胴開始部からシリコン融液１２の液面までの距離である（図７）。）がＨ2（Ｈ2は冷却完了時のシリコン単結晶１４の直胴開始部からシリコン融液１２の液面までの距離である。）に達するまで、即ちシリコン単結晶１４の冷却が完了するまで行われる。なお、シリコン単結晶１４がシリコン融液１２から切離された後は、シリコン単結晶１４の引上げ高さＨ1にδ（例えば５０ｍｍ）だけ加え、上記と同様に上記第１ステップから第６ステップまでを繰返す。
【００２０】
シリコン単結晶１４の引上げ高さＨ1がＨ2に達すると、第８ステップに移行する。第８ステップでは、シリコン単結晶１４をシリコン融液１２から成長させて引上げ始めたときｔ0から、シリコン単結晶１４をシリコン融液１２から切離して更にシリコン単結晶１４を引上げ、その冷却が完了したときｔ1までの時間を、所定の間隔Δｔ秒（微小時間間隔）毎に区切る。このときシリコン単結晶１４内の格子間シリコン及び空孔の拡散係数及び境界条件のみならず、後述するボイド及び高温酸素析出物の密度分布及びサイズ分布を求めるための式に用いられる定数をそれぞれコンピュータに入力する。上記区切られた時間間隔Δｔ秒毎に、第７ステップで求めたシリコン単結晶１４のメッシュの座標及び温度のデータから、シリコン単結晶１４の引上げ長Ｌ1及び引上げ高さＨ1と、シリコン単結晶１４内の温度分布とを求める。
【００２１】
即ち、第１〜第７ステップでシリコン単結晶のメッシュの座標及び温度を引上げ長δ毎に求め、シリコン単結晶を例えば５０ｍｍ引上げるのに数十分要するため、この数十分間でのシリコン単結晶のメッシュの温度変化を時間の関数として微分することにより、時刻ｔ0からΔｔ秒後におけるシリコン単結晶１４の引上げ長Ｌ1及び引上げ高さＨ1とシリコン単結晶１４内の温度分布を算出する。次にシリコン単結晶１４内の空孔及び格子間シリコンの拡散係数及び境界条件に基づいて拡散方程式を解くことにより、Δｔ秒経過後の空孔及び格子間シリコンの密度分布を求める（第９ステップ）。
具体的には、空孔の密度Ｃvの計算式が次の式（７）で示され、格子間シリコンの密度Ｃiの計算式が次の式（８）で示される。式（７）及び式（８）において、密度Ｃv及び密度Ｃiの経時的進展を計算するために、空孔と格子間シリコンの熱平衡がシリコン単結晶の全表面で維持されると仮定する。
【００２２】
【数５】
上記平衡式は時間で微分され、空孔及び格子間シリコンに対してそれぞれ次の式（９）及び式（１０）になる。
【００２３】
【数６】
ここで、Θ(ｘ)はヘビサイド関数（Heaviside function）である。即ち、ｘ＜０のときΘ(ｘ)＝０であり、かつｘ＞０のときΘ(ｘ)＝１である。またＴnは高温酸素析出物の形成開始温度Ｔpとボイドの形成開始温度Ｔv0とを比較したときの高い方の温度である。更に式（９）及び（１０）のそれぞれ右側の第１項はフィックの拡散式であり、右側の第１項中のＤv及びＤiは、次の式（１１）及び（１２）に表される拡散係数である。
【００２５】
【数７】
ここで△Ｅv及び△Ｅiはそれぞれ空孔及び格子間シリコンの活性化エネルギーであり、ｄv及びｄiはそれぞれ定数である。また式（９）及び式（１０）のそれぞれ右側の第２項中の
【数８】
【００２７】
一方、上記式（９）が成立つのは、空孔が析出するための流速が十分大きく、ＳｉマトリックスとＳｉＯ2との単位質量当りの体積差を埋められる場合、即ちＤv（Ｃv−Ｃv e）≧γＤ0Ｃ0の場合である。上記以外の場合には、次式（１３）が成立つ。
【００２８】
【数９】
【００２９】
【数１０】
【００３０】
第１１ステップとして、シリコン単結晶１４内のそれぞれのメッシュの格子点における温度が次第に低下してボイドの形成開始温度Ｔv0になったときに、次の近似式（１５）を用いてボイドの密度Ｎv0を求める。
【００３１】
【数１１】
【００３２】
【数１２】
【００３３】
一方、上記シリコン単結晶引上げ機１１におけるシリコン単結晶１４内の高温酸素析出物の密度分布及びサイズ分布のシミュレーション方法を図４〜図７に基づいて説明する。
第１〜第９ステップまでは、上記ボイドの密度分布及びサイズ分布を求める場合と同一である。
次いで第１０’ステップとして、拡散方程式（９）及び（１０）を解くことにより求めた空孔の密度Ｃv分布に基づいて、高温酸素析出物の形成開始温度Ｔpを次の式（１７）から求める。
【００３４】
【数１３】
【００３５】
次に第１１’ステップとして、シリコン単結晶１４内のそれぞれのメッシュの格子点における温度が次第に低下して高温酸素析出物の形成開始温度Ｔpになったときに、次の近似式（１８）を用いて高温酸素析出物の密度Ｎpを求める。なお、式（１８）において、ａ2は定数である。
【００３６】
【数１４】
【００３７】
【数１５】
【００３８】
一方、上記式（１９）が成立つのは、空孔が析出するための流速が十分大きく、ＳｉマトリックスとＳｉＯ2との単位質量当りの体積差を埋められる場合、即ちＤv（Ｃv−Ｃv e）≧γＤ0Ｃ0の場合である。上記以外の場合には、次式（２０）が成立つ。
【００３９】
【数１６】
【００４０】
上述のように、シリコン融液１２の対流を考慮してシリコン融液１２から成長するシリコン単結晶１４内の温度分布を求めた後に、シリコン融液１２から切離されたシリコン単結晶１４の冷却過程、即ちシリコン単結晶１４がシリコン融液１２から切離された後のシリコン単結晶１４の引上げ速度を考慮し、シリコン単結晶１４の徐冷及び急冷の効果を結果に反映して解析することにより、シリコン単結晶１４内のボイド及び高温酸素析出物の密度分布及びサイズ分布を正確に予測できる。この結果、シリコン融液１２から引上げられるシリコン単結晶１４内の欠陥の密度分布及びサイズ分布を所望の分布にするために、引上げ機１１の設計段階で引上げ機１１の構造を検討することができる。
なお、この実施の形態では、単結晶としてシリコン単結晶を挙げたが、ＧａＡｓ単結晶，ＩｎＰ単結晶，ＺｎＳ単結晶若しくはＺｎＳｅ単結晶でもよい。
【００４１】
【実施例】
次に本発明の実施例を比較例とともに詳しく説明する。
＜実施例１＞
図７に示すように、石英るつぼ１３に貯留されたシリコン融液１２から直径６インチのシリコン単結晶１４を引上げる場合の、シリコン単結晶１４内のボイドの密度分布及びサイズ陥分布を、図１〜図３のフローチャートに基づくシミュレーション方法により求め、更にシリコン単結晶１４内の高温酸素析出物の密度分布及びサイズ陥分布を、図４〜図６のフローチャートに基づくシミュレーション方法により求めた。
【００４２】
即ち、シリコン単結晶引上げ機１１のホットゾーンをメッシュ構造でモデル化した。ここで、シリコン融液１２のシリコン単結晶１４直下のシリコン単結晶１４の径方向のメッシュを０．７５ｍｍに設定し、シリコン融液１２のシリコン単結晶１４直下以外のシリコン単結晶１４の径方向のメッシュを１〜５ｍｍに設定した。またシリコン融液１２のシリコン単結晶１４の長手方向のメッシュを０．２５〜５ｍｍに設定し、乱流モデル式の乱流パラメータＣとして０．４５を用いた。更にシリコン単結晶１４の引上げ開始時ｔ0から冷却完了時ｔ1までの引上げ長及び引上げ高さの段階的な変更を５０ｍｍずつとした。このような条件下で、シリコン融液１２の対流を考慮してシリコン単結晶１４内の温度分布を求めた。
【００４３】
次にシリコン単結晶１４の引上げ開始時ｔ0から冷却完了時ｔ1までの時間を所定の間隔Δｔ秒（微小な時間の間隔）毎に区切り、この時間間隔Δｔ秒毎に上記シリコン単結晶１４のメッシュの座標及び温度のデータから、シリコン単結晶１４の引上げ長及び引上げ高さとシリコン単結晶内１４の温度分布を求め、更にボイドの形成開始温度を求めてボイドの密度分布及びサイズ分布を求めた。即ち、シリコン単結晶１４をシリコン融液１２から切離した後のシリコン単結晶１４の徐冷及び急冷を考慮して、シリコン単結晶１４内のボイドの密度分布及びサイズ陥分布をコンピュータを用いてそれぞれ求めた。その結果を図９（ａ）及び（ｂ）に示す。
【００４４】
一方、上記と同様に、シリコン単結晶１４の引上げ開始時ｔ0から冷却完了時ｔ1までの時間を所定の間隔Δｔ秒（微小な時間の間隔）毎に区切り、この時間間隔Δｔ秒毎に上記シリコン単結晶１４のメッシュの座標及び温度のデータから、シリコン単結晶１４の引上げ長及び引上げ高さとシリコン単結晶内１４の温度分布を求め、更に高温酸素析出物の形成開始温度を求めて高温酸素析出物の密度分布及びサイズ分布を求めた。即ち、シリコン単結晶１４をシリコン融液１２から切離した後のシリコン単結晶１４の徐冷及び急冷を考慮して、シリコン単結晶１４内の高温酸素析出物の密度分布及びサイズ陥分布をコンピュータを用いてそれぞれ求めた。その結果を図９（ｃ）及び（ｄ）に示す。
【００４５】
＜比較例１＞
実施例１と同一形状のシリコン単結晶を引上げたときにライフタイムを実験により求めた。その結果を図９（ｅ）に示す。
【００４６】
図９（ａ）及び図９（ｂ）から明らかなように、シリコン単結晶内のボイドの密度分布は下方に向うに従って多くなり、シリコン単結晶内のボイドのサイズ分布は上方に向うに従って大きくなった。
図９（ｃ）及び図９（ｄ）から明らかなように、シリコン単結晶内の高温酸素析出物の密度分布は略Ｖ字状に最も多くなる領域が現れ、シリコン単結晶内の高温酸素析出物の密度分布も略Ｖ字状に最も多くなる領域が現れた。
上記結果を図９（ｅ）の実際に測定した例と比較すると、実測例中のボイドなどが現れる領域と、図９（ａ）及び図９（ｂ）におけるボイドが存在する領域とが一致し、実測例中のＯＳＦリングが観察される領域と、図９（ｃ）及び図９（ｄ）でのＶ字状の領域とが一致することが判った。
【００４７】
【発明の効果】
以上述べたように、本発明によれば、融液の対流を考慮して融液から成長する単結晶内の温度分布をコンピュータを用いて求めるだけでなく、更に冷却過程における単結晶内の温度分布までも求めることによって、即ち融液から切離された単結晶の冷却過程における単結晶の徐冷及び急冷の効果を考慮することによって、単結晶内の欠陥の密度分布及びサイズ分布をコンピュータを用いて求める。この結果、単結晶内のボイド及び酸素析出物の密度分布及びサイズ分布を正確に予測できる。
【図面の簡単な説明】
【図１】本発明実施形態シリコン単結晶内のボイドの密度分布及びサイズ分布のシミュレーション方法の前半を示すフローチャート。
【図２】そのボイドの密度分布及びサイズ分布のシミュレーション方法の中盤を示すフローチャート。
【図３】そのボイドの密度分布及びサイズ分布のシミュレーション方法の後半を示すフローチャート。
【図４】シリコン単結晶内の高温酸素析出物の密度分布及びサイズ分布のシミュレーション方法の前半を示すフローチャート。
【図５】その高温酸素析出物の密度分布及びサイズ分布のシミュレーション方法の中盤を示すフローチャート。
【図６】その高温酸素析出物の密度分布及びサイズ分布のシミュレーション方法の後半を示すフローチャート。
【図７】本発明のシリコン融液をメッシュ構造としたシリコン単結晶の引上げ機の要部断面図。
【図８】従来例のシリコン融液をメッシュ構造としたシリコン単結晶の引上げ機の要部断面図。
【図９】実施例１のボイド及び高温酸素析出物の密度分布及びサイズ分布と、比較例１の実際に測定した例を示す図。
【符号の説明】
１１ シリコン単結晶引上げ機
１２ シリコン融液
１４ シリコン単結晶
Ｓ シリコンの三重点
Claims (2)
- 引上げ機(11)による単結晶(14)の融液(12)からの引上げ開始時から前記単結晶(14)の冷却完了時までの前記引上げ機(11)のホットゾーンをメッシュ構造でモデル化する第１ステップと、
前記ホットゾーンの各部材毎にメッシュをまとめかつこのまとめられたメッシュに対する前記各部材の物性値とともに前記単結晶(14)の引上げ長及びこの引上げ長に対応する前記単結晶(14)の引上げ速度をそれぞれコンピュータに入力する第２ステップと、
前記各部材の表面温度分布をヒータの発熱量及び前記各部材の輻射率に基づいて求める第３ステップと、
前記各部材の表面温度分布及び熱伝導率に基づいて熱伝導方程式を解くことにより前記各部材の内部温度分布を求めた後に融液(12)が乱流であると仮定して得られた乱流モデル式及びナビエ・ストークスの方程式を連結して解くことにより対流を考慮した前記融液(12)の内部温度分布を更に求める第４ステップと、
前記単結晶(14)及び前記融液(12)の固液界面形状を前記単結晶の三重点(S)を含む等温線に合せて求める第５ステップと、
前記第３ステップから前記第５ステップを前記三重点(S)が前記単結晶(14)の融点になるまで繰返し前記引上げ機(11)内の温度分布を計算して前記単結晶(14)のメッシュの座標及び温度を求めこれらのデータをそれぞれ前記コンピュータに入力する第６ステップと、
前記単結晶(14)の引上げ長及び引上げ高さを段階的に変えて前記第１ステップから前記第６ステップまでを繰返し前記引上げ機(11)内の温度分布を計算して前記単結晶(14)のメッシュの座標及び温度を求めこれらのデータをそれぞれ前記コンピュータに入力する第７ステップと、
前記単結晶(14)の前記融液(12)からの引上げ開始時から前記単結晶(14)の冷却完了時までの時間を所定の間隔毎に区切り前記区切られた時間間隔毎に第７ステップで求めた前記単結晶(14)のメッシュの座標及び温度のデータから前記単結晶(14)の引上げ長及び引上げ高さと前記単結晶(14)内の温度分布とを求める第８ステップと、
前記単結晶(14)内の空孔及び格子間原子の拡散係数及び境界条件に基づいて拡散方程式を解くことにより前記所定の時間間隔の経過した後の空孔及び格子間原子の密度分布を求める第９ステップと、
前記空孔の密度分布に基づいてボイドの形成開始温度を求める第１０ステップと、
前記単結晶(14)内のそれぞれのメッシュの格子点における温度が次第に低下して前記ボイドの形成開始温度になったときの前記ボイドの密度を求める第１１ステップと、
前記単結晶(14)内のそれぞれのメッシュの格子点における温度が前記ボイドの形成開始温度より低いときの前記ボイドの半径を求める第１２ステップと、
第８ステップから第１２ステップを前記単結晶(14)の冷却が完了するまで繰返す第１３ステップと
を含むコンピュータを用いて単結晶内欠陥の密度分布及びサイズ分布のシミュレーションを行う方法であって、
前記乱流モデル式が次の式（２）で表されるｋｌ−モデル式であり、このモデル式の乱流パラメータＣとして０．４〜０．６の範囲内の任意の値が用いられたことを特徴とする、コンピュータを用いて単結晶内欠陥の密度分布及びサイズ分布のシミュレーションを行う方法。
- 引上げ機(11)による単結晶(14)の融液(12)からの引上げ開始時から前記単結晶(14)の冷却完了時までの前記引上げ機(11)のホットゾーンをメッシュ構造でモデル化する第１ステップと、
前記ホットゾーンの各部材毎にメッシュをまとめかつこのまとめられたメッシュに対する前記各部材の物性値とともに前記単結晶(14)の引上げ長及びこの引上げ長に対応する前記単結晶(14)の引上げ速度をそれぞれコンピュータに入力する第２ステップと、
前記各部材の表面温度分布をヒータの発熱量及び前記各部材の輻射率に基づいて求める第３ステップと、
前記各部材の表面温度分布及び熱伝導率に基づいて熱伝導方程式を解くことにより前記各部材の内部温度分布を求めた後に融液(12)が乱流であると仮定して得られた乱流モデル式及びナビエ・ストークスの方程式を連結して解くことにより対流を考慮した前記融液(12)の内部温度分布を更に求める第４ステップと、
前記単結晶(14)及び前記融液(12)の固液界面形状を前記単結晶の三重点(S)を含む等温線に合せて求める第５ステップと、
前記第３ステップから前記第５ステップを前記三重点(S)が前記単結晶(14)の融点になるまで繰返し前記引上げ機(11)内の温度分布を計算して前記単結晶(14)のメッシュの座標及び温度を求めこれらのデータをそれぞれ前記コンピュータに入力する第６ステップと、
前記単結晶(14)の引上げ長及び引上げ高さを段階的に変えて前記第１ステップから前記第６ステップまでを繰返し前記引上げ機(11)内の温度分布を計算して前記単結晶(14)のメッシュの座標及び温度を求めこれらのデータをそれぞれ前記コンピュータに入力する第７ステップと、
前記単結晶(14)の前記融液(12)からの引上げ開始時から前記単結晶(14)の冷却完了時までの時間を所定の間隔毎に区切り前記区切られた時間間隔毎に第７ステップで求めた前記単結晶(14)のメッシュの座標及び温度のデータから前記単結晶(14)の引上げ長及び引上げ高さと前記単結晶(14)内の温度分布とを求める第８ステップと、
前記単結晶(14)内の空孔及び格子間原子の拡散係数及び境界条件に基づいて拡散方程式を解くことにより前記所定の時間間隔の経過した後の空孔及び格子間原子の密度分布を求める第９ステップと、
前記空孔の密度分布に基づいて高温酸素析出物の形成開始温度を求める第１０’ステップと、
前記単結晶(14)内のそれぞれのメッシュの格子点における温度が次第に低下して前記高温酸素析出物の形成開始温度になったときの前記高温酸素析出物の密度を求める第１１’ステップと、
前記単結晶(14)内のそれぞれのメッシュの格子点における温度が前記高温酸素析出物の形成開始温度より低いときの前記高温酸素析出物の半径を求める第１２’ステップと、
第８、第９、第１０’、第１１’及び第１２’ステップを前記単結晶(14)の冷却が完了するまで繰返す第１３’ステップと
を含むコンピュータを用いて単結晶内欠陥の密度分布及びサイズ分布のシミュレーションを行う方法であって、
前記乱流モデル式が次の式（２）で表されるｋｌ−モデル式であり、このモデル式の乱流パラメータＣとして０．４〜０．６の範囲内の任意の値が用いられたことを特徴とする、コンピュータを用いて単結晶内欠陥の密度分布及びサイズ分布のシミュレーションを行う方法。
Priority Applications (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2001227192A JP4106880B2 (ja) | 2001-07-27 | 2001-07-27 | 単結晶内欠陥の密度分布及びサイズ分布のシミュレーション方法 |
Applications Claiming Priority (1)
| Application Number | Priority Date | Filing Date | Title |
|---|---|---|---|
| JP2001227192A JP4106880B2 (ja) | 2001-07-27 | 2001-07-27 | 単結晶内欠陥の密度分布及びサイズ分布のシミュレーション方法 |
Publications (2)
| Publication Number | Publication Date |
|---|---|
| JP2003040695A JP2003040695A (ja) | 2003-02-13 |
| JP4106880B2 true JP4106880B2 (ja) | 2008-06-25 |
Family
ID=19059894
Family Applications (1)
| Application Number | Title | Priority Date | Filing Date |
|---|---|---|---|
| JP2001227192A Expired - Lifetime JP4106880B2 (ja) | 2001-07-27 | 2001-07-27 | 単結晶内欠陥の密度分布及びサイズ分布のシミュレーション方法 |
Country Status (1)
| Country | Link |
|---|---|
| JP (1) | JP4106880B2 (ja) |
Families Citing this family (4)
| Publication number | Priority date | Publication date | Assignee | Title |
|---|---|---|---|---|
| US7282094B2 (en) | 2003-05-28 | 2007-10-16 | Sumco Corporation | Method of simulation with respect to density distribution and size distribution of void defect within single crystal and oxygen precipitation nucleus within single crystal |
| JP4604462B2 (ja) * | 2003-06-06 | 2011-01-05 | 株式会社Ｓｕｍｃｏ | 単結晶内酸素析出核の密度分布及びサイズ分布のシミュレーション方法 |
| JP6052188B2 (ja) | 2014-01-16 | 2016-12-27 | 信越半導体株式会社 | シリコン単結晶ウェーハの熱処理方法 |
| JP6135611B2 (ja) * | 2014-07-03 | 2017-05-31 | 信越半導体株式会社 | 点欠陥濃度計算方法、Ｇｒｏｗｎ−ｉｎ欠陥計算方法、Ｇｒｏｗｎ−ｉｎ欠陥面内分布計算方法及びこれらを用いたシリコン単結晶製造方法 |
-
2001
- 2001-07-27 JP JP2001227192A patent/JP4106880B2/ja not_active Expired - Lifetime
Also Published As
| Publication number | Publication date |
|---|---|
| JP2003040695A (ja) | 2003-02-13 |
Similar Documents
| Publication | Publication Date | Title |
|---|---|---|
| JP6604338B2 (ja) | シリコン単結晶の引き上げ条件演算プログラム、シリコン単結晶のホットゾーンの改良方法、およびシリコン単結晶の育成方法 | |
| Kinney et al. | Application of turbulence modeling to the integrated hydrodynamic thermal-capillary model of Czochralski crystal growth of silicon | |
| JP4380537B2 (ja) | シリコン単結晶を製造する方法 | |
| JP2009190926A (ja) | シリコン単結晶製造における数値解析方法 | |
| KR100411553B1 (ko) | 단결정과 용융액 사이의 고체-액체 계면 형상 및 단결정의점결함 분포의 시뮬레이션 방법 | |
| JP4106880B2 (ja) | 単結晶内欠陥の密度分布及びサイズ分布のシミュレーション方法 | |
| Ding et al. | The influence mechanism of melt flow instability on the temperature fluctuation on the crystal/melt interface during Czochralski silicon crystal growth | |
| Dezfoli | Design and optimization of complex single heater for vertical gradient freeze (VGF) grower | |
| KR100719207B1 (ko) | 단결정 내 보이드 결함과 단결정 내 산소 석출핵의 밀도분포 및 사이즈 분포의 시뮬레이션 방법 | |
| JP4604462B2 (ja) | 単結晶内酸素析出核の密度分布及びサイズ分布のシミュレーション方法 | |
| JP4403722B2 (ja) | シリコン単結晶内ボイド欠陥の密度分布及びサイズ分布のシミュレーション方法 | |
| JP3846155B2 (ja) | 単結晶及び融液の固液界面形状のシミュレーション方法 | |
| JP4096499B2 (ja) | 単結晶の点欠陥分布のシミュレーション方法 | |
| JP4449347B2 (ja) | シミュレーションによるｏｓｆリングの分布予測方法 | |
| Derby | Theoretical modeling of Czochralski crystal growth | |
| KR100665683B1 (ko) | 실리콘 단결정 제조방법 | |
| Kitamura et al. | Radial distribution of temperature gradients in growing CZ-Si crystals and its application to the prediction of microdefect distribution | |
| Podkopaev et al. | Multiphysical simulation analysis of the dislocation structure in germanium single crystals | |
| Voigt et al. | Transient and quasi‐stationary simulation of heat and mass transfer in Czochralski silicon crystal growth | |
| JP4192704B2 (ja) | 単結晶の無欠陥領域を最大化するシミュレーション方法 | |
| Dupret et al. | Mathematical modeling of the growth of large diameter Czochralski silicon crystals considering melt dynamics | |
| Voigt et al. | Controlling Point Defects in Single Silicon Crystals Grown by the Czochralski Method |
Legal Events
| Date | Code | Title | Description |
|---|---|---|---|
| A625 | Written request for application examination (by other person) |
Free format text: JAPANESE INTERMEDIATE CODE: A625 Effective date: 20050408 |
|
| A977 | Report on retrieval |
Free format text: JAPANESE INTERMEDIATE CODE: A971007 Effective date: 20071217 |
|
| A131 | Notification of reasons for refusal |
Free format text: JAPANESE INTERMEDIATE CODE: A131 Effective date: 20071225 |
|
| A521 | Written amendment |
Free format text: JAPANESE INTERMEDIATE CODE: A523 Effective date: 20080208 |
|
| TRDD | Decision of grant or rejection written | ||
| A01 | Written decision to grant a patent or to grant a registration (utility model) |
Free format text: JAPANESE INTERMEDIATE CODE: A01 Effective date: 20080311 |
|
| A61 | First payment of annual fees (during grant procedure) |
Free format text: JAPANESE INTERMEDIATE CODE: A61 Effective date: 20080324 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20110411 Year of fee payment: 3 |
|
| R150 | Certificate of patent or registration of utility model |
Ref document number: 4106880 Country of ref document: JP Free format text: JAPANESE INTERMEDIATE CODE: R150 Free format text: JAPANESE INTERMEDIATE CODE: R150 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20120411 Year of fee payment: 4 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20130411 Year of fee payment: 5 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| FPAY | Renewal fee payment (event date is renewal date of database) |
Free format text: PAYMENT UNTIL: 20140411 Year of fee payment: 6 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| R250 | Receipt of annual fees |
Free format text: JAPANESE INTERMEDIATE CODE: R250 |
|
| EXPY | Cancellation because of completion of term |