Particle-PLUS コラム - プラズマモデリング(連載) - 第2回

ウェーブフロントCAEソリューションParticle-PULSコラムプラズマモデリング(連載) > 第2回

ウェーブフロントでは「粒子法」でプラズマをモデル化したシミュレーションソフトウェア『Particle-PLUS』や 「流体法」のプラズマシミュレーションソフトウェア『VizGlow』を取り扱っており, これらを適宜使い分けながら皆様方へ長年ソリューションを提供しております. 本コラムを読んでプラズマシミュレーションに少しでも興味をお持ちになりましたら 資料請求無料セミナーなど承りますので, いつでもお気軽にお問い合わせください.

第1回では, プラズマの運動論方程式としてBoltzmann方程式を導出しました. 第2回では,そのままでは一部の状況下でしか解くことが困難なBoltzmann方程式の速度モーメントを計算して, 現実的な計算コストでプラズマの運動を記述できる流体方程式を導出していきます. (長くなったので,コラム第3回と2回に分けて導出します. 途中の式変形に興味がなければ「1.4. 連続方程式(質量保存則)」は読み飛ばしても大丈夫です.)

1.3. 分布関数のモーメントMoments of distribution function

プラズマの流体方程式は,力が全く電磁的な場合のBoltzmann方程式 \begin{align} \frac{\partial f_{s}}{\partial t} + \vec{v}\cdot\left(\frac{\partial f_{s}}{\partial\vec{x}}\right) + \frac{1}{m_{s}}\left(\vec{E}+\vec{v}\times\vec{B}\right)\cdot\left(\frac{\partial f_{s}}{\partial\vec{v}}\right) = \left(\frac{\delta f_{s}}{\delta t}\right)_{C} \label{eq_boltzmann} \end{align} に関する速度モーメントを計算することで得られます.

まずは準備として,速度モーメントに関連するいくつかの巨視的物理量を下記のように定義しておきます. また,任意の物理量 $X$(スカラー・ベクトル・テンソルを問わない)に関する速度モーメント量を, \begin{align} \langle X\rangle = \frac{1}{n_{s}}\int X f_{s} \mathrm{d}\vec{v} \end{align} と表すことにします.

1.3.1. モーメント量の定義

[巨視的物理量] \begin{alignat}{2} \text{粒子数密度: } && n_{s} & = \int f_{s}\mathrm{d}\vec{v} \label{eq_density}\\ \text{平均速度(流速): } && \vec{u}_{s} & \equiv \langle\vec{v}\rangle \label{eq_mean_velocity}\\ \text{質量密度: } && \rho_{s} & \equiv m_{s}n_{s} \label{eq_mass_density}\\ \text{熱速度: } && \vec{w}_{s} & \equiv \vec{v} - \vec{u}_{s} \label{eq_thermal_velocity}\\ \text{粒子数フラックス: } && \vec{\varGamma}_{s} & \equiv n_{s}\vec{u}_{s} \label{eq_particle_flux}\\ \text{平均エネルギー: } && \varepsilon_{s} & \equiv \frac{1}{2}m_{s}\left(\vec{u}_{s}\cdot\vec{u}_{s}\right) \label{eq_mean_energy}\\ \text{応力: } && {\bf{P}}_{s} & \equiv m_{s}n_{s}\langle\vec{w}_{s}\otimes\vec{w}_{s}\rangle = m_{s}n_{s}\left(\langle\vec{v}\otimes\vec{v}\rangle - \vec{u}_{s}\otimes\vec{u}_{s}\right) \label{eq_stress}\\ \text{温度: } && T_{s} & \equiv \frac{m_{s}}{3k_{\mathrm{B}}}\langle\vec{w}_{s}\cdot\vec{w}_{s}\rangle \qquad \text{($k_{\mathrm{B}}$ はBoltzmann定数)} \label{eq_temperature}\\ \text{熱エネルギー密度: } && \theta_{s} & \equiv \frac{3}{2}n_{s}k_{\mathrm{B}} T_{s} = \frac{1}{2}m_{s}n_{s}\langle\vec{w}_{s}\cdot\vec{w}_{s}\rangle \label{eq_thermal_energy}\\ \text{熱流束: } && \vec{Q}_{s} & \equiv \frac{1}{2}m_{s}n_{s} \langle\left(\vec{w}_{s}\cdot\vec{w}_{s}\right)\vec{w}_{s}\rangle \label{eq_thermal_flux} \end{alignat} [散乱項] \begin{alignat}{2} \text{生成項(ソース項): } && S_{s} & \equiv \int \left(\frac{\delta f_{s}}{\delta t}\right)_{C} \mathrm{d}\vec{v} \label{eq_source_term}\\ \text{運動量交換項: } && \vec{K}_{s} & \equiv \int m_{s}\vec{v} \left(\frac{\delta f_{s}}{\delta t}\right)_{C} \mathrm{d}\vec{v} \label{eq_friction_term}\\ \text{エネルギー交換項: } && H_{s} & \equiv \int \frac{1}{2}m_{s} \left(\vec{v}\cdot\vec{v}\right)\left(\frac{\delta f_{s}}{\delta t}\right)_{C} \mathrm{d}\vec{v} \label{eq_energy_trans} \end{alignat}

1.3.2. 記法について(ベクトルの二項積・ベクトルとテンソル)

応力テンソルの定義に含まれるようなベクトル$\vec{A}$と$\vec{B}$の直積$\vec{A}\otimes\vec{B}$の表記に関して, 流体力学の分野内ではしばしば直積記号$\otimes$を省略して$\vec{A}\vec{B}$と書くことがあります. しかしながら,本コラムにおいては分野外の方にとっての数学的な分かりやすさを保てるように,その旨を明示しない限り直積記号を省略しないことにします. ちなみに,このような二項表示$\vec{A}\vec{B}$はダイヤディック(dyadic)表示と呼ばれ, この演算自身やその結果は二項積(dyadic)や二項テンソル(dyadic tensor)と呼ばれます. また,本コラム内においては,ベクトルを矢印付き $\vec{A}$ で,テンソルを太字立体 ${\bf{A}}$ で書き分けることにします.

1.4. 連続方程式(質量保存則)Continuity equation (mass conservation)

Boltzmann方程式の0次の速度モーメントからは連続方程式が導けます. \begin{align*} \int \left\{\frac{\partial f_{s}}{\partial t} + \vec{v}\cdot\left(\frac{\partial f_{s}}{\partial \vec{x}}\right) + \frac{q_{s}}{m_{s}}\left(\vec{E}+\vec{v}\times\vec{B}\right)\cdot\left(\frac{\partial f_{s}}{\partial\vec{v}}\right) \right\} \mathrm{d}\vec{v} = \int \left(\frac{\delta f_{s}}{\delta t}\right)_{C} \mathrm{d}\vec{v} \end{align*}

[左辺第1項] \begin{align*} \int \frac{\partial f_{s}}{\partial t} \mathrm{d}\vec{v} = \frac{\partial}{\partial t}\int f_{s} \mathrm{d}\vec{v} = \frac{\partial n_{s}}{\partial t} \end{align*}

[左辺第2項] \begin{alignat*}{2} \int \vec{v}\cdot\left(\frac{\partial f_{s}}{\partial\vec{x}}\right) \mathrm{d}\vec{v} & = \frac{\partial}{\partial\vec{x}} \cdot \int \vec{v}f_{s} \mathrm{d}\vec{v} && \text{ (演算子交換関係 $(\partial / \partial\vec{x})\cdot\vec{v}=\vec{v}\cdot (\partial /\partial\vec{x})$ より)} \\ & = \frac{\partial}{\partial\vec{x}} \cdot \left(n_{s}\vec{u}_{s}\right) && \\ & = \frac{\partial}{\partial\vec{x}} \cdot \vec{\varGamma}_{s} && \end{alignat*}

[左辺第3項] \begin{alignat*}{2} \int \vec{E}\cdot\left(\frac{\partial f_{s}}{\partial\vec{v}}\right) \mathrm{d}\vec{v} & = \int \frac{\partial}{\partial\vec{v}}\cdot\left(\vec{E}f_{s}\right) \mathrm{d}\vec{v} && \text{ (演算子交換関係 $(\partial /\partial\vec{v})\cdot\vec{E}=\vec{E}\cdot(\partial /\partial\vec{v})$ より)} \\ & = \int_{S} \left(\vec{E}f_{s}\right) \cdot \mathrm{d}\vec{S} && \text{ (Gaussの発散定理より)} \\ & = 0 && \text{ (下記参照)} \end{alignat*} 積分領域 $S$ は速度空間における $|\vec{v}|\to\infty$ での閉曲面で,$\vec{S}$はその面積ベクトルです. この積分は,$|\vec{v}|\to\infty$ 極限での $\vec{E}f_{s}$ の値を計算しています. 速度空間において全粒子は積分領域の内側にあるので, およそ $|\vec{E}|\propto|\vec{v}|^{-2}$ となるはずです*1. $f_{s}$ が有限の速度モーメント持つような分布関数であれば(言い換えれば,考えている粒子系が有限エネルギーであれば) $|\vec{v}|\to\infty$ において $f_{s}|\vec{v}|^{-2}\to0$ が成り立ちますので,上式の結果が得られます.

[左辺第4項] \begin{alignat*}{2} \int \left(\vec{v}\times\vec{B}\right)&\cdot\left(\frac{\partial f_{s}}{\partial\vec{v}}\right) \mathrm{d}\vec{v} && \\ & = \int \frac{\partial}{\partial\vec{v}}\cdot\left\{\left(\vec{v}\times\vec{B}\right)f_{s}\right\} \mathrm{d}\vec{v} - \int \left\{\frac{\partial}{\partial\vec{v}}\cdot\left(\vec{v}\times\vec{B}\right)\right\}f_{s} \mathrm{d}\vec{v} && \text{ (部分積分)} \\ & = \int \frac{\partial}{\partial\vec{v}}\cdot\left\{\left(\vec{v}\times\vec{B}\right)f_{s}\right\} \mathrm{d}\vec{v} && \text{ (下記参照)} \\ & = \int_{S} \left(\vec{v}\times\vec{B}\right)f_{s} \cdot \mathrm{d}\vec{S} && \text{ (Gaussの発散定理より)} \\ & = 0 && \text{ (下記参照)} \end{alignat*} $(\vec{v}\times\vec{B})\perp\vec{v}$ より,明らかに $(\partial /\partial\vec{v})\cdot(\vec{v}\times\vec{B})=0$ です. また,先ほどの左辺第3項の場合と同様に, $f_{s}$ が有限の速度モーメント持つような分布関数であれば $|\vec{v}|\to\infty$ において $(\vec{v}\times\vec{B})f_{s}\to0$ が成り立ちますので, 最後の等式が示せます.

[右辺] \begin{align*} \int \left(\frac{\delta f_{s}}{\delta t}\right)_{C} \mathrm{d}\vec{v} = S_{s} \end{align*} 散乱過程や輻射過程による電離や結合が起きなければ位相空間全域において粒子数は増減しないので,この項は $0$ になります. 逆に言えば,$S_{s}$ こそが散乱過程や輻射過程の結果による粒子の生成項(ソース項)を意味します.

[連続方程式] \begin{align} \frac{\partial n_{s}}{\partial t} + \vec{\nabla}\cdot\vec{\varGamma}_{s} = S_{s} \label{eq_continuity_eq} \end{align}


速度に関するモーメント積分を行ったことにより,連続方程式に現れる物理量が持つ独立変数は $(\vec{x},t)$ になりました. 即ち,$\vec{v}$ に依存しなくなりましたので,ベクトル微分演算子は空間微分のみが残ります. 以上のことから,今後はベクトル微分演算子の作用素が明らかに空間微分である場合には $(\partial/\partial\vec{x})\to\nabla$ と書くことにします.

導出過程より,密度分布関数 $f_{s}$ が一定であるなどの速度モーメントを計算できないような系に対しては $|\vec{v}|\to\infty$ における表面項が非零となり,連続方程式が成り立たないことが分かります. また,電磁気力以外の外力がある場合,その形によっても連続方程式が成り立たない(零にできていた左辺第3,4項が非零になるかもしれない)ことが分かります. しかしながら,有限エネルギー系におけるプラズマを考える限りにおいては,上記のような状況に出会うことはないでしょう.


*1 典型的な電磁気学の演習問題「一様に帯電した有限半径の球の周りの電場を求めよ」と似た状況であることを考えると明らかです.

第2回まとめ

今回は,いくつかの巨視的物理量の定義を与えました. これらの定義は今後も使いますので覚えておいてください. また,Boltzmann方程式の0次モーメントを計算して連続方程式を導きました. 次回は,1次モーメントから運動量方程式を,2次モーメントからエネルギー方程式を導きます. (ほとんど今回と同様の手順です.)

ウェーブフロントでは「粒子法」でプラズマをモデル化したシミュレーションソフトウェア『Particle-PLUS』や 「流体法」のプラズマシミュレーションソフトウェア『VizGlow』を取り扱っており, これらを適宜使い分けながら皆様方へ長年ソリューションを提供しております. 本コラムを読んでプラズマシミュレーションに少しでも興味をお持ちになりましたら 資料請求無料セミナーなど承りますので, いつでもお気軽にお問い合わせください.

著者プロフィール
  上野 崚一郎 | 博士(理学)
1991年 広島県生まれ
2019年 ウェーブフロント入社
2019年 広島大学大学院理学研究科 博士後期課程修了

学生時代は数値シミュレーションを使った素粒子論(格子ゲージ理論)の研究に従事. 入社後は,専門職(エンジニア)としてプラズマ解析ソフトウェアの開発をはじめとして, 解析コンサルティング業務や国内外のユーザー向けの技術サポート・トレーニングなどを担当.

最後までお読みいただきありがとうございます. お気づきの点や扱ってほしい話題がございましたらお気軽にお問い合わせください.