Transcript 第18回プラズマエレクトロニクス講習会資料 (ppt:4.8MB)

第18回プラズマエレクトロニクス講習会
プラズマプロセスモデリングと将来技術の最前線
実用シミュレーション技術
2007年11月2日
ペガサスソフトウェア株式会社
松永史彦、中舘博
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
1
シミュレーションの必要性
プロセス開発の現状
ブラックボックス
装置の運転条件
幾何形状、ガス種、
ガス圧、パワー等
装置内
基板
プラズマやラジカ
ルの状態
薄膜や微細構
造の状態
最終的な製品の状態を評価して、運転条件にフィードバッ
クする。装置内の状態を測定するのは困難で、ブラック
ボックスとして扱われる。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
2
シミュレーションの必要性
シミュレーションの有効利用
数値シミュレーションにより
装置内の状態を評価する。
装置の運転条件
幾何形状、ガス種、
ガス圧、パワー等
基板に到達するイオン、ラ
ジカル、スパッタ粒子の
•フラックス
•エネルギー分布
•入射角度分布 など
装置内
基板
プラズマやラジカ
ルの状態
薄膜や微細構
造の状態
装置内の状態を可視化することで、プロセ
ス開発の効率を向上させることを目指す。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
3
プラズマ装置内に生じている現象
方程式の非線形性
ポアソン方程式
荷電粒子の運動方程式
荷電粒子の生成
Ra = （電子エネルギー分布の関数）
•荷電粒子の運動は電位分布（電界）の影響を受ける
•荷電粒子が運動すれば空間電荷分布が変化し、電位分布も変化する
•荷電粒子の生成率は電子エネルギー分布の関数であり、電子エネル
ギー分布は電位分布の関数である
非線形性が非常に強い。短いタイムステップ
で時間を追跡していくしかない。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
4
プラズマ装置内に生じている現象
気相反応 （CF4の場合）
もっと大きい分子(例えば
C4F8)であれば、考慮すべき反
応式の数はさらに増大す
る！！
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
5
プラズマ装置内に生じている現象
表面反応
•エッチング
•デポジション
•スパッタリング
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
6
PIC/MCC法
Particle-In-Cell /
Monte Carlo Collision Method
プラズマの粒子シミュレーション法
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
7
PIC/MCC法
基本的な考え(1)
装置内に存在する例えば1015個程度の
電子・イオンを、106個程度の
“超粒子”(super particle, サンプル粒
子)
で置き換えてシミュレートする。
1015/106=109 超粒子の“重み”という。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
8
PIC/MCC法
基本的な考え
1. ポアソン方程式を解く
2. 荷電粒子の運動方程式
3. 境界条件
4. 荷電粒子の衝突
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
9
PIC/MCC法
実行のイメージ
イオン
電子
中性気体をバックグラウンドとして
扱い、電場、磁場中の電子・イオ
ンの超粒子の運動を追跡する。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
10
PIC/MCC法
クーロン相互作用
クーロン力は遠距離力である。
周りの全ての荷電粒子との相互作
用を計算するのは計算時間がかか
りすぎ非現実的。
空間をメッシュ分割し、超粒子の電
荷をグリッド上に集約して、ポアソン
方程式を解く。
q３
q４
q1
q2
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
11
PIC/MCC法
荷電粒子の運動方程式
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
12
PIC/MCC法
境界条件

電位の境界条件




ノイマン境界条件（法線方向の電位勾配が０）
ディリクレ境界条件（境界の電位を与える）
周期境界条件
荷電粒子の境界条件
 境界で消滅
 材料表面に表面電荷が蓄積
 2次電子放出
 周期境界条件
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
13
PIC/MCC法
荷電粒子の壁面での境界条件
+
壁面
イオンが壁面に入射すると、
•中性化して反射する
•固体の内部に進入する
イオンは壁面で消滅するとして取り扱う
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
14
PIC/MCC法
荷電粒子の壁面での境界条件
+
誘電体表面
- - -
-
誘電体表面に荷電粒子が入射
→表面に電荷が蓄積する
ポアソン方程式を解くときにこ
の表面電荷も考慮する。
イオンと電子の入射フラックス
が等しくなるように電位が変化
する。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
15
PIC/MCC法
荷電粒子の壁面での境界条件
+
壁面
-
イオンが壁面に入射すると、電子が
放出される場合がある。入射者イオ
ン1個あたり平均何個の電子が放出
されるか→2次電子放出係数
金属材料では 0.1前後の値をとるも
のが多い。MgO のように 3.0 とい
う大きい値をとる材料もある。
DCプラズマの場合には2次電子放出
がプラズマを維持するため、2次電子
放出係数は特に重要。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
16
DSMC法と
テスト粒子モンテカルロ法
希薄気体流れ場の
粒子シミュレーション方法
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
17
DSMC法




Direct Simulation Monte Carlo Method
直接シミュレーションモンテカルロ法、モンテカルロ直接法
1960年代に G.A. Bird が最初に始めた。
希薄気体の流れ解析に適する方法であるが、原理的には
（計算量の問題を無視すれば）連続流の領域も解析可能。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
18
DSMC法
希薄気体



ナビエストークス方程式が適用できない。
ボルツマン方程式に基づいた解析が必要。
Kn数> 0.01 という条件が目安。
 Knudsen number, クヌーセン数、クヌッセン数
 気体の流れ場の希薄度を表す無次元数。
 Kn = l/d l; 平均自由行程, d; 代表長さ
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
19
DSMC法
基本的な考え(1)
装置内に存在する例えば1023個程度の
気体分子を、106個程度の
“超粒子”(super particle, サンプル粒
子)
で置き換えてシミュレートする。
1023/106=1017 超粒子の“重み”という。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
20
DSMC法
基本的な考え(2)
超粒子の運動

並進移動

粒子同士の衝突

壁との衝突、流入、流出
を短い時間ステップで追跡する。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
21
DSMC法
基本的な考え(3)

計算領域を平均自由行程の大き
さのセルに分割する。

セル内にある粒子同士が、ある
確率で衝突を起こすものとして
取り扱う。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
22
DSMC法
実行イメージ(1)
超粒子が速度を持って
運動している。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
23
DSMC法
実行イメージ(2)
1タイムステップ並進
移動後の状態。
1タイムステップの並進
移動中は粒子は衝突し
ないとの仮定
“分離の原理”
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
24
DSMC法
実行イメージ(3)
セル内の粒子を2
個選び、衝突試行。
選ばれた2粒子の相対速
度の大きさから実衝突が
起こるかどうか判定される。
“最大衝突数法”
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
25
DSMC法
実行イメージ(4)
衝突が起こったと判定さ
れれば、2粒子の速度を
変える。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
26
DSMC法
サンプリング(1) 密度分布
各セルに超粒子が平均何個
存在しているかを統計処理し
て密度分布を求める。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
27
DSMC法
サンプリング(2) 速度分布
超粒子の速度を平均して、各
セルの流速を求める。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
28
DSMC法
衝突する粒子対の選択方法
いくつかの方法が提案されてきた。



南部法
変形南部法
最大衝突数法
最大衝突数法がもっとも計算効率が良
いとされている。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
29
DSMC法
最大衝突数法(1)


粒子 i と粒子 j が衝突する確率は、i と j の相対速度
の関数（剛体球の場合は比例）である。
セル内の全ての2粒子の組み合わせの相対速度を評
価するのは、計算コストがもったいない。例えば、セル
内に粒子が 20個あれば、20x19=380対の相対速度
を評価しなくてはならない。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
30
DSMC法
最大衝突数法(2)


多くても衝突数は最大これくらい、という最大衝突数を
見積もる。
その最大衝突数だけランダムに2対を選ぶ。これを仮
衝突
(衝突試行)と呼ぶ。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
31
DSMC法
最大衝突数法(3)




2対の相対速度から、衝突確率を求める。
乱数を振って、衝突するかどうか判定する。
衝突したと判定されれば、2粒子の速度を変更す
る。
これを実衝突とよぶ。
衝突の散乱角はまた別の乱数を振って決める。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
32
DSMC法
壁面での境界条件(1)



鏡面反射
拡散反射
Maxwell型反射
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
33
DSMC法
壁面での境界条件(2)
鏡面反射
壁面
•
高温、高真空で脱気した
平坦な金属面
•
壁の原子量が大きく、気
体の原子量が小さい場
合
光が鏡で反射するように、粒子
が壁面で反射する境界条件
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
34
DSMC法
壁面での境界条件(3)
拡散反射
壁面 温度T[K]
入射速度によらず、T で決
まる速度分布に従う確率
で反射する。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
35
DSMC法
壁面での境界条件(4)
Maxwell型反射条件
確率αで拡散反射を起こし、
確率(1-α) で鏡面反射する。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
36
DSMC法
流入/流出境界(1)
計算領域外
計算領域内
？
？
流入境界 どのようなフラックスと速度
分布で粒子を入れるか？
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
37
DSMC法
流入/流出境界(2)
計算領域外
1/4nv
ある温度、圧力の平
衡状態の気体が存在
する。
計算領域内
単位時間あたりに境
界を横切って入ってく
る粒子の速度分布
流入境界
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
38
DSMC法
流入/流出境界(3)
計算領域外
計算領域内
流出境界
計算領域外に出た粒
子は単に取り除く
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
39
DSMC法
流入/流出境界(4)
計算領域外
計算領域内
流出境界の計算領域外にも
ある温度、圧力の平衡状態の
気体が存在すると考えるのが
自然。
流出境界
流出境界からも粒子は
入ってくる。この効果を反
射率で与えたり、流入境界
と同様の条件で与えたりす
る。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
40
DSMC法
流入/流出境界(5)



粒子は流入境界と流出境界のいずれからも入ってくるし、
また、出て行く。
境界から出て行くフラックスは、計算結果として求まる量
である。境界から入ってくるフラックスは境界条件として
与える量である。
DSMC法では流入境界と流出境界に区別はない。粒子
の出入りの正味のバランスの計算結果として、流入境界
か流出境界かが決まる。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
41
DSMC法
分離の原理
1タイムステップ
時間
並進移動
衝突
並進移動
衝突
衝突
本当は1タイムステップの間のばらばらな
時刻に衝突が起こっているはず。衝突と並
進移動を分離して取り扱ってよいための条
件を分離の原理という。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
42
DSMC法
分子モデル


2分子間の距離と働く力の関係をモデル化したもの。
多くの場合、剛体球モデルで十分。
分子間ポテンシャル
分子間ポテンシャル
距離
距離
d
剛体球モデル
逆べき分子モデル
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
43
テスト粒子モンテカルロ法
密度の大きい成分をバック
グラウンドとして扱う。その
中で密度の小さい成分のみ
を超粒子として取り扱い、そ
の運動を追跡する。
密度の小さいラジカル種や、
スパッタ粒子の挙動に興味
がある場合などに使われる。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
44
連続体流体と希薄気体の
振舞の違いの例(1)
クエット流れ
速度Vで壁が移動
流体
速度分布は？
静止
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
45
連続体流体と希薄気体の
振舞の違いの例(2)
連続体流体の場合
流体の速度は、壁面
上で壁の速度と等しく
なる。
希薄気体の場合
速度の滑りが
生じる。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
46
連続体流体と希薄気体の
振舞の違いの例(3)
Kn数無限大の極限では..
vy<0の粒子のvx速度分布
y
x
vy>0の粒子のvx速度分布
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
47
2体衝突近似に基づく
モンテカルロシミュレーション
Monte Carlo simulations using the binary collision
approximation
(MC-BCA)
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
48
MC-BCA
実行イメージ（１）
固体表面に高エネル
ギーのイオンが入射
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
49
MC-BCA
実行イメージ（２）
入射粒子はターゲット原
子と衝突する。
周りの原子からの相互作
用は無視して2体衝突と
して取り扱う。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
50
MC-BCA
実行イメージ（３）
ターゲット原子が受け取っ
たエネルギーが結合エネ
ルギーより大きければ、は
じき出される。反跳原子
はじき出された原子がまた
さらに別の原子をはじき出
す。反跳カスケード
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
51
MC-BCA
実行イメージ（４）
表面結合エネルギーよ
り大きいエネルギーで
はじき出された原子は
固体表面を飛び出す→
スパッタリング
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
52
MC-BCA
実行イメージ（５）
衝突によって入射粒子はエ
ネルギーを失う。
全ての反跳原子と入射粒
子のエネルギーが十分小さ
くなるまで追跡を続ける。
入射粒子が停止した位置
(深さ)をサンプリング
→Depth Profile
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
53
BC-MCA
計算例 入射エネルギー依存性
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
54
BC-MCA
計算例 入射角依存性
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
55
BC-MCA
計算例 Depth Profile
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
56
BC-MCA
計算例 放出角度分布
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
57
マグネトロンスパッタ装置
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
58
マグネトロンスパッタ装置
解析モデル
基板
ターゲット
磁石
２次元円柱座標系
でモデル化
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
59
マグネトロンスパッタ装置
１．静磁場解析
磁束密度分布
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
磁力線
60
マグネトロンスパッタ装置
２．PIC/MCC法 プラズマ解析
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
61
マグネトロンスパッタ装置
２．PIC/MCC法 プラズマ解析
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
62
マグネトロンスパッタ装置
２．PIC/MCC法 プラズマ解析
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
63
マグネトロンスパッタ装置
２．PIC/MCC法 プラズマ解析
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
64
マグネトロンスパッタ装置
２．PIC/MCC法 プラズマ解析
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
65
マグネトロンスパッタ装置
２．PIC/MCC法 プラズマ解析
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
66
マグネトロンスパッタ装置
３．MC-BCA スパッタリング解析
スパッタリングフラックス分布 スパッタリング粒子放出角度分布
（ターゲットの損耗分布）
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
67
マグネトロンスパッタ装置
４．スパッタ粒子の輸送
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
68
マグネトロンスパッタ装置
４．スパッタ粒子の輸送
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
69
マグネトロンスパッタ装置
４．スパッタ粒子の輸送
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
70
流体モデル
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
71
流体モデル
前提
平均自由行程 λ ∝ 1/p （p : ガス圧）
装置の代表長さ L
P 大
λ 小
粒子モデルでは計算がきつくなる。
λ/L ≪ 1 なら電子・イオンの集団をそれぞれ連続体と
して取り扱うことができる。
例 アルゴンプラズマの場合
電子温度 1 [eV] として λe～8.4×10-5/pTorr
λe/L ≪1
pTorr≫8.4×10-5/L
Ar+ 温度 室温 として λi～4.8×10-5/pTorr
λi/L ≪1
pTorr≫4.8×10-5/L
p > 10 ～100 mTorr
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
72
流体モデル
マクロ量
数密度： nα(r,t)=∫fα (r,cα,t)d3cα
流速：
vα(r,t)=∫ cα fα (r,cα,t)d3cα/ nα(r,t)=< cα>
エネルギー ：<εα>(r,t)=< mαc2α/2>=
∫ (mαc2α/2) fα (r,cα,t)d3cα/ nα(r,t)
r：位置座標, t：時刻, c：粒子の速度ベクトル
m：粒子の質量
f(r,c,t)：速度分布関数
Boltzmann 方程式
添字α：粒子の種類(e,ion)
これらの量をもちいて運動を記述する。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
73
流体モデル
基礎方程式のセット
連続の式
∂ne
∂t
+ ▽・Γe = Re
運動方程式
Γe = - neμe E - De▽ne
∂nion
∂t
+ ▽・Γion = Rion
Γion = ± nionμion E - Dion▽nion
電子エネルギーバランス式
5
5
∂
neDe▽<εe> +eΓe・E = Reng
Γ
<ε
>
e
e
(ne<εe>) +▽・ 3
3
∂t
電界
- ▽・( εrε0 ▽φ) = e ( ΣionnionZion – ne )
, E = - ▽φ
Γα:粒子フラックス (≡nαvα)
μα:移動度
Rα:粒子の生成レート
Dα:拡散係数
Reng:電子エネルギーの損失レート
E:電界
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
74
流体モデル
計算手順
これらの基礎方程式を離散化（差分法）して解く
E
電位に関するポアソン方程式
コイル電流から
生ずる誘導電界
・パワー
これが一番時間がかかる
Γe、Γion を計算
連続の式
ne, nion
電子エネルギーバランス式
t=t+Δt
<εe>
この間、中性粒子種の
密度は変化しないと考
える
モンテカルロ法
による電子エネ
ルギー分布関数
の計算
ハイブリッドモ
デルの場合
中性粒子種に
関する流れの
計算
数周期毎に１回やり取りする
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
75
流体モデル
基礎方程式の導出
基礎方程式（粒子種を表す添字は省略）
∫
1
c
eZ
∂f
+ c・▽rf +
E・ ▽cf
m
∂t
c2/2
=
1
c
∫
c2/2
∂f
∂t
( )
col
d 3c
d3c
E ：電界
∂f
∂t
( )
col
：衝突項
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
76
流体モデル
基礎方程式
連続の式
∂n
+ ▽・nv = R
∂t
運動方程式
∂mnv
∂t
+ ▽・mnvv + ▽・P - neZE= R
mom
エネルギーバランスの式
∂n<ε>
∂t
+ ▽・nv<ε> + ▽・P ・v + ▽・Q + env・E= R
eng
R:電離等の反応による粒子の生成レート
Rmom:衝突による運動量の交換レート
Reng:衝突によるエネルギーの交換レート
P:応力テンソル=mn<(c-v)(c-v)>=mn<c’c’>
Q:熱流ベクトル=mn<c’2c’>/2
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
77
流体モデル
生成レート（電子の場合）
R= (νi – νa) n
Rmom= - mnνm v
Reng= - nνm(2m/M) <ε>
- nνiEi –nνexEex - ……
νk/N = ∫σk(c)f(c)c d3c
νi
νa
νm
νex
Ei
Eex
N
M
σk
：電離の衝突周波数
：電子付着の衝突周波数
：運動量移行衝突周波数
：励起の衝突周波数
：電離エネルギー
：励起エネルギー
：衝突相手分子の数密度
：衝突相手分子の質量
：k衝突の衝突断面積
（分子は電子に対して静止しているという近似）
・各種衝突について、衝突断面積のデータが必要
・電子の速度分布関数を知る必要
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
78
流体モデル
運動方程式の簡略化（ドリフトー拡散モデル）
・応力テンソルPは等方的  スカラー圧力 p
p=mn<c’2>/3 = 2nεT/3 = nkBT
1
mnνm
T:電子温度
∂mnv
+ ▽・mnvv + ▽nkBT + neE = - v
∂t
左辺第１項のオーダー : v/νmτv=(τm/τv)v
左辺第２項のオーダー: v2/νmL = (τm/L)v2 = (λm/L)(v/vT) v
τm :運動量緩和に関する平均衝突時間（1/νm）
τv : マクロな速度が変化する時間尺度
λm : 運動量緩和に関する平均自由行程
L :マクロな速度が変化する空間尺度（例えば電極間距離）
vT: 電子の熱（ランダム）運動の速度スケール √<c’2>
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
79
流体モデル
運動方程式の簡略化（ドリフトー拡散モデル）
・τm/τv ≪ 1 の場合は第１項（時間微分）は無視できる。
アルゴン中の電子の場合 νm ～ 1010ptorr [/s]
τv=107 [s]（10MHzに相当）とすると
τm/τv～10-3/ptorr
p > 10 [mTorr] なら条件が成り立つ。
・v/vT ≪ 1 の場合は第２項（対流項）は無視できる。
プロセスプラズマではだいたいこの条件は
満たされていると考えてよい。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
80
流体モデル
運動方程式の簡略化（ドリフトー拡散モデル）
圧力勾配項
1
nevT2
1
▽nekBTe
meneνm = ne ▽( 3νm ) = ne ▽(
1
De
= ne ▽(neDe)～
ne ▽ne
eE
meνm
外力(電界）項
neτmvT2
)
3
= μeE
μe =
e
meνm
拡散係数
De =
τmvT2
3
フラックス
Γe = - neμe E - De▽ne
移動度
=
kBTe
meνm
=
μekBTe
e
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
81
流体モデル
エネルギーバランス（電子）
∂n<ε>
∂t
+ ▽・nv<ε> + ▽・P ・v + ▽・Q + env・E= R
eng
・熱流ベクトルは温度勾配に比例
Q= - K▽Te = -
2K
▽εT = 3kB
5
neDe▽εT
3
(K=5De/2)
・圧力はスカラー
▽・P・ve = ▽・(pve) = 2 ▽・(neεTve) =
3
2
▽・(ΓeεT)
3
・平均速度の運動エネルギーは、熱運動のエネルギーに対して無視
<ε> ≒ εT ≫ meve2/2
・エネルギー損失
Reng = -ne(νm
2me ν
+ i
M
Ei
Eex
)<ε> = -neνε<ε>
+νex
<ε>
<ε>
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
82
流体モデル
電子に関する基礎方程式のセット
連続の式
∂ne
∂t
運動方程式
Γe = - neμe E - De▽ne
+ ▽・Γe = (νi -νa) ne
エネルギーバランス式
∂
∂t
5
(ne<ε>) +▽・ 3 Γe<ε>
-
5
neDe▽<ε>
3
+eΓe・E
= -neνε<ε>
・νε≪ νm であるのでエネルギーバランス式においては
一般に時間微分項を無視することはできない。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
83
流体モデル
イオンに関する基礎方程式のセット
・運動方程式に関しては、電子と同様ドリフトー拡散モデル
・イオンのエネルギーは電離や電子付着にほとんど影響しない
ので、エネルギーバランス式は解かない。イオン温度は一定
と仮定。
連続の式
∂nion
+ ▽・Γion =
∂t
νine ( for positive ion)
νane ( for negative ion)
運動方程式
Γion = ± nionμion E - Dion▽nion
μion =
(mion+M)e|Z|
mionMνm,ion
Dion =
μionkBTion
e
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
84
流体モデル
電界についての基礎方程式
・電子／イオンの運動を計算するために電界Eを求める必要
・電気的に中性とはいえ、シースなどでは電気的中性
が破られる。( nion > ne )
・電子／イオン密度の空間分布の変化により電界も変化
電位に関するポアソン方程式を解く
- ▽・( εrε0 ▽φ) = e ΣαnαZα
E = - ▽φ
φ: 電位
ε0 : 真空の誘電率
εr : 比誘電率 （プラズマ領域では 1）
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
85
流体モデル
電子の速度分布関数に対する仮定
衝突周波数： νk/N = ∫σk(c)f(c)c d3c
（１） 局所マクスウェル分布を仮定
f(c)=fLM(c;Te)= 4π
me
2πkBTe
3
2
2
m
c
e
c2 exp 2kBTe
νk/N はTeの関数
色々なTeの値に対してνk /N を計算し、Te～ νk /Nの関係を
テーブルとして用意しておく。
連続の式
運動方程式（粒子フラックス）
電子エネルギーバランス式
電位のポアソン方程式
Te（セル毎）
νk
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
μe, De
R, Rmom,Reng
86
流体モデル
電子の速度分布関数に対する仮定
（２） 一様電界中の電子スウォームの分布関数を用いる
一様な電界Eの中を移動する電子群に対する２項近似の
ボルツマン方程式を解き分布関数を求める。
Nμe、NDe、νk/N および<ε>はそれぞれ換算電界E/Nの関数
色々なE/Nの値に対して計算を行い、 Nμe、NDe、νk/N ～E/N
またはNμe、NDe、νk/N ～ <ε>関係のテーブルを作成しておく。
連続の式
運動方程式（粒子フラックス）
電位のポアソン方程式
E/N（セル毎）
μe, De
E/Nの関数
R, Rmom,Reng (局所電界近似)
または
電子エネルギーバランス式
<ε>（セル毎）
μe, De
<ε>の関数
R, Rmom,Reng
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
87
流体モデル
ハイブリッドモデル（粒子モデル＋連続体モデル）
（３） 電子のモンテカルロシミュレーションを用いる
プラズマの連続体モデルの計算でえられた電子密度、分子
密度、および(一般には時間依存の）電界E(r,t)をもちいて、
同じ計算体系（形状）で電子のモンテカルロシミュレーション
をおこなう。
電子の速度分布関数f(r,c)をえる。（セル毎）
（ 実際にはエネルギー分布関数f(r,ε) ）
νk および Te がセル毎に求まる。
連続の式
運動方程式（粒子フラックス）
μe, De
R, Rmom,Reng
ne(r)（周期平均）
電界E(r,t)
電位のポアソン方程式
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
88
流体モデル
イオンの輸送パラメータ
電子と異なる点
・質量 mion ～ M
・分子速度が無視できない
νm/N = ∫∫σm(cr)f(cr)cr fM(cn)d3c d3cn
cr=c-cn , cnは分子の速度
・分子は多様なので、衝突断面積のデータがみつから
ない場合も多い
ドリフト速度 vd は換算電界E/N に依存する。
vd ∝ E/N (E/Nの小さい領域 ≪ 5.6 Td )
(1Td=10-17 volt・cm2 )
∝√(E/Nの大きい領域) (E/N  ∞)
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
理論が難しい
89
流体モデル
イオンの輸送パラメータ
（１） 実測データを使用する
・vd, Nμ,ND ～E/N の関係をテーブル表にしておく。
・文献等ではμそのものではなく、μ0 (reduced
mobility )の値をE/N またはE/p0に対してプロット
していることが多い。
μ0=(p0/760)μ, E/N[Td] =2.828×E/p0[volt/cm/Torr]
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
90
流体モデル
イオンの輸送パラメータ
（２） 次の近似式からμを求める。
3e|Z| π(mion+M)
Nμ=
8
2mionMkBTeff
1/2
1
Ω(1,1)(Teff)
kBTeff=kBTn +(1/3)Mvd2
Ω(1,1)(Teff)=
1 (k T )-3∫ε∞2σ (ε)exp(-ε/k T )dε
B eff
m
B eff
0
2
D=
μkBTeff
e
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
91
流体モデル
電磁場とのカップリング
ICPの場合
▽2Eθ – iωμ0σEθ= 0 を解く
Eθ：誘導電界
ω：コイル電流の角振動数
μ0:真空の透磁率
σ：プラズマの導電率
連続の式
運動方程式
=
e2ne
me(νm+iω)
ne(r)（周期平均）
νm
電位のポアソン方程式
Eθ
σ|Eθ|2(Power)
電子モンテカルロ計算
（ハイブリッドの場合）
電子エネルギーバランス式
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
92
流体モデル
中性ガス流れとのカップリング
衝突による励起、解離
ガスの流動
原子、励起種生成
原料ガスの消費
中性粒子種の密度変化
電離、電子付着などに影響
・中性ガスの流れシミュレーション
各種成分の質量保存式（対流－拡散形）
混合ガスの運動方程式
vgas
混合ガスのエネルギー保存式
Tgas
nα
（α:中性粒子種）
・電子・イオン運動の時間スケールとガス流れの時間
スケールはオーダーが異なるので、分離しておこなう。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
93
流体モデル
境界条件
・境界面における粒子フラックスを規定する。
Γdrift =
(ξnμE・n)^n ^（ξE・n >^0 のとき）
0
（ξE・n^> 0 のとき）
Γdiff = (-D ▽n・n^) n^
（n境界=0）
Γout = aΓdrift + bΓdiff
(a,b=0 または1)
ξ=±1
あるいは
Γout = (1/4)nvT n^, vT=√(8kBT/πm)
n: 境界面における外向き単位法ベクトル
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
94
流体モデル
境界条件
通常、電子に関しては
Γe,out = (1/4)neveT n^
Γe,in = - Σαγse,αΓα,out
γseは正イオンの衝突による２次電子放出係数。
総和は正イオン種αに関して取る。
Γe,境界= Γe,out + Γe,in
イオンに関しては
Γion,境界= Γion,out = Γion,drift + Γion,diff
電子エネルギーフラックスに関しては
qe,境界=2kBTeΓe,out + 2kBTe,in Γe,in （Te,inは仮定）
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
95
流体モデル
離散化
Γz i,j+1/2
Ez i,j+1/2
・２次元直交メッシュ
軸対称
平面２次元
・変数配置
staggered
Γr i-1/2,j
Er i-1/2,j
・コントロールボリューム法
nij
Teij
Dij
φij
セル(i,j)
μij
Γr i+1/2,j
Er i+1/2,j
Γz i,j-1/2
Ez i,j-1/2
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
96
流体モデル
Sharfetter – Gummel スキーム(指数法)
ni
Γi+1/2
vi+1/2
ni+1
Δx
Γi+1/2 = vi+1/2
exp(αΔx)ni - ni+1
exp(αΔx) - 1
α= (v/D)i+1/2 = (±μED)i+1/2
・Γ=nv-D(dn/dx) として、Γi+1/2を評価するのに、１次元定常問題
dΓ/dx=0 の解析解を局所的に使用。
・|α|Δx ≫1 では１次の風上法、 |α|Δx ≪1では中心差分になる。
計算が安定
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
97
流体モデル
ポアソン方程式のsemi-implicit 解法
∂ρ
= - ▽・Σe ZαΓα
∂t
= - ▽・Σe |Zα|nαμαE
～ - e(Σ|Zα| nαμα) ▽・ E
= - e(Σ|Zα | nαμα)
ρ
ε0
アルゴン
p=100 [mTorr]
ne=1012 [/cm3]
τM=1.5×10-12 [s]
ε0
誘電緩和時間: τM =
eΣ|Zα| nαμα
陽解法をもちいる場合 Δt < τM とする必要
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
98
流体モデル
ポアソン方程式のsemi-implicit 解法
- ▽・(ε0 ▽φ) = e ΣαnαZα
nα= nα,old + Δt (-▽・Γα+Rα)
= nα,old - Δt (|Zα|/Zα)μα nα ▽・E
+ Δt (Dα ▽2nα,old +Rα)
Δt
1+
τM
▽2φ= -
1
eΣαZα nα,old + Δt (-▽・Γα,diff + Rα)
ε0
Δt > τM としても安定に計算できる。
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
99
流体モデル
時間刻みΔt
Δt < min Δr/(μeEr), Δz/(μeEz)
・電界ドリフトによる移動が１ステップで１メッシュを
こえない。
・電子のドリフト速度で評価
・生成レートが大きい場合には、この評価よりもかなり
小さくしなければならない場合も有る。
・だいたいΔt =10-10～10-12 程度
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
100
流体モデル
その他
・(周期）定常の判定
目視による。
密度のピーク値や、モニター点における各物理量
の時間変化の出力が不可欠
・リスタート計算機能をつけるべき
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
101
ドライエッチング装置のシミュレーション
ハイブリッドモデル
ＤＳＭＣ
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
102
ドライエッチング装置
解析モデル
1 KW, 13.56 MHz
ガス圧 3 Pa
流量 200 sccm
SF6
150 V, 13.56 MHz
ハイブリッドモデル
+
DSMC法(ガス流れ）
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
103
ドライエッチング装置
反応式
e+F
e+F2
e+F2
e+F2
→
→
→
→
F(+)+2e
e+F2
F(-)+F
F+F+e
e+SF6 → e+SF6
e+SF6 → SF6 (-)
e+SF6 → SF5(+)+F+2e
e+SF6 → SF4(+)+2F+2e
e+SF6 → SF3(+)+3F+2e
e+SF6 → SF2(+)+4F+2e
e+SF6 → SF(+)+5F+2e
e+SF6 → SF5(-)+F
e+SF6 → SF3+3F+e
e+SF6 → SF2+4F+e
e+SF6 → SF5+F(-)
e+SF5 → SF4(+)+F+2e
e+SF5 → SF5(+)+2e
e+SF4 → SF4(+)+2e
e+SF3 → SF3(+)+2e
e+SF2 → SF2(+)+2e
e+SF → SF (+)+2e
X(+)+Y(-) →X+Y (X,Y∈SFn,F)
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
壁面反応
すべての壁面
F → 1/2F2
F2, SFn(n=1,5)は壁面
で0.1%損失
104
ドライエッチング装置
誘導電界・パワー
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
105
ドライエッチング装置
電子密度・生成レート
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
106
ドライエッチング装置
イオン密度
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
107
ドライエッチング装置
イオン密度・Ｆ原子密度
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
108
ドライエッチング装置
入射イオンのエネルギー分布関数
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
109
ＤＣホローカソード放電プラズマシミュレーション
ハイブリッドモデル
ＤＳＭＣ
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
110
ＤＣホローカソード
解析モデル
10
単位：mm
解析条件
2
・Ａｒ：133Pa
・カソード印加電圧およびγ+
1) -240[V], 0.035
2) -280[V], 0.04
30
陰極
・ハイブリッドモデル
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
111
ＤＣホローカソード
プラズマ密度
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
112
ＤＣホローカソード
電離レート・電子温度・電位分布
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
113
ＤＣホローカソード
入射イオンエネルギー分布
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
114
PDP cell DBD グロー放電シミュレーション
(1) RF電極
(2) DC電極
・ 流体モデル
・ 電子エネルギーバランス式
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
115
Plasma Display Panel (1)
誘電体
εr = 8
320μm
80μm
320μm
ガス : Ar
圧力 : 200 [Torr]
電圧 : Vpp=400 [V]
freq.=10 [MHz]
位相差180°
γ+= 0.04
40μm
600μm
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
116
電子密度・イオン密度
周期平均値
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
117
電離レート・電子温度
周期平均値
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
118
Plasma Display Panel (2)
ガス : He
圧力 : 500 [Torr]
50μm
400μm
400μm
50μm
誘電体
εr=3.8
600[V]
γ+= 0.35
0[V]
250μm
誘電体
εr=10
80μm
1050μm
300[V]
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
119
電子密度・イオン密度
10[ns]
20[ns]
30[ns]
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
1000[ns]
120
電位・電子温度
10[ns]
20[ns]
30[ns]
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
1000[ns]
121
電離レート
10[ns]
20[ns]
30[ns]
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
1000[ns]
122
プラズマ密度・ポテンシャルの時間変化(0～50 ns)
0～50[ns]
電子密度
He(+)密度
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
123
プラズマ密度・ポテンシャルの時間変化(0～1000 ns)
0～1000[ns]
電子密度
He(+)密度
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
124
おわりに
・熱、構造、流体、電磁気分野ではシミュレーションが定着
・プラズマ分野ではシミュレーションの普及はまだこれから
・計算時間の短縮
・断面積等のデータ整備
・実験とシミュレーション結果の照らし合わせ
・３次元化
・使い勝手
・メーカー・大学等研究機関・フトウェア開発者
の連携
装置設計により役立つシミュレーションツールを目指す
Copyright 2002-2007 PEGASUS Software Inc., All rights reserved.
125