【文献】
宮下大 外5名,3次元ハイプリッドPIC/MC法によるRPD装置内プラズマシミュレーション,第59回応用物理学関係連合講演会講演予稿集,2012年,p.08-053
【文献】
田口俊弘,高速点火核融合の統合シミュレーション −多階層プラズマシミュレーションシステム"FI3"− 3.計算コードと解析例 3.3 高速電子の伝播(ハイブリッドコード),プラズマ・核融合学会誌,社団法人プラズマ・核融合学会,2006年 3月25日,第82巻,第3号,pp.150−157
(58)【調査した分野】(Int.Cl.,DB名)
前記第2ステップにおける計算領域は、前記第1又は第3ステップにおける計算領域を、各座標軸において等距離拡大するように設定されてなる、請求項1又は2記載のプラズマシミュレーション方法。
前記第2ステップにおける計算領域は、前記第1又は第3ステップにおける計算領域と重複する領域と、前記第1又は第3ステップにおける計算領域と重複しない領域と、を含み、
前記重複しない領域の計算点の密度は、前記重複する領域の計算点の密度よりも小さい、請求項1〜3の何れか一項記載のプラズマシミュレーション方法。
【発明の概要】
【発明が解決しようとする課題】
【0004】
ここで、上記プラズマシミュレーション方法では、流体モデルを用いて計算を行う際、自己誘導磁場の計算時において、無限遠の境界条件を再現するために計算領域を大きくする場合がある。しかし、このように計算領域を大きくすると、例えば膨大な規模の格子について電子流体方程式を取り扱うことになることから、計算が冗長となり計算時間かかり過ぎる虞がある。
【0005】
本発明は上記実情に鑑みてなされてものであり、計算コストを低減することができるプラズマシミュレーション方法及びプラズマシミュレーションプログラムを提供することを課題とする。
【課題を解決するための手段】
【0006】
上記課題を解決するために、本発明に係るプラズマシミュレーション方法は、粒子モデル及び流体モデルを用いたプラズマシミュレーション方法であって、粒子モデルを用いて、イオン及び中性粒子に関する計算を行う粒子モデル計算ステップと、流体モデルを用いて、電子に関する計算を行う流体モデル計算ステップと、を備え、流体モデル計算ステップは、電子流体方程式を用いて電子温度及び電位を求める第1ステップと、自己誘導磁場を求める第2ステップと、移動度を更新する第3ステップと、を有し、第2ステップにおける計算領域は、第1及び第3ステップにおける計算領域よりも大きい。
【0007】
このプラズマシミュレーション方法では、電子温度及び電位を求める場合及び移動度を更新する場合には、例えば所定領域を計算領域として電子流体方程式を計算し、自己誘導磁場を求める場合には、当該所定領域よりも大きい領域を計算領域として計算することが可能となる。すなわち、自己誘導磁場の計算領域を大きくして計算精度を維持しつつ、当該計算領域よりも小さい計算領域にて電子温度及び電位の算出並びに移動度の更新を行い、その計算を簡潔化することができる。従って、本発明によれば、計算コストを低減することが可能となる。
【0008】
また、上記作用効果を好適に奏する場合として、具体的には、第1及び第3ステップにおける計算領域は、シミュレートの対象となる所定領域とされている場合がある。
【0009】
また、第2ステップにおける計算領域は、第1又は第3ステップにおける計算領域を、各座標軸において等距離拡大するように設定されてなることが好ましい。この場合、シミュレーションの精度を高めることが可能となる。
【0010】
また、第2ステップにおける計算領域は、第1又は第3ステップにおける計算領域と重複する領域と、第1又は第3ステップにおける計算領域と重複しない領域と、を含み、重複しない領域の計算点の密度は、重複する領域の計算点の密度よりも小さいことが好ましい。このように、重複しない領域の計算点の密度を重複する領域の計算点の密度よりも小さくすることにより、シミュレーションの精度を維持しながら、計算コストを一層低減することができる。
【0011】
また、本発明に係るプラズマシミュレーションプログラムは、上記プラズマシミュレーション方法をコンピュータに実行させる。このプラズマシミュレーションプログラムにおいても、上記作用効果、すなわち、計算コストを低減することが可能となるという作用効果が奏される。
【発明の効果】
【0012】
本発明によれば、計算コストを低減することが可能となる。
【発明を実施するための形態】
【0014】
以下、添付図面を参照して、本発明の好適な実施形態について詳細に説明する。なお、以下の説明において同一又は相当要素には同一符号を付し、重複する説明を省略する。
【0015】
図1は、一実施形態に係るプラズマシミュレーション方法の対象となるRPD装置を示す概略図である。
図1に示すように、本実施形態に係るプラズマシミュレーション方法及びプラズマシミュレーションプログラムは、一例として、RPD(Reactive Plasma Deposition)装置1をシミュレートの対象としている。そこで、まず、このRPD装置1について説明する。
【0016】
RPD装置1は、例えば、材料蒸発源Maを低基板温度条件下(〜200℃)で成膜用基板2に成膜するイオンプレーティング装置である。材料蒸発源Maとしては、ITOやZnO等の透明導電材料や、SiON等の絶縁封止材料が例示され、成膜用基板2としては、例えばガラス基板やプラスチック基板等の板状部材が例示される。
【0017】
このRPD装置1では、真空チャンバ3において、プラズマ源4から発生したプラズマPが、ステアリングコイル5とハースコイル6a及びハース磁石6bから構成されるプラズマビームコントローラ6とで形成される磁場によって、材料蒸発源Maへ導かれる。プラズマP中を流れる電流が例えば150Aと大きいことから、材料蒸発源Maが熱蒸発され、一部がプラズマPと相互作用を起こしてイオン化・励起される。そして、イオンはプラズマP中の電位勾配によって加速され、成膜粒子Mbとして広がり、成膜用基板2にて成膜が行われる。
【0018】
特に、ここでのRPD装置1は、プラズマP中を流れる電流が大きいことから、当該電流により生じる自己誘導磁場の影響が強いものとされており、この自己誘導磁場によってプラズマPが中心軸からズレ、ねじれが発生する場合があるものとされている。
【0019】
図2は、一実施形態に係るプラズマシミュレーションプログラムの構成を示す概略図である。
図2に示すように、プラズマシミュレーションプログラム10は、プラズマシミュレーション方法をコンピュータに実行させるためのプログラムであって、コンピュータが備える記録媒体50に形成されたプログラム格納領域50a内に格納されている。
【0020】
このプラズマシミュレーションプログラム10は、初期条件の入力処理を行う初期条件入力モジュール11と、粒子モデルを用いてイオン及び中性粒子に関する計算処理を行う粒子モデル計算モジュール12と、流体モデルを用いて電子に関する計算処理を行う流体モデル計算モジュール13と、収束判定処理を行う収束判定モジュール14と、を備えている。
【0021】
粒子モデル計算モジュール12は、運動方程式の計算処理を行う運動方程式計算部12aと、粒子境界条件の入力処理を行う粒子境界条件入力部12bと、衝突計算処理を行う衝突計算部12cと、を含んでいる。流体モデル計算モジュール13は、電子流体方程式の計算処理を行う電子流体方程式計算部13aと、自己誘導磁場の計算処理を行う自己誘導磁場計算部13bと、移動度の更新処理を行う移動度更新部13cと、を含んでいる。なお、上記の各処理については、後述のプラズマシミュレーション方法の説明において詳説する。
【0022】
図3は、一実施形態に係るプラズマシミュレーション方法を示すフローチャートである。本実施形態のプラズマシミュレーション方法は、粒子モデル及び流体モデルを用いて計算を行う、いわゆるハイブリッドPIC(Particle In Cell)法を利用するものである。
【0023】
すなわち、本実施形態のプラズマシミュレーション方法は、速度が大きく(〜10
6m/s)且つ粒子モデルでは計算コストのかかる電子について流体モデルで取り扱う一方で、速度分布を求めたいイオンや中性粒子については粒子モデルで取り扱う、電気的準中性を仮定している。
【0024】
なお、電子のエネルギ分布は、1温度のマックスウェル・ボルツマン分布に近似できることが見出されることからも、流体モデルを適用することは好ましい。また、本実施形態では、3次元的にプラズマシミュレーションを実施しており、シミュレートの対象となる所定領域(求めたい領域)として真空チャンバ3内を設定している。
【0025】
図3に示すように、本実施形態のプラズマシミュレーション方法では、まず、初期条件を入力する(S1)。初期条件としては、例えば、イオン及び中性粒子の初期位置及び初期速度、プラズマP中の初期電位、初期電子密度分布、初期磁場分布、及び、流体モデル計算における計算領域を少なくとも含んでいる。各初期値は、任意の値とすることもできるし、例えば初期磁場分布は、プラズマビームコントローラ6による磁場分布とすることもできる。
【0026】
次に、粒子モデルを用いて、イオン及び中性粒子に関する計算を行う粒子モデル計算パート(粒子モデル計算ステップ)を実施する。粒子モデル計算パートでは、運動後の粒子がセルを跨がず且つ1つの粒子が1時間刻みの間に2回衝突しないように、十分に時間刻みを短く設定し、粒子の運動と衝突現象とを分離して計算させている。なお、粒子の所属するセルを簡単に計算するために、xyz軸が直交するカルテシアン座標系を用いて数値計算を行っている(以下の数値計算において同じ)。
【0027】
具体的には、下式(1)を用いて、イオン及び中性粒子の速度及び位置を計算する(S2)。下式(1)の時間積分は、時間刻みを大きくとった際にも軌道計算の誤差が小さくなるように、高次精度の古典的4次のルンゲ・クッタ法を用いている。外力としては、電場及び磁場のローレンツ力を考慮し、重力は無視している。
【数1】
m
α:粒子αの質量、V
α:粒子αの速度、q:電化、φ:電位、
B(=B
self+B
ex):磁束密度、B
self:コイル等による磁束密度、
B
ex:プラズマP中の電流による磁束密度。
【0028】
ここでは、粒子の所属するセルの8隅にある格子点に電場・磁束密度の値を持たせ、粒子の位置へ線形補間を行っている。また、補間方法に合わせて粒子が数値計算のための力を受けないようにするために、1次の形状関数を与えて数密度分布や流速分布等の統計処理を行っている。
【0029】
続いて、以下の粒子境界条件を与える(S3)。流入粒子数は、実際の蒸発レートからの推測値を用いている。なお、粒子の種類や衝突の種類は、支配的な現象及び再現したい現象を絞って選定することができる。
アルゴン,酸素分子 : 鏡面反射
その他の粒子 : 消失
【0030】
続いて、粒子の衝突計算を実施する(S4)。ここでの衝突計算では、以下に示すように、少なくとも、電子−中性粒子衝突モデルによる計算と、重粒子同士衝突モデルによる計算と、を行っている。
【0031】
[電子−中性粒子衝突モデル]
電気的準中性の仮定からイオン密度と電子密度は等しく、また、後述の電子流体方程式の計算結果より、あるセルにおける電子温度が求まり反応確率がわかることから、反応確率に応じてイオン化等の中性粒子の反応処理を行う。
【0032】
[重粒子同士衝突モデル]
重粒子同士の衝突を計算し、運動量輸送を考慮する。本実施形態では、演算の簡略化のために剛体球衝突を仮定する。このように、重粒子同士の衝突を考慮することで、イオンが材料蒸発源Ma付近の背景アルゴンガスを排除する現象を再現できる。
【0033】
次に、流体モデルを用いて、電子に関する計算を行う流体モデル計算パート(流体モデル計算ステップ)を実施する。この流体モデル計算パートでは、まず、電子流体方程式を用いて、電子温度及び電位を求める(S5:第1ステップ)。具体的には、以下に示すように、運動量保存則、質量保存則及びエネルギ保存則から、電子温度及び電位を算出する。
【0034】
[運動量保存則]
RPD装置1中のプラズマPの振動周期(〜10
−12sec)は、イオンや中性粒子の運動を計算する時間刻み(〜10
−6sec)では、電子の運動が定常状態に達していると仮定できる。定常状態を仮定すると、運動量保存則から電子流束Γ
eを下式(2)で表すことができる。磁束密度の方向に対して複数の拡散・移流現象が考えられるが、下式(2)では、電子の流れとして、磁力線に対して平行方向及び垂直方向ともに古典拡散(μ
||)を仮定した。但し、磁力線を横切る方向(μ
⊥)に関しては、衝突周波数とラーマー周波数との比で決まる補正係数で補正している。なお、[μ
e]は、後述の移動度テンソルである。
【数2】
n
e:電子密度、T
e:電子温度。
【0035】
[質量保存則]
上記の運動方程式から導かれる電子流束の保存則を下式(3)式に示す。ここでは、電位φに対して線形な方程式として、有限体積法を用いて離散化して解いている。
【数3】
γ
ion:単位時間・体積当たりのイオン化による電子の発生数密度。
【0036】
境界条件としては、下式(4)に示すように、熱流束に壁とプラズマ電位とのボルツマンファクターを乗ずることでシースの影響を考慮する。
【数4】
α:壁面上の法線ベクトルと磁束密度の成す角、
△φ:プラズマ端と壁面との電位差。
【0037】
[エネルギ保存則]
エネルギ保存則の式は、下式(5)で表される。下式(5)では、電子温度T
eに対して線形に近似して解いており、また、一項目の対流項と二項目の拡散項との比率に応じて、風上差分及び中心差分を使い分けるハイブリッド法により離散化を評価している。なお、下式(6)に示すように、エネルギ流束の境界条件としては、流束同様に熱流束を仮定している。
【数5】
【数6】
【0038】
続いて、自己誘導磁場を算出する(S6:第2ステップ)。自己誘導磁場の計算では、イオン及び電子の運動が定常状態に達しているとの仮定から、磁場も定常状態に達していると考えられるため、プラズマP中の電流が作るベクトルポテンシャル(下式(7))の回転として自己誘導磁場B
selfを計算する(下式(8))。なお、境界条件は、A
x=A
y=A
z=0としている。
【数7】
【数8】
μ
0:真空中の誘電率。
【0039】
このとき、電子流体計算全体を通して、計算を安定化させて行列の優対角性を保証するために、対角項と同符号の項を湧出し項へ移項して1反復毎に係数を修正すると共に、不足緩和係数をかけて物理量を更新する。また、各物理量に対して非現実的な解を除去するために、各物理量に最大値及び最小値を設けて制限する。なお、計算を安定化させるために、上述の安定化手法に限らず、他の安定化手法を用いてもよい。
【0040】
続いて、移動度を更新する(S7:第3ステップ)。下式(9)に示すように、移動度テンソル[μ
e]は、カルテシアン座標系において、磁束密度ベクトルの座標系から座標変換Λで変換される。
【数9】
【0041】
次に、計算結果の場における値の変動が十分に(所定値よりも)小さく、収束しているか否かを判定する(S8)。上記S8にてNoの場合、上記S2に戻り、イオン及び中性粒子の位置及び速度が更新されると共に、電子温度及び電位が更新される。一方、上記S8にてYesの場合、シミュレートが終了し、結果が出力される。当該出力結果には、例えば、シミュレートの対象となる所定領域における、イオン及び中性粒子の位置分布、速度分布、電位分布、電子密度分布、及び磁場分布が少なくとも含まれている。
【0042】
ここで、流体モデル計算パート(上記S5〜上記S7)において、自己誘導磁場の計算時の境界条件として無限遠を再現するために、計算領域を大きくする場合がある。しかしこの場合、例えば膨大な規模の格子について電子流体方程式を取り扱うことになるために計算が冗長となり、計算時間かかり過ぎるという虞がある。
【0043】
そこで、
図4に示すように、本実施形態では、自己誘導磁場の演算(上記S6)における計算領域R1を、電子流体方程式の演算(上記S5)及び移動度更新(上記S7)における計算領域R2よりも大きく設定している。つまり、上記S5,S7の計算領域R2を大きくせずに、上記S6の計算領域R1については、無限遠の境界条件を再現するべく大きい領域として計算点を生成している。
【0044】
具体的には、計算領域R2は、シミュレートの対象となる所定領域であって、求めたい計算領域である真空チャンバ3内に対応している。一方で、計算領域R1については、無限遠を再現するために好ましいとして、計算領域R2に対し各座標軸において所定倍以上拡大するように設定されている。
【0045】
所定倍としては、例えば計算結果の誤差が許容範囲(10%等)に収まる値とされており、好ましいとして5倍以上とされている。図示するように、計算領域R1,R2を立方体としたとき、ここでの計算領域R1は、各座標軸の正方向及び負方向それぞれにおいて、計算領域R2よりも5L大きくなるように設定されている。
L:計算領域R2の一辺の長さ。
【0046】
また、計算領域R1は、計算領域R2を含むように設定されており、計算領域R2に対し重複する重複領域を中央部分に有していると共に、計算領域R2よりも外側は重複しない非重複領域となっている。計算領域R1において外側の非重複領域では、その計算点密度(単位体積当たりの計算点の数)が重複領域の計算点密度よりも小さくされており、重複領域の格子よりも粗い格子となるように設定されている。
【0047】
以上のように、本実施形態において、上記S5,S7では、シミュレートしたい領域を計算領域R2とする一方、自己誘導磁場を求める上記S6では、当該計算領域R2よりも大きい領域を計算領域R1として計算することが可能となる。よって、自己誘導磁場の計算領域R1を、境界条件として無限遠を再現するように大きくして計算精度を維持しつつ、当該計算領域R1よりも小さい計算領域R2にて電子温度及び電位の算出並びに移動度の更新を行い、その計算を簡潔化できる。従って、計算コストを低減することが可能となる。
【0048】
すなわち、上記S5〜上記S7のそれぞれは、その計算成分に下式(10)をx軸、y軸及びz軸に関して有している。そのため、計算点Nの数に計算コストが大きく依存するところ、上記S6ではメッシュが大きい一方で上記S5,S7ではメッシュが必要最低限で小さく維持されていることから、上記S6のみで計算点Nの数が大きくなることとなる。従って、計算コストを低減でき、例えば計算時間を1/100まで短縮化することが可能となる。
【数10】
S:行列、
a,b:ベクトル。
【0049】
また、本実施形態の計算領域R1は、上述したように、計算領域R2を各座標軸において等距離拡大するように設定されている。これにより、例えば、計算領域R2を偏重して拡大するように計算領域R1が設定される場合に比べ、シミュレーションの精度を高めることが可能となる。
【0050】
また、本実施形態では、上述したように、計算領域R1において、計算領域R2と重複しない領域の計算点の密度が、重複する領域の計算点の密度よりも小さくされている。これにより、シミュレーションの精度を維持しながら、計算コストを一層低減することができる。
【0051】
以上、本発明の好適な実施形態について説明したが、本発明は上記実施形態に限られるものではなく、各請求項に記載した要旨を変更しない範囲で変形し、又は他のものに適用したものであってもよい。
【0052】
例えば、本発明は、RPD装置1のシミュレーションに適用できるだけでなく、核融合のシミュレーションやイオンエンジンのシミュレーション等にも適用することができる。特に、本発明は、プラズマP中を流れる電流が大きく当該電流により生じる自己誘導磁場の影響が強い場合のシミュレーションに好適に適用することができる。
【0053】
また、上記計算領域R1,R2の形状は限定されるものではなく、例えば球形状であってもよい。なお、本発明は、上記プラズマシミュレーション方法を実施するプラズマシミュレータとして捉えることもでき、また、上記プラズマシミュレーションプログラム10を記憶する記憶媒体として捉えることもできる。