非一様空間格子幅・時間ステップ幅 および衝突効果 寺田 直樹 名古屋大学太陽地球環境研究所

１．Introduction 非一様な格子分布の利用は、数値流体力学の分野では ２０年以上にわたり盛んに研究が行われている。 複雑な形状の物体周りの場を解く （適切な境界条件） 解像が必要な領域に計算格子を集中 例１．飛翔体周りの電磁環境の解析 例２．衝撃波の反射解析 非一様な格子分布の利用は、数値流体力学の分野では ２０年以上にわたり盛んに研究が行われている。 ⇒　しかし、粒子・ハイブリッドモデルでは 　　 直交格子を用いるのが一般的

講義の流れ ・ 非一様な格子分布 （テキスト p.233～） ・ 粒子数の動的変更 （テキスト p.237） を粒子・ハイブリッドモデルで利用する為の 手法を中心に説明。 　・ 粒子数の動的変更　（テキスト p.237） 　・ 非一様な Δｔ　（テキスト p.238～） 　・ 衝突効果　（テキスト p.240～） 非一様な格子分布と併用されることがしばしばある

２．非一様な格子分布 （流体計算） 格子の種類 ・ 構造格子 (b) ・ 非構造格子 (a) ・ 構造格子と非構造格子の組み合わせ ・ 解適合格子

非構造格子 ・ 格子点の並びに規則性をもたない ・ 格子の自由度が高い ・ 有限要素法や有限体積法で使用される ・ 三角形（三次元では四面体）が最も多く用い られるが、四角形（三次元では六面体）や その他のセル形状やその混在もしばしば 用いられる 例．非構造格子を用いた デルタ翼周りの流れ場の計算

構造格子 ・ 格子点の並びに規則性をもつ ・ 差分法や有限体積法で使用される 境界適合格子 などがある マルチブロック 重なり格子 格子自由度の 低さを補う 例．境界適合格子を用いた 電極間プラズマの計算

非構造格子と構造格子の組み合わせ 解適合格子 三角形（三次元では四面体）セルを用いると 流速計算が多くなる ⇒　高レイノルズ数流れの計算では、高い 　　 空間分解能が必要となる物体近傍に 　　 のみ構造格子（またはプリズム形状の 　　 非構造格子）を用いる　→ 例 解適合格子 流れの計算結果をもとに格子点を移動、 または格子点を追加・削除

非構造格子 [Togashi et al., 2001]

[Togashi et al., 2001]

CFD (Computational Fluid Dynamics) 　　　　航空宇宙分野などが牽引役 　　　　・ 1970年代後半～80年代 　　　　　　　効率の良い計算が最重要課題 　　　　　　　⇒　境界適合格子 　　　　・ 1990年代～ 　　　　　　　実形状周り流れの計算を目指す 　　　　　　　⇒　非構造格子

CFD　　　　　　　　　　　　　実形状　　　　　　　　　　　非構造格子 　　　　　　　　　　　　　（形状の微妙な差で 　　　　　　　　　　　　　　空力性能が決まる） Space　plasma　　　　　　　　実形状　　　　　　　　　　構造（or 非構造） の大部分　　　　　　　　プラズマの粒子性　　　　　　　　 格子 　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　＋ 　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　　粒子コード

２.２ 粒子コードへの適用 ・ 非構造格子 ・ 構造格子（境界適合格子） ・ 解適合格子 ⇒ e.g., Lipatov [2002] ２.２　粒子コードへの適用 ・ 非構造格子 　　　　格子形状が三角形（四面体）　⇒　e.g., Kazeminezhad et al. [1995] 　　　　　　プログラミングがやや複雑 　　　　　　物体形状が複雑な場合に使用 （飛翔体環境など） 　　　　格子形状が四角形（六面体）　⇒　e.g., Kallio and Janjunen [2001] 　　　　　　三次元で六面体を用いた場合は、+n レベルの格子では格子体積が 1/8n に。 　　　　　　　　　　　→　隣接するグリッド間内の粒子数の変化が激しい。 　　　　　　　　　　　　　 粒子分割・合体（後述）が必要となる場合が多い。 ・ 構造格子（境界適合格子） 　　　　　　プログラミングが比較的容易で汎用性に富む 　　　　　　物体形状が相当複雑でない限り十分 ・ 解適合格子 　⇒　e.g., Lipatov [2002] 粒子分割・合体

２.３　境界適合格子 ＋ 粒子コード ・ 境界適合格子における変換式 ・ 粒子コードへの適用 　　- 補間法 　　- 粒子局在化法 ・ 適用例

適用例 （１） デバイス内での電子軌道計算 [Westermann, 1994]

適用例 （２） 金星電離圏 金星電磁圏のグローバルハイブリッド シミュレーション [Terada et al., 2002]

球面上への準一様な格子分布 [Tsugawa et al., 2000]

・ 局所的なグリッド集中による計算の効率化 が行われる。 非一様な格子分布の導入により、 　・ 複雑な形状場での適切な境界条件 　・ 局所的なグリッド集中による計算の効率化 が行われる。 さらに計算の効率化を行う、 　・ 粒子数の動的変更　（section 3) 　・ 非一様時間ステップ幅　(section 4) を説明していく。 非一様な格子分布と併用されることがしばしばある

３．粒子数の動的変更 非一様格子を用いる際には、粒子数の動的変更が 必要となる場合がある 　　　　　　　　 ・ 粒子分割 　　　　　　　　 ・ 粒子合体

粒子分割・合体の例 [Lapenta and Brackbill, 1994] ・ 粒子分割 　　　同一速度の２粒子に分割。 　　　ソース項が変化しないように粒子を分割する。 ・ 粒子合体 　　　運動量保存とエネルギー保存を共に満たすことはほぼ不可能。 　　　　　　　　　　　　　のものを合体。 合体により粒子数を１/２にした分布のテスト

Lapenta and Brackbill [1994] では計算効率が２倍ほど改善 しかし、 　　・ 合体は近似でしかない。 　　・ 分割において、細かいグリッドでの速度分布に偏りが 　　 出ないようにする為には、上流の粗いグリッド内にも 　　 十分な数の粒子が必要。 　　・ 分割・合体によって系の発展が変わる。 などの問題点がある。

グリッド幅 　→　 粒子数 グリッド幅 　←　 粒子数 イオン慣性長 　　グリッド幅 　　　　　　　　　　　　　　　　→　粒子が多い所では 　　　　　　　　　　　　　　　　　　 グリッドを細かくする

４．非一様時間ステップ幅 時間ステップ幅 Δｔ は、対象とするプラズマ現象を解くのに 十分小さくなければならない。 例えば一般的な陽解法ハイブリッドコードにおいては、 を共に満たす必要がある。

このことから、 などの場合は、小さな Δｔ が必要な領域でのみ小さな Δｔ を用いることによって、計算効率を大幅に改善する ことが期待できる。

非一様時間ステップ幅 １．背景場が空間変化する場合 例． Wu [1985] 地球磁気圏グローバルシミュレーションでは内部磁気圏領域で Δｔ を小さくしなければならない 　　・ 強い固有磁場のために　　　　　　　が小さくなる 　　・ 地球中心から 1RE 地点でのアルベン速度は ニ領域で異なる Δｔ を用いたスキーム

非一様時間ステップ幅 １．背景場が空間変化する場合 小さな時間ステップの間に速いグリッド側から 遅いグリッド側に伝搬していく短周期擾乱 ⇒　有限の伝搬速度の双曲型システムの場合、 　　 緩衝層（通常２～３グリッド）で対処

非一様時間ステップ幅 ２．グリッド幅が空間変化する場合 例． Lipatov [2002] ・ 格子幅は （n は格子分割のレベル） 　・ 時間ステップ幅は 　・ 先程と同様に緩衝層を用いて遅いグリッド側の場を補正 　・ 緩衝層では（計算の安定性のために）一時的に粒子分割を行うことも 非一様 Δｘ と非一様 Δｔ の併用

５．衝突効果 ・ 特定のエネルギー帯で衝突頻度が高い衝突過程 ・ 衝突以外の物理過程の介在 （例．惑星電離圏） など･･･ ⇒　プラズマ速度分布が非マクスウェル分布に 流体モデル　＋　衝突効果 粒子モデル　＋　衝突効果

[cf. Vahedi and Surendra, 1995] ＰＩＣコード　＋　モンテカルロ衝突 [cf. Vahedi and Surendra, 1995] １タイムステップ（Δｔ）ごとに衝突粒子を決定する方式を採用 （Δｔ 間の複数回衝突についてはテキスト 5.2.2 参照）

Δｔ 間の衝突確率 Ar+-Ar 衝突断面積

モンテカルロ衝突のアルゴリズム

モンテカルロ衝突のアルゴリズム （「空衝突」導入） ・ 全ての粒子についての計算は 3. のみに 下　図 ・ 全ての粒子についての計算は　 3. のみに ・ 4., 5. は N×Pnull 個の粒子に 　 ついて計算

境界適合格子 ・ 写像の概念を用いて、等間隔直交格子上（計算空間）で計算 を行う （詳細は後述） ・ 格子にそこそこの柔軟性はあるが、特定の領域に格子点を 極端に集中させるのは困難な場合が多い

複雑な形状の場合には、複数の格子を組み合わせて 流れ場を覆う方法が一般的に用いられる マルチブロック、重なり格子 マルチブロック法 重なり格子法 複雑な形状の場合には、複数の格子を組み合わせて 流れ場を覆う方法が一般的に用いられる

境界適合格子における変換式 (Appendix 1) という形の方程式は、計算空間 (x, h) では、 となる。 ここで である。 ここで 　　　　　　　　　　　　　　　である。 例えば、点 (i, j) での の差分式は などで数値的に評価しておく。

粒子コードへの適用 粒子コードに境界適合格子を適用するためには、 　　多数個の粒子量（ラグランジュ量）　と 　　非一様な格子点上の場の量（オイラー量） を効率的に結びつける方法が必要である。 例えば、粒子量（ 　　　　）からソース項（　　　）を 格子点上でいったん求めれば、場の量に関する 方程式については Appendix 1 の変換式がその まま使用できる。

粒子コードへの適用（２） 境界適合格子 ＋ 粒子コード の具体的な計算手順は、 １．粒子量（ ）を格子点上に集めてソース項（ ）を計算する 境界適合格子　＋　粒子コード　の具体的な計算手順は、 １．粒子量（　　　 ）を格子点上に集めてソース項（　　 ）を計算する ２．格子点上で場の量（　　　 ）を進める　（計算空間、 Appendix 1） ３．格子点上の場の量（　　　 ）から、粒子に働く力　　　を計算する ４．粒子量（　　　 ）を進める　（物理空間） 下線部分の粒子量と格子点上の　　　　　　　　 「補間法」 場の量とのやりとりを行う方法　　　　　　　　　　「粒子局在化法」

補間法（１）： 補間係数の求め方 通常の等間隔直交格子上での area-weighting method では、 粒子に働く場は格子点上の場 補間法（１）：　補間係数の求め方 通常の等間隔直交格子上での area-weighting method では、 粒子に働く場は格子点上の場 より、 で与えられる。

補間法（２）：　補間係数の求め方

粒子局在化法 「粒子追跡アルゴリズム」

金星電磁圏シミュレーション 金星電離圏 金星電磁圏のグローバルハイブリッド シミュレーション [Terada et al., 2002]