第 6 章 DC 型セルの二次元放電シミュレーション 101
6.2 二次元シミュレーションの解析方法
y x
Anode(A)
Cathode(K) Wall
Resistor( ) Cell center
y
/2
図 6.6 二次元シミュレーションの放電セル構造（表 6.1のCase 1）
のために電極間距離dを50の格子，セル幅の半分Wyを20の格子に等間隔で分割 した。
解析手法についても一次元解析と同様，位置に対応する瞬時電界に対応する定 常状態のボルツマン方程式により求めた電離係数や励起係数を用いる局所場近似 による二次元流体モデルを採用する。この中で発光特性に大きく関与する共鳴放 射の閉じ込めは，二次元以上で解くことは非常に困難なため，一次元解析で用い た手法 [7]をy軸方向にも適用した。この扱いは共鳴放射の閉じ込めを実効的な寿 命 [8]として近似して解析する例が多い中，電子衝突による共鳴準位原子の階段励 起などを正確に取り扱うことができるために高い精度が期待できる。
s粒子(s=1∼12)の基本的な粒子密度連続の式は次式で与えられる。
∂Ns
∂t + ∂Φsx
∂x + ∂Φsy
∂y =Gs(x, y, t)−Ls(x, y, t) (6.12) ここで，Nsはs粒子の密度，ΦsxとΦsyはそれら粒子束のx軸とy軸方向の成分，
GsとLsはs粒子と他の粒子との反応による生成項と消滅項である。
粒子束ΦsxとΦsyは次式で与えられる。
Φsx =Nsvsx−Ds∂Ns
∂x (6.13a)
Φsy =Nsvsy −Ds∂Ns
∂y (6.13b)
ここで，vsxとvsyはそれぞれドリフト速度のx軸とy軸方向の成分，Dsは拡散 係数である。
セル内の電位V はポアソンの式で計算される。
∂2V
∂x2 +∂2V
∂y2 =−ρ
ε (6.14)
Ex =−∂V
∂x (6.15a)
Ey =−∂V
∂y (6.15b)
ここで，ρはイオンと電子の密度差による空間電荷密度，εはガスの誘電率，Ex
とEyは電界のx軸とy軸方向の成分である。
6.2.3 境界条件と拡散寿命によるセル壁面損失の取り扱い
境界条件として陽極とセル壁面で，電子，イオンおよび励起粒子の密度を零と する。また陰極面の励起粒子の密度も零とするが，陰極直前のイオン束密度は一 定とした。
4種のイオンと2種の準安定原子および5種の光子が陰極に入射したときの陰極 での二次電子放出は，次式によって与える。
Φe=ΣγiΦi+ΣγmΦm+ΣγpΦp (6.16) ここで，Φe, Φi, ΦmおよびΦpは，二次電子束密度，イオン束密度，準安定原子束 密度および光子束密度，そしてγi, γmおよびγpは，それぞれ3種類の下付文字に 対応するイオン，準安定原子および光子による二次電子放出係数である。
電界の計算に必要な積分定数として，次式による時間(t)と位置(y)について，
セルの印加電圧を境界条件として与える。
Vg(t) =− d
0 Ex(x, y, t)dx=Va(t)−IdR (6.17) ここで，dは電極間距離，Idは放電電流，そしてRは電流制限のための直列抵抗 の抵抗値である。セルの隔壁の電位は，電極の境界条件とセルの中心線上の電位 分布を境界条件として与えたポアソンの式から見積もるが，このセルの中心線上 の電位分布は一次元解析により計算する。
電子とイオンによる放電電流Idは次式で与えられる。
Id(t) = 2Wzq d
Wy
0
d
0 (vexNe−ΣvixNi)dxdy (6.18) ここで，WyとWzは放電セル幅の半分の長さとセルの厚みであり，vexとvixは電 子とイオンのドリフト速度のx成分，NeとNiは電子とイオンの密度である。
x–y座標の二次元モデルでは，z軸方向については境界の設定もなく，粒子密度な ど全ての量を一様として扱うために損失は生じない。しかしながら，実際のパネル のセル壁での粒子の損失を近似して取り入れるために，励起粒子について式(6.12) のLsの中に，次式で求めた実効寿命 [9]から推定した損失項を加えている。
τ =λ2/D (6.19)
ここで，実効寿命τは拡散長λ (= WZ/π)と拡散係数Dにより計算される。この 拡散係数は電子とイオンの作用による両極性拡散を考える必要があるが，電子と 4種のイオンによる両極性拡散係数は複雑な計算になる。そこで電子と4種イオン の中で二桁以上も密度が大きいXe+を対としてその係数を見積もり使用した。
6.2.4 計算方法および数値解析の時間短縮と安定化手法
微分形式の粒子密度連続の式について，時間進展を陽的解法とした有限体積法 によって数値解析を行った。この計算での時間刻みΔtは，放電セルの一つの格子 Δxを最も速い電子が通過する時間よりも短くする必要があり，それはどの位置に おいてもクーラン条件Δt ≤ Min[Δx/ve(x, y)] を満足しなければならない。そこ で，一般的な条件Δt= 0.5·Min[Δx/ve(x, y)]で計算をしたが，Δtは2∼10 ps程度 の大きさであった。このように陽的解法ではΔtを充分小さくする必要があるが，
陰的解法のように逆マトリクスを解く必要がないので，ここで扱う数値解析では 総合的に計算時間を短くすることができる。多次元を扱うシミュレーションでは マトリクスの次元が増加するために計算時間が膨大になることが大きな相違点で ある。そこで，連立する非線形の方程式は，一つ前の時間刻みの値を使用するこ とによって，線形の方程式に変換した。
ここで，計算手順を以下に示す。
1 x軸方向のs粒子の密度の差分値を，一つ前の時間（t−Δt）で計算した生 成項と消滅項を使用して式 (6.12)の連続の式により時刻tにおける値を計算 する。
2 y軸方向の粒子密度の差分値も同様に計算する。
3 電界は前記の境界条件と式(6.14)によって計算する。
4 時間をΔtだけ増やしt+ Δtで同様の手順を繰り返す。
このようにして放電を維持させて繰り返して得られる計算結果は，計算の最初 に与える粒子密度分布の影響を受けることはないが，その分布は収斂の速さに影 響する。このため，初期密度は可能な限り計算結果に近い値にすることが望まれ る。例えば，同じ計算条件の下での一次元シミュレーションの結果から密度の初 期値を見積もるなどの方法も考えられる。
計算時間の短縮については一つの手法を提案してその効果を確認した。電子の ドリフト速度や拡散係数の値は，イオンや励起粒子に比べて二桁程度大きい。そ こで，動きの速い電子の密度は時間刻みΔt毎に計算するが，動きの遅い他の粒 子の密度は20倍の時間刻み20Δtで計算する。この手法を導入して計算した場合，
t=1.0μs時の最大放電電流値で0.002 %以内の誤差であった。また，数値解析の安 定性向上のために，放電電流値に時間刻みΔt毎の瞬時値を使用するのではなく，
表 6.1 DC型セルの二次元シミュレーションの計算条件
Case 1 Case 2 Case 3
Filling gas He(90 %)-Xe(10 %) ← ←
Gas pressure 27 kPa (200 Torr) ← ←
Cathode width 0.18 mm 0.4 mm ←
Anode width 0.4 mm 0.18 mm ←
Electrode distance (d) 0.35 mm ← 0.2 mm
length (WZ) 0.1 mm ← ←
Pulse Voltage (VS) 165 V 200 V ←
period 4μs ← ←
duration 1μs ← ←
Bias voltage (VB) 150 V ← ←
Resistance (R) 1 MΩ ← ←
数Δtの平均値を使用した。この放電電流の取り扱いについての妥当性も一次元シ ミュレーションを用いて確認している。
この二つの手法を取り入れて開発したシミュレーションコードにより，高い安定 性を維持しつつ計算時間を短縮することに成功した。一例であるが，この二次元 シミュレーションの計算時間は，中央演算装置（CPU: Central Processing Unit）
にクロック周波数200 MHzのUltra SPARK CPUを搭載したワークステーション を使用した場合，Case 1（以下に述べる）の条件で印加パルス1周期（4.0μs）の 計算時間が約3時間であった。この処理時間は，電子とイオンの時間刻みを同じ にして解析する通常の二次元シミュレーションの所要時間の約1/4である。