目次

第１章
序論
１
１．１ 電気推進機 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．１
１．２ マイクロ波放電型イオンエンジン ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．２
１．３ マイクロ波放電型イオンエンジンの課題 ．．．．．．．．．．．．．．．．．．．．．．．．．３
１．４ 本研究の目的 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．３
図 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．５
参考文献．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．１０

第２章
計算手法
１１
２．１
緒言 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．１１
２．２
プラズマ解析手法 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．１１
２．２．１ 粒子シミュレーション ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．１１
２．２．２ PIC 法 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．１３
２．２．３ 荷電粒子と中性粒子の衝突 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．１５
図．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．２１
２．３マイクロ波解析手法 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．２３
２．３．１ FDTD 法 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．２３
２．３．２ 安定条件 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．２５
２．３．３ 境界条件 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．２５
２．３．４ 給電方法 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．２６
２．３．５ FDTD コードの妥当性の検証 ．．．．．．．．．．．．．．．．．．．．．．．．．．．２７
図．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．２８
２．４
PIC 法と FDTD 法のカップリング ．．．．．．．．．．．．．．．．．．．．．．．．．．．．３１
２．４．１ カップリングコードの妥当性の検証 ．．．．．．．．．．．．．．．．．．．．．３１
２．４．２ カップリングコードを用いた計算 ．．．．．．．．．．．．．．．．．．．．．．．３３
図 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．３５
参考文献 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．４０

第３章

１０ｃｍ級マイクロ波放電型イオンスラスタ（MUSES-C model）に関する解析
４１
３．１ 緒言 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．４１
３．２ 計算体系 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．４１
３．３ 計算条件 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．４２
i３．４ 結果と考察 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．４２
３．４．１ グリッド（Z = 0.3 m）の境界条件を変化させた場合の計算結果の比較 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．４２
３．４．２ 電子の挙動解析 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．４３
３．５ 第 3 章のまとめ ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．４５
図 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．４７
参考文献 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．５８
第４章
マイクロ波放電型中和器に関する解析
５９
４．１ 緒言 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．５９
４．２ 計算体系 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．６０
４．３ 計算条件 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．６０
４．４ 計算結果と考察 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．６２
４．４．１ 電場強度の空間分布 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．６２
４．４．２ 電子エネルギーの空間分布 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．６２
４．４．３ 電子エネルギー分布関数（EEDF） ．．．．．．．．．．．．．．．．．．．．．．．６３
４．４．４ 電子の軌道とエネルギーの時間変化 ．．．．．．．．．．．．．．．．．．．．６３
４．４．５ 電離衝突発生位置 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．６３
４．４．６ 考察 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．６３
４．５ 第 4 章のまとめ ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．６４
図 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．６６
参考文献 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．７８
第５章
アンテナを利用したマイクロ波放電型プラズマ源に関する解析
７９
５．１ 緒言 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．７９
図 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．８１
５．２ 一様磁場中での解析 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．８２
５．２．１ アンテナと磁場が平行な場合 ．．．．．．．．．．．．．．．．．．．．．．．．．．．８２
５．２．２ アンテナと磁場が垂直な場合 ．．．．．．．．．．．．．．．．．．．．．．．．．．．８３
５．２．３ マイクロ波電場強度依存性 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．８３
５．２．４ 磁場強度依存性 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．８４
図 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．８５
５．３ 永久磁石の磁場中での解析 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．９２
５．３．１ 電子エネルギーの空間分布 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．９２
５．３．２ 電子エネルギー分布関数（EEDF） ．．．．．．．．．．．．．．．．．．．．．．．９３
５．３．３ 電離衝突発生位置 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．９３
５．４ 第 5 章のまとめ ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．９３
図 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．９５
ii

参考文献
第６章

．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．１０３

結論

１０４

謝辞

１０７

付録 A 小型マイクロ波放電型イオンエンジンに関する実験
１０８
A．１ 緒言 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．１０８
A．２ 実験装置 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．１０８
A．３ 実験結果と考察 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．１１０
A．４ 付録 A のまとめ ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．１１１
図 ．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．１１２
参考文献．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．．１１９

iii

第１章

序論

１．１ 電気推進機
１８９５年に Tsiolkowskii が「宇宙飛行とロケットについて」のなかで有名なロケット方程式を導出し，その約７０年後の１９６９年にはアポロ１１号が月面に降り立った．それからさらに３６年が経過した現在，人類の夢であった火星への有人探査 １−１）は，今まさに達成可能な目標となりつつある ．
アポロ 計画では全ての行程が化学推進によるものであったが，人類のより遠くへという願望を達成する手段として電気推進機が開発されてきた．化学推進では化学反応で得られる反応熱以上の排出速度を得る事は出来ないために，使用する推進剤で排出速度が決定してしまう．電気推進機では排出速度
は電磁気的に決定される為に，使用する推進剤 ではなく，使用できる電力に依存する．推進機の推進
剤利用効率は排出速度に比例するために電気推進機 は化学推進と比較してその効率が高い．推進機での推進剤利用効率は比推力 Isp（排出速度の 1/g，g : 重力加速度で定義され，sec の単位を持つ）という値で評価され，化学推進では Isp は最高で 500 sec，電気推進では Isp は 1000 ~ 10000 sec となる．高比推力であることは同一の速度増分（⊿V）を必要とするミッションにおいては化学推進と比較して推進剤が少量ですみ ，結果としてペイロード（目的地に到達できる燃料以外の衛星の質量）の増大が見込める．また，燃料を多く搭載する事でより遠方に到達できる事を意味している．これらの特徴から電気推進機は４０年程前から開発が進められ，現在使用されている代表的 な電気推進機はアークジェッ
トスラスタ ，MPD（Magnet-Plasma -Dynamic）スラスタ ，ホールスラスタ ，イオンエンジン 等である ．
これら 電気推進機はミッション によって使い分けがなされており，なかでも イオンエンジンは最も比推力が高く，長期間の人工衛星の姿勢制御や惑星探査に適している．
イオンエンジン は，アーク 放電やマイクロ波等で推進剤 を加熱・電離させて プラズマをつくり，２
枚ないし３枚からなる多孔状の電極に 1000 ~ 2000 Vの電圧を印加してイオンを加速するという静電加速型の推進機である．Fig. 1-1 にイオンエンジンの概念図を示す １−２）．イオンエンジンの物理過程はプ
ラズマ 生成，イオンビーム 引き出し，イオンビーム の中和に大別することができ ，プラズマの生成方法によって直流放電型，RF 放電型 １−３），マイクロ波放電型 １−４，５）に分類される．直流放電型はさらに Kaufman 型 １−６）とカスプ型に分けられる．この中でもっとも作動実績があるのは 直流放電型のカスプ型である．カスプ型の概念図 を Fig. 1-2 に示す．カスプ型におけるプラズマの生成はホロカソードと呼ばれる 電子源から電子を加速，抽出する事により行われる．抽出された 電子は放電室内に配置された永久磁石が形成する磁力線 に沿って運動し，ミラー 磁場の効果によって閉じ込められ往復運動をする過程で推進剤（おもに キセノン）を電離し，アノードへと到達する．カスプ 型の中で特筆すべき成果を残しているのは１９９９年に打ち上げられた「Deep Space-1」 １−７）である．このイオンエンジ

1

ンはこれまでに運用されている宇宙推進機の中で最長の稼働時間 16,000 時間 １−８）（２００１年１２月まで）を達成している．

１．２ マイクロ波放電型イオンエンジン
直流放電型の推進性能は既にほぼ 実用の域に達している．しかし ，エンジンの寿命，システムの信頼性といった観点からは改良すべき 点があるのも事実である．エンジンの寿命に関しては上述の直流放電型の寿命は 16,000~18,000 時間程度 であると考えられている ．これはホロカソードの損耗によるものであり，放電室で生成された 高エネルギーイオン の衝突によるグリッドや放電室壁の損耗も寿命を短くする要因である．システムの信頼性では，イオン引き出し以外にホロカソードの使用にはヒーター用の電源，電子抽出の為のキーパー用の電源が必要であり，これがシステム構成の複雑化 を招いている．また，ホロカソード の加熱に要するイオンエンジンの起動時間の長さも問題の一つとして挙げられる ．イオンエンジンは長期間作動させる 事によって高い⊿V を達成するのであるから，上述の問題点は極力排除する必要がある．
１９９４年に日本では宇宙科学研究所（現宇宙航空研究開発機構）より小惑星からのサンプルリターンが提案された．これは MUSES-C（Mu Space Engineering Spacecraft C（C は３番目のプロジェクトである事を意味する））計画と呼ばれ，火星の外側の小惑星帯へ探査機を送り，小惑星の表面の岩石（サンプル ）を採取して地球に帰還するという計画である．この計画は２００３年５月に実行に移され，探査機 「HAYABUSA」は２年間の宇宙航行の後に，２００５年１０月に小惑星「ITOKAWA」に到達している．MUSES -C 計画に求められる推進機の性能の一つとして往復約４年間（20,000 時間以上 ）の稼働時間が求められる．上述の直流放電型ではこの稼働時間の条件は保障されていない．そこで 宇宙科学研究所では，直流放電型に起因する故障要因を排除するために，MUSES-C 計画においては無電極
放電のマイクロ波放電型イオンエンジンを採用する事となった．マイクロ波放電型イオンエンジンでは ECR（電子サイクロトロン共鳴）を利用する為に無電極放電が可能である．ECR を利用したプラズ
マ源は材料プロセスや加速器等に利用されているが，イオンエンジンへの応用は MUSES-C 計画において初めて成された．マイクロ波放電型イオンエンジンの概念図を Fig. 1-3 に示す．マイクロ波放電型イ
オンエンジンはマイクロ波伝送系，磁気回路系，イオン 加速系，中和器系から構成される．マイクロ波はアンテナから放射され導波管を伝播し，放電室に到達する．放電室内には永久磁石によって ECRの条件を満たす磁場が生成されている．放電室内に存在する電子がマイクロ波によって共鳴的に加速されキセノンを電離し，電離によって新たに生成された 電子が再び電離を起こしプラズマが生成される．マイクロ波放電型イオンエンジンの推進性能は実用のレベルに達し，２００３年から上述のように実際の運用に至っている．また，耐久性能としては地上の耐久テストで 18,000 時間 １−９）を達成している．
また，イオンエンジンには衛星本体の電気的中性を保つために中和器と呼ばれる 電子源が必要であ2

る．直流放電型 ではこの中和器 にもホロカソードが使用されている．マイクロ波放電型では中和器にもマイクロ波放電型が採用されており，イオンスラスタ と同様に耐久性 には問題が無い事が確認されている．マイクロ波放電型中和器の概念図を Fig. 1-4 に示す．マイクロ波放電型中和器はマイクロ波を導入する L 型アンテナ，磁気回路から構成され，電子の抽出はイオンスラスタ本体と中和器の間の電位差によって行われる．イオンスラスタ本体とはプラズマの生成機構が異なり，アンテナが ECR 領域に直接挿入されており，アンテナとプラズマが接している．

１．３ マイクロ波放電型イオンエンジンの課題
マイクロ波放電型イオンエンジン は長寿命化が達成され実用レベルにあるが，推進性能では直流放電型に劣る事が現在の課題である．Table 1-1 に航空宇宙 技術研究所（現宇宙航空研究開発機構）で開発された直流放電型と MUSES-C で使用されているマイクロ波放電型の性能パラメータ １−１０）を比較して示す．この中でイオンエンジンの性能として 重要な項目は推力密度（推力をイオンビーム 引き出し面積で割った値），イオン生成コスト，推進剤利用効率である．推進剤利用効率ではその差は２％と僅かであるが ，推力密度ではマイクロ放電型は直流放電型 に比べて７０％と低い．また，引き出された
イオン １個がどれだけのエネルギー で作られたかを 示すイオン生成コストでは２倍と高い．マイクロ波放電型イオンスラスタでのプラズマの生成分布が Fig. 1-5 に示すように ECR領域に局在化しており，プラズマの生成が均一でない事が船木らの放電室プラズマの計測結果 １−１１）から分かっている．この事が直流放電型 と比べて性能が低い原因と考えられる．マイクロ波放電型イオンスラスタでは放電室
中央に導波管が配置されている 為に放電室中央でのプラズマ生成は不可能であり ，導波管を使用する
以上この問題は回避できないと考えられる．
この問題を解決する方法として ，上述したマイクロ波放電型中和器の特性が注目された．マイクロ波放電型中和器 でのプラズマ計測の結果からは ，カットオフ密度を超えるプラズマの生成が示唆されている １−１２）．アンテナを用いる事で導波管を排除することが可能となり，放電室内での均一なプラズ
マ生成の可能性が考えられる． また，中和器を小型マイクロ波放電型イオンスラスタとして作動させる研究も行われている １−１３，１４）．Fig. 1-6 に小型マイクロ波放電型イオンスラスタの概念図を示す．
小型マイクロ波放電型イオンスラスタでは中和器で採用されていた L 型アンテナではなく，イオン引
き出し面から見て対称なアンテナを採用している．

１．４ 本研究の目的
マイクロ 波放電型 イオンスラスタの 性能向上 や改良は 現在も 行われている １ − １５ ， １ ６ ）． 特 にMUSES -C に搭載されているマイクロ波放電型イオンスラスタに関しては，プラズマの特性を把握する為の実験が多数行われてきた．しかし，マイクロ波放電型イオンスラスタおいて 重要となる 電子とマ
イクロ 波の相互作用に注目した研究は実験により行う事が難しく積極的 に行われていない．これに対3し，数値解析による 試み １−１７，１８，１９）もなされているが，マイクロ波放電はマイクロ波，永久磁石 の作る磁場，プラズマが組み合わさった複雑な体系であるため，多くの仮定無 しに解析することが 困難となっていた ．そこで 本研究で は，電子 とマ イ ク ロ波の挙 動を PIC（ particle-in -cell） 法と FDTD
（finite-difference-time-domain : 有限差分時間領域）法により解析するコード（カップリングコード）を開発した．PIC 法を採用する利点は衝突を取り扱える事とエネルギー分布関数が得られることにある．
また，FDTD 法はマイクロ波のプラズマ中での減衰を再現できる事とアンテナの取り扱いが簡単であるという利点を有している．これによりマイクロ波放電型イオンスラスタにおける電子とマイクロ波の挙動を仮定を取り入れずに解析を行うことが可能となる．そこで ，本研究 ではカップリングコード
を利用し，
１．

マイクロ波放電型イオンスラスタおよび 中和器における電子とマイクロ波の相互作用を数値解析により明らかにする，

２．

アンテナを用いたマイクロ波放電型イオンスラスタの設計指針を得る，

３．

マイクロ波とプラズマの解析を簡便に行える数値解析コード を開発する，ことを目的として研究を行った．
本論文は以下のように構成される
第１章の序論では電気推進機の概要について説明し，マイクロ波放電型イオンエンジンの特徴と課題について述べた．そして本研究において行う数値解析の意義について述べた．
第２章の計算手法では本研究で使用した PIC 法と FDTD 法について説明し，両方式のカップリング
方法について述べる．また，コードの妥当性を検証する．
第３章ではカップリングコードを用い て 10cm 級 マイクロ波放電型イオンスラスタ（ MUSES -C
model）に関する解析を行った．
第４章ではカップリングコードを用いてマイクロ波放電型中和器に関する解析を行った．
第５章ではカップリングコード を用いてアンテナを利用したマイクロ波放電型プラズマ源に関する解析を行った．
第６章において結論を述べ，本論文についてのまとめを行う．

4

Fig. 1-1 Schematic of ion engine system

5

Fig. 1-2 Schematic of DC ion engine system

Fig. 1-3 Schematic of microwave discharge
ion engine system

6

Fig. 1-4 Schematic of microwave discharge neutralizer

Table 1-1 Comparison of thrust performance
between DC discharge and microwave discharge

7Fig. 1-5 Distribution of plasma density in the microwave discharge ion thruster
adopted in MUSES-C

Fig. 1-6 Schematic of miniature microwave ion engine

8

参考文献

１−１）

NASA ホームページ（http://www.nasa.gov/missions/solarsystem/bush_vision.html）.

１−２）

R. G. Jahn, Physics of Electric Propulsion, McGraw-Hill, New York (1968)

１−３）

G. Saccoccia, “European Electric Propulsion Activities and Progra mmes” IEPC-01-002 (2001)

１−４）

H. Kuninaka, K. Nishiyama, I. Funakai, Y. Shimizu and J. Kawaguchi, “Asteroid Rendezvous of
HAYABUSA Explorer Using Microwave Discharge Ion Engines”, IEPC-2005-10 (2005)

１−５）

J. E. Foster, T. Haag, H. Kamhawi, M. Pattercon, S. Malone, G. J. Williams, J. S. Sovey and
C. Carpenter, “The High Power Electric Propulsion (HiPEP) Ion Thruster”, AIAA 2004-3812(2004)

１−６）

H. R. Kaufman, “An Ion Rocket with an Electron Bombardment Ion Source”, NASA,TN D-585 (1961)

１−７）

Deep Space 1 ホームページ （http://nmp.jpl.nasa.gov/ds1/tech/sep.html）

１−８）

J. M. Sankovic and J. W. Dunning, Jr. , “An Overview of Current NASA Activities in Electric
Propulsion Technology”, IEPC-99-003 (1999)

１−９）

K. Toki, H. Kuninaka, I. Funaki, K. Nishiyama, Y. Shimizu, “Development Status of a Microwave
Ion Engine System for the MUSES-C Mission”, IEPC-99-139 (1999)

１−１０） 栗木恭一，荒川義博，電気推進ロケット入門，東京大学出版(2003)

１−１１） I. Funaki, H. Kuninaka, S. Satori and K. Toki, “Plasma Diagnostics and Numerical Modeling of a
Microwave Ion Engine”, AIAA-98-3478 AIAA-98-3341 (1998)
１−１２） I. Funaki and H. Kuninaka, “Overdense Plasma Production in a Low-power Microwave Discharge
Electron Source”, Jpn. J. Appl. Phys., 40, 2495 (2001)
１−１３） Y. Takao, H. Masui, T. Miyamoto, H. Kataharada, H. Ijiri and H. Nakashima, “Development of
Small-Scale Microwave Discharge Ion Thruster”, Vacuum, 73, 449 (2004)
１−１４） 中山宜典，船木一幸，國中均 ，中島秀紀， 小型マイクロ波イオン 推進機の実験的研究 ，
日本航空宇宙学会誌論文集 , 53, 461 (2005)
１−１５） K. Nishiyama, Y. Shimizu, T. Miyamoto, M. Fukuda, T. Nakai and H. Kuninaka, “Research andDevelopment Status Microwave Discharge Ion Thruster µ20”, IEPC-2005-055 (2005)
１−１６） Y. Nakayama, H. Hayashi, M. Cho, I. Funaki, H. Kuninaka, “Feasibility Study on High Specific
Impulse Ion Thruster with C/C Grid System”, IEPC-2005-233 (2005)
１−１７） 西山和孝，

電子サイクロトロン共鳴型イオンスラスタの放電プラズマ ，東京大学大学

院工学研究科博士論文(1997)

9

１−１８） 平山美晴，中北正人，鳴のシミュレーション

マイクロ波放電型イオンスラスタにおける電子サイクロトロン共，日本航空宇宙学会論文集，47，267(1999)

１−１９） 吉永睦生，佐鳥新，船木一幸，小野寺範義，石川芳男， マイクロ波中和器のプラズマ解析, 平成１０年度宇宙輸送シンポジウム (1999)

10

第２章

２．１

計算手法

緒言

この章では，本計算で使用したカップリングコード の詳細について述べる．プラズマを解析する手法としては，プラズマを電磁流体として取り扱う方法とプラズマを粒子として取り扱う方法に大別される．本計算では特に電子の放電室内での軌道，中性粒子との衝突などの電子の粒子としての挙動を詳しく解析する為に，プラズマを粒子として取り扱う粒子法を採用した．本計算 では粒子法 の中でも一般的 に良く使われている PIC（particle -in-cell）法を採用した．また，マイクロ波の解析には FDTD
（finite-difference-time-domain）法と呼ばれる電磁波の解析に最適な手法を採用した．FDTD 法はアルゴ
リズム が単純である事，境界条件の取り扱いが容易などの特徴を持つ．以下にプラズマ解析手法，マイクロ 波解析手法の順に説明を行い，両手法のカップリング方式を述べ，本計算 で使用するカップリ
ングコードについて妥当性を含めて説明する．

２．２

プラズマ解析手法

２．２．１

粒子シミュレーション

２．２．１．１

基礎式２ −１）プラズマの挙動を解析する際の基礎式として運動方程式を用いる．式を以下に示す．

mdv= q( E + v × B )
dt

(2.1)

dx
=v
dt

(2.2)

ここで， m は粒子質量， v は粒子の速度， x は粒子の位置， t は時間， q は粒子の持つ電荷， E は電場， B は磁束密度を表わす．一般的な実験室プラズマは 1 cm3 当りに 1010 〜1014 個が含まれ，この粒子全部を取り扱うことは最新のスーパーコンピュータ を用いても不可能であり ，実用的でない ．そこで粒子シミュレーション で取り扱う粒子として，実際の粒子の電荷質量比 を一定に保ったまま多数の粒子の電荷と質量を１つにまとめた超粒子を使用する事とする．超粒子数は，１セルに含まれる粒子数が 10 個以上になる事と取り扱う粒子数よりも超粒子の重みが少ない事（統計性から）から決定す る．

２．２．１．２

計算方法（2.2）式に leap-frog 法による時間中心差分を施すと，

11

x n+1 − x n= v n+1 2
∆t

(2.3)

となる．ここで，各変数の上付き添え字は時間に関するもので， n はタイムステップである．
（2.1）式に時間中心差分を施すと，

v n+1 2 − v n-1 2 q  n 1 n+1 2= E + v+ v n-1 2 × B n 
∆t
m
2


(

)

(2.4)

である．n- 1 2 はある時刻 n から ∆t 2 だけ以前の値であることを示す．同様に，n + 1 2 はある 時刻 nから ∆t 2 だけ進んだ値であることを示す．
プラズマ 粒子の運動を求めるために，（2.4）式において未知量である v

n +1 2について 解き，求めた

v n+1 2 と（2.3）式より新たな粒子の位置 x n + 1 を求めることができる．（2.4）式を変形すると，
Av n+1 2 = S(2.5)と書ける．ここで， α ， β ， γ ， δ を次のように定義する．

α=q∆ t n
Bx
2m

(2.6)

β=

q∆ t n
By
2m

(2.7)

q∆t n
Bz
2m
δ =1+α 2 + β 2 +γ 2
γ =

(2.8)
(2.9)

なお，ここでの E ， B は粒子位置におけるそれぞれの値を表わしている ．これらを用いて行列 A とベクトル S は，以下のように構成される．

 1

A= γ
− β


−γ
1
α

β 

−α 
1 

 v xm + cE xn + γv ym − β vzm 
 m

n
m
m
S =  v y + cE y + αv z − γv x 
 m
n
m
m
 vz + cEz + βv x − αv y (2.10)(2.11)

ここで， c = q∆t m ， m = n − 1 2 である．さらに（2.5）式を次のようにして変形し，速度を求めることができる．

v n +1 2 = A −1 S 1+ α 21−1
A =  αβ − γδ αγ + βαγ − β 
1+ ββγ + α βγ − α 1 + γ 2 (2.12)

αβ + γ212(2.13)

メッシュサイズ，タイムステップ２− ２）

２．２．１．３

計算に用いるメッシュサイズ ∆ は，電子温度 Te，ボルツマン定数 kb ，電子密度 n e から計算される
デバイ長λDe，

ε 0 k b Te
ne e 2λDe =(2.14)に対して，∆ ≤ λDe(2.15)であることが望ましい．また，タイムステップはプラズマ周波数 ? pe，nee2
me ε 0ω pe =(2.16)に対して，∆t ≤ ω pe−1

(2.17)とすれば，数値的な不安定性は生じない．時間積分の精度まで考えると，通常は∆t ≤ 0.2ω pe−1

(2.18)が用いられる．また，磁場中の荷電粒子のサイクロトロン運動にともなう現象をシミュレートするた
めには，電子サイクロトロン周波数 ? c

ωc =

eB
me

(2.19)

に対して，∆t << ω c

−1

(2.20)

である必要がある．

２．２．２

PIC 法２−３）

２．２．２．１

内挿計算

シミュレーションでは，電磁場量を空間格子点上で求める．そのため，計算粒子の電荷，電流を空間の格子点上に分配し，マクスウェル方程式を用いて電磁場を計算しなければならない ．また，粒子が感じる電磁場 を計算する際，格子点上に割り振られた 電磁場を内挿する必要がある．このような計算を行うため，粒子シミュレーションにおいて標準的な手法である PIC 法（Particle In Cell method）を

13

用いる．

(i, j, k ) を 求 め る ． 次 に ， Fig. 2-1 に 示 すまず，粒子の位置から近傍のセル座標

dx1,dx 2 ,dy1, dy 2, dz1,dz 2 を求める．その値を用いて，図中にあるような直方体 V1〜V8 の体積を求める．この体積比から粒子が感じる電磁場を内挿することができる．

(

)

(

)

具体的には，座標 x , y , z 上にある粒子が感じる電場 E x , y , z は，計算セルの体積を V ，図に示すような粒子によって仕切られた空間体積をそれぞれ V1〜V8 とすると，格子点上の重みは

An =

Vn
V

(n = 1 〜 8)

(2.21)となり，E (x , y , z ) = A1E (i , j , k ) + A2 E(i + 1, j , k ) + A3 E(i + 1, j , k + 1) + A4 E (i , j ,k + 1)
+ A5 E (i , j + 1,k ) + A6 E(i + 1, j + 1, k ) + A7 E (i + 1, j + 1, k + 1) + A8 E (i , j + 1, k + 1)
(2.22)と表わされる．

２．２．２．２

電荷密度および電流密度計算

電荷密度，電流密度はそれぞれ Fig. 2-1 に示されるセルの格子点上で定義され，セル中に N 個の粒()

子を代表した１個の超粒子が存在したとき，格子点 i, j , k 上に割り振られる電荷は ，

ρ (i , j , k ) =

V1
qN
V2

(2.23)となる．すべての粒子について和をとれば，格子点状での電荷密度が得られる．また，電流密度は

j = qnv(2.24)

より，

j (i , j , k ) =

V1
qNv
V2

(2.25)

PIC 法では(2.23)式から求まる電荷密度を元にポアソン方程式を解くことにより，静電場の効果をシ
ミュレート可能となる．静電場を考慮する為には(2.15)式を満たす事が必要である．(2.15)式の意味は，
実際のプラズマ中での電場はデバイ長で遮蔽される為，この効果を PIC 法において再現する場合は，デバイ 長以内のセルの電場の影響だけを考慮すればよい 事を示している．一般のプラズマの場合，デバイ長は λDe ~ 10 −4 m である ．このオーダーのメッシュサイズを採用すると汎用の計算機では限られた領域のみしか解析できず，本計算 で対象とする直径が 0.02 ~ 0.1 mの放電室長 でのマイクロ波の電場中での電子の挙動を解析することは不可能となる ．そこで本計算においてはメッシュサイズを(2.15)式より大きく設定し，静電場を考慮しない事とした ．それによりイオン は背景とし電子のみを取り扱うこ
ととする．

14

２．２．３ 荷電粒子と中性粒子の衝突
荷電粒子は電磁場 の影響を受けて，中性粒子と衝突する．プラズマは電子，イオン，中性粒子からなるとする．電子とイオンの数密度は中性粒子の数密度より非常に小さいため，電子-電子，電子-イオ
ン，イオン-イオンの各衝突は無視できるものとする．電子と中性粒子は弾性，励起，電離衝突をすると仮定する．なお，以下に用いる乱数 ２−４）は（0,1）間に一様に分布する乱数である．

２．２．３．１

平均自由行程と透過度

面 α ，厚さ l ，その物質の数密度 N の板を考える．エネルギー E の粒子ビームが，面 α に垂直に当るものとする ．衝突をせずに板を透過してくる粒子数は，入射粒子数より減少する．その減少度は σ を衝突の断面積とすれば，と表わされ，これよりdn = −nNσdl(2.26)

n = n0 exp(− Nσl )(2.27)となる．ここで， n 0 は面 α （ l = 0 ）でのビーム中の粒子である．それゆえ，（ l ,l + dl ）間で最初の衝突が起こる確率は，

Pinitial(l )dl = exp(− Nσl )Nσdl(2.28)となり， l 以下で最初に衝突が起こる確率は，P(l ) = ∫0 Pinitial(l )dll= 1 − exp(− Nσl )(2.29)()となる．最初の衝突までに粒子が走る平均距離 λ は P l の一次のモーメントである．
∞

λ = ∫ lP (l )dl0(2.30)

1=
Nσ
これを 平均自由行程という ．よって モンテカルロ法によって，任意の点から出発し最初の衝突が起こるまでの距離 l を決めるには，乱数 r を用いて，

r = P (l ) l 
= 1 − exp −  λ

で表わされ，これを変形して，

l = −λln (1 − r )(2.31)(2.32)となる．ここで， r は一様乱数であるから 1 − r を r とすれば，

l = −λln( r )15(2.33)

と表わされる．また，速さ v の粒子が気体分子と衝突する際の断面積を σ とすれば，この粒子が１秒間に気体分子と衝突する回数 υ は，

v
λ
= σNv

υ=

(2.34)となり，これを衝突周波数という．

２．２．３．２

衝突確率，衝突判定

すべての衝突の断面積は，荷電粒子の全エネルギーの関数として表わされる．

σ = σ (E )(2.35)

荷電粒子速度がわかれば，荷電粒子と中性粒子の各衝突 における衝突周波数 がわかり，全衝突周波数

υ total がわかる．したがって，ある粒子 m の衝突確率 Pcollision,m は以下の式で表わすことができる．
υ total = ngas σ total v m(2.36)

Pcollision,m = 1 − exp(− ngas σ total vm ∆t )= 1 − exp(− υ total ∆t )(2.37)

ここで， ngas は中性粒子の数密度 ， v m は粒子 m の速度， σ total は全衝突断面積 である ．ある乱数 r1（ 0 < r1 < 1 ）を用いて， Pcollision,m > r1 のときに粒子 m はタイムステップ ∆t の間に衝突が起こるとす
る．

２．２．３．３

衝突計算

乱数 r2 を用いて以下のように衝突の種類を決める．ここでは例として，電子と中性粒子の衝突について示す．弾性衝突断面積を σ elastic ，励起衝突断面積を σ excitation ，電離衝突断面積を σ ionizationとする．

∗ 0 ≤ r2 ≤σ elastic
σ totalのとき弾性衝突∗

σ elasticσ+ σ excitation≤ r2 ≤ elastic
σ totalσ totalのとき励起衝突∗

σ elastic + σ excitationσ+ σ excitation + σ ionization≤ r2 ≤ elastic
σ totalσ totalのとき電離衝突

ここで， σ total = σ elastic + σ excitation + σ ionization である．このように断面積の比と乱数を対応させ，衝突の種類を決める．計算では中性粒子としてキセノンを用い，その断面積データ ２−５，６，７）を Fig. 2-2 に示す．

16

２．２．３．４

エネルギーと散乱角

粒子は衝突の前後でエネルギー の受け渡しを行い散乱される．以下にそれぞれの 衝突におけるエネ
ルギーと散乱角について述べる．
①弾性衝突

(

)

(

)

まず，分子同士の弾性衝突を考える．速度 v1 ,v 2 を持つ分子対 の衝突後の速度を v1′ ,v ′2 とすると，

v1′ =

1
(v1 + v 2 − u ′)2

(2.38)

v ′2 =

1(v1 + v 2 + u′)2(2.39)と表わすことができる．ここで， u ′ = v ′2 − v1′ （相対速度）である．電子と中性粒子の衝突も同様に考える．衝突前の電子，中性粒子の速度をそれぞれ v e ， v n ，衝突後 の速度をそれぞれ v ′e ， v ′n とする．
電子と中性粒子の衝突では，電子は中性粒子に比べ質量が非常に小さいため，電子のみが散乱される．
中性粒子の速度 v n を無視すれば，相対速度は u = v n − v e ≅ −v e としてよい ．よって，電子の質量を m ，
中性粒子の質量を M とすれば，次式が成立する ２−８）．

v ′e =

m + M cos χM sin χ
ve +
h
m+ M
m+M

(2.40)

h の直交座標成分は
hy = −
hz = −

h x = ver cos φv ex v ey cos φ − ve vez sin φverv ex v ez cos φ − v e v ey sin φ(2.41)(2.42)(2.43)ver

ここで，v er = v ey2 + v ez2(2.44)

v e = vex2 + vey2 + v ez2(2.45)

であり， χ は散乱角， φ は任意の角である．（2.40）式はエネルギー損失を含んでいる．
次に，衝突前後のエネルギー収支を考える．衝突前の電子のエネルギーを Eincident,e とし，衝突後の電子のエネルギーを Escattered,e とすれば，エネルギー収支は以下のように表わされる．

Escattered,e = Eincident,e(2.46)

電子は，エネルギー Eincident,e が小さいとき等方散乱され，大きいときは主に前方散乱となる．これを考慮した散乱角 χ の確率密度を g (χ ) として 以下に示す ２−９）．ここで， Eincident,e の単位を eV とする ．

17

g (χ ) =E incident,e sin χχ
2 1 + E incident,e sin 2 ln 1 + Eincident,e 2(0 ≤ χ ≤ π )(2.47)

これを解けば，{

2 1 − (1 + Eincident,e ) r3
cos χ = 1 +Eincident,e}(2.48)となる．この値を（2.40）式に用いて衝突後の速度を求める．このとき，任意の角は

φ = 2πr4(2.49)

として与える．

②励起衝突
基底状態にある中性粒子をある準位に励起したとき，電子が失う励起エネルギーを Eexcitationとし，入射電子のエネルギーを Eincident,e ，励起後のエネルギーを Escattered,e ，励起後の速度を v~ とする．エネルギ
ー収支は，

Escattered,e = Eincident,e − Eexcitation(2.50)E~v = v e 1 − excitation
Eincident,e(2.51)と表わされ，励起後の速度はとなる．励起衝突を励起と弾性衝突に分けて考えれば，励起衝突でエネルギーを失った後に（2.51）式
中 v~ の速度で弾性衝突をすると考えることができる．
（2.40）式中 v e に v~ を用いて，衝突後の速度 v ′e を求める．このときの散乱角には（2.48）式を用いるが（2.48）式中の Eincident,e には（2.50）式から求められる Escattered,e を用いる．励起エネルギー Eexcitation は，キセノンの場合 8.34 eVである．

③電離衝突
入 射 電 子のエネルギーを Eincident,e ，散乱電子のエネルギー を Escattered,e ，生成電子のエネルギーを

Ecreated,e ，電離エネルギーを Eionizationとする．エネルギー収支はEscattered,e + Ecreated,e = Eincident,e − Eionization(2.52)と表わされる．
電離後の散乱電子のエネルギーは以下のように求められる ２−１０）．

E− Eionization  -1  incident,e 
Escattered,e = B tanr5 tan 
2
B



18

(2.53)

ここで，B はガスに固有の値でキセノンでは 8.7 e V となる．
また，入射電子が電離によって失うエネルギーを ∆ E とすると，∆E = Eincident,e − Escattered,e(2.54)と表わされる．
（2.51）式から電離後の電子の速度は，~v = ve 1 −∆E(2.55)E incident,eとなり，励起衝突の場合と同様に，この速度で弾性衝突をすると考える．
（2.40）式中 v e に v~ を用いて，衝突後の速度 v ′e を求める．散乱角には（2.48）式を用いるが，（2.48）式中の Eincident,e には（2.53）式から求められる Escattered,e を用いる．電離 エネルギー E ionization は，キセノンの場 合 12.13 e V で あ る ．

電離後の生成電子の持つエネルギーは，Ecreated,e = ∆E − Eionization(2.56)となる．また，入射電子の速度 v e を用いて生成電子の速度はv created,e =veve

2E created,em(2.57)と表わされる．この速度で弾性衝突をすると考えれば，v created,e を v e の代わりに（2.40）式に代入して，衝突後の速度を求める．このときの散乱角は生成電子のエネルギー Ecreated,e を用いて（2.48）式から求める．

２．２．３．５

Null-collision 法 ２ −１ １， １２）

計算時間短縮のため，Null-collision 法を用いる．これには，架空の衝突断面積 σ fake を導入する．すべての電子のエネルギー E に対して

σ fake < σ total(2.58)であるとし，例えば電子と中性粒子の衝突においてはσ total = σ elastic + σ excitation + σ ionization + σ fake(2.59)とする．ここで，σ totalv =υ collision= constantn gas(2.60)

を満たす σ fake を導入すれば ， υ collision が一定，すなわち Pcollision,m が一定となりエネルギーに依存しない．したがって，粒子毎 に Pcollision,m を計算する必要が無くなり，全粒子に対して衝突計算を行わず，衝突を起こす粒子のみについて 衝突の種類を決めればよい．これにより ，計算時間を大幅に短縮でき

19

る．衝突の種類を決めるときは，２．２．3．３節において，σ total = σ elastic + σ excitation + σ ionization + σ fakeとし，∗

σ elastic + σ excitation + σ ionization≤ r2 ≤ 1
σ totalのとき Null-Collision Processとする ．Null-Collision Process が選ばれたときには，粒子の運動は何ら変化しないものとする．Fig. 2-3に電子と中性粒子の衝突におけるキセノンガスの中性粒子密度 で規格化 した衝突周波数を示す．本計
算では Null-Collision 法で用いた υ collision n gas の値は 2. 00 × 10 -12 m 3 / sec である．

20

Fig. 2-1 Construction of cell in the PIC method

21Fig. 2-2 Collision cross section of electron with Xenon

Fig. 2-3 Collision frequency of electron with Xenon normalized by using Xenon density

22

２．３

マイクロ波解析手法２−１３，１４）

２．３．１

FDTD 法

２．３．１．１

基礎式

マイクロ波伝播解析には，FDTD 法（Finite Difference Time Domain method）を用いる．FDTD 法は数値電磁場解析の分野において，モーメント法や有限要素法と同様に，アンテナや平面回路など数多くの複雑な問題に適用され，主要な解析手法として多く用いられている．計算には，マクスウェル 方程式を基礎式として用いる．以下に基礎式として用いる微分形式のマクスウェル方程式を示す．

∂H
∂t
∂E
∇ × H = ε0
+J∂t
∇ × E = − µ0(2.60)(2.61)

ここで， E は電場， H は磁場， µ 0 は真空の透磁率， ε 0 は真空の誘電率 ， J は電流密度を表わす．

２．３．１．２

差分式

三次元の FDTD 法で用いられる電場，磁場の空間配置を Fig. 2-4 に示す．このように，電場，磁場ともに空間的 に交互に配置され，電場はセルの各辺に沿って，磁場はセル面の中心で与えられる．これは Yee によって提案され ２−１５），電場の回転（ ∇ × E ）が磁場を，磁場の回転（ ∇ × H ）が電場の方
程式を満たすように配置されている．
この電磁場の配置のもと ，（2.60）式，（2.61）式を空間，時間とも中心差分した式を以下に示す（（2.62）()

〜（2.67）式）．ここで， ∆x , ∆y , ∆z は空間幅， ∆t は時間幅，添え字 i, j , k はセル中での座標，添え字 n は時間ステップを表わす．
n+

1

n−

1

H x 2 (i , j + 1 2 , k + 1 2) = H x 2 (i , j + 1 2 ,k + 1 2 ) +∆tµ 0 (i , j + 1 2 , k + 1 2 ) E yn (i , j + 1 2 ,k + 1) − E yn (i , j + 1 2 ,k ) Ezn (i , j , k + 1 2) − Ezn (i , j + 1,k + 1 2) 
+


∆z
∆y


(2.62)
n+

Hy

1
2

n−

1

(i + 1 2 , j , k + 1 2) = H y 2 (i + 1 2 , j ,k + 1 2 ) +∆tµ 0 (i + 1 2 , j , k + 1 2 ) E zn (i + 1, j , k + 1 2 ) − Ezn (i , j , k + 1 2 ) E xn (i + 1 2 , j ,k ) − E xn (i + 1 2 , j , k + 1) 
+


∆x
∆z


(2.63)

23

n+

Hz

1
2

n−

1
2

(i + 1 2 , j + 1 2 , k ) = H z (i + 1 2 , j + 1 2 , k ) +∆tµ 0 (i + 1 2 , j + 1 2 , k ) E xn (i + 1 2 , j + 1, k ) − E xn (i + 1 2 , j ,k ) E yn (i , j + 1 2 , k ) − Eyn (i + 1, j + 1 2 , k ) 
+


∆y
∆x


(2.64)

E xn+1 (i + 1 2 , j ,k ) = E xn (i + 1 2 , j ,k ) +∆t
ε 0 (i + 1 2 , j , k )1
 n+ 12n+
 H z (i + 1 2 , j + 1 2 ,k ) − H z 2 (i + 1 2 , j − 1 2 , k )

∆y


n+

+

Hy

1
2

n+

1
2

(i + 1 2 , j , k − 1 2 ) − H y (i + 1 2 , j , k + 1 2 )∆z− J x (i + 1 2 , j , k )

(2.65)

E yn+1 (i , j + 1 2 ,k ) = E yn (i , j + 1 2 ,k ) +∆t
ε 0 (i , j + 1 2 , k )1
 n+ 12
n+
2(
(
)
H
i
,
j
+
1
2
,
k
+
1
2
−
H
i , j + 1 2 , k − 1 2)
 x
x

∆z


1
1
n+
n+

H z 2 (i − 1 2 , j + 1 2 ,k ) − H z 2 (i + 1 2 , j + 1 2 , k )+− J y (i , j + 1 2 , k )∆x


(2.66)

Ezn+1 (i , j ,k + 1 2 ) = E zn (i , j ,k + 1 2 ) +∆t
ε 0 (i , j , k + 1 2 )1n+
 n+ 12
 H y (i + 1 2 , j , k + 1 2) − H y 2 (i − 1 2 , j , k + 1 2)

∆x


n+

+

Hx

1
2

n+

1
2

(i , j − 1 2 , k + 1 2) − H x (i, j + 1 2 ,k + 1 2)∆y
− J z (i , j , k + 1 2 )

(2.67)（2.62）〜（2.67）式から，時間ステップが n + 1 2 での磁場 H

24

1
n+
2 は，その半ステップ前での電場

Enと１ステップ前の磁場 H
前の磁場 Hn+n−

12 から求めることができ，時間ステップが1
2 と１ステップ前の電場

n + 1 での電場 E n +1 は，半ステップ

E n から求めることができる．初期値については， H

n+

1
2は

Enか

ら計算されるので，アンテナの根元に E n だけを与える．

２．３．２

安定条件

計算結果の正確さを保つため，空間の増分値 ∆x , ∆y , ∆z は波長 λ に比べて十分小 さくなければならない．（ λ 10 以下）
時間ステップ ∆t は，Courant の安定化条件より

∆t ≤

1
1
1
1
v
+
+
2
2
(∆x ) (∆y ) (∆z )2(2.68)

を満たさなければならない ．この条件は極めて厳しく，わずかでも満足しなければ計算が不安定 にな
る．

２．３．３
２．３．３．１

境界条件
PEC(Perfect Electric Conductor )境界条件

境界面が金属のような完全導体である場合，境界に平行な電場成分を 0 とする PEC 境界条件を用いる．放電室内壁は導体を仮定しているので，この境界条件を適用する．

２．３．３．２

Mur 吸収境界条件 ２−１６）

開放領域を取り扱う場合，計算領域を仮想的 な境界で閉じなければならない．この境界を吸収境界といい ，この吸収境界が完全でなければ，反射波が解析領域に戻ってしまい ，誤差の原因となる ．こ
こでは，定式化が容易で広く使用されている Mur の吸収境界条件を用いる．

x = 0 の吸収境界に x の正の方向から Ez 成分を持つ平面波が入射する場合を考える．伝播速度を vとすると， x の負の方向に進む波は

Ez = Ez ( x + vt )(2.69)∂E z 1 ∂Ez−=0∂x v ∂t(2.70)と表わされ，次の微分方程式を満足する．

x = 0 の境界で完全に反射がないならば ，電場は（2.69）式の形を保ったまま伝播するはずである．すなわち，境界でも（2.70）式が満足される．これを FDTD 法の形式に定式化するため，時間について差分化すると25

∂E z Ezn − E zn −1
=
∂t
∆t
n−

=v

∂E z
∂x

(2.71)

1
2となる ．右辺の x についての微分は，電場が x = 0, ∆x に割り当てられていることから， x = ∆x 2 で行い，（2.70）式に代入すると

E (1 2 ) − E zn−1 (1 2 )
n
z

∆t
n

( )

n −1

となる ．
（2.72）式の Ez 1 2 , E z

n−

=v

Ez

1
2

n−

1

(1) − Ez 2 (0)∆x(2.72)(1 2 ) も FDTD 法では直接求めることはできないため，前後の値で平均して求める．そして， Ez についてまとめるとEzn (0, j ,k + 1 2 ) = E zn −1 (1, j , k + 1 2 ) +v∆t − ∆x n
Ez (1, j , k + 1 2) − Ezn−1 (0 , j , k + 1 2 )v∆t + ∆x{}(2.73)となる．これを Mur の１次吸収境界条件という．上式は吸収境界に平行な電場成分について適用され
る．

２．３．４

給電方法

本計算 における給電方法は同軸ケーブルを模擬し，Fig. 2-5 に示すように動径方向の４つの電場を中心導体に与える同軸線路給電を用いる．また，アンテナ導入部には Fig. 2-6 に示すように，Mur の一次吸収境界条件を与え，中心導体には PEC 境界条件を与えた．
中心導体に与える電場の x 方向成分 Ex ，y 方向成分 Ey は近似的に

Ex =V
 ∆x 
rln  
 a 

Ey =

V
 ∆y 
rln   a 

(2.74)(2.75)

で与えられる．ただし ， V は中心導体と外部導体間 の電圧， r は動径方向の距離， a はアンテナの半径である．
また，同軸ケーブルのインピーダンス Z は以下の式で表される ．

Z=

µ 1
R
ln
ε 2π r

(2.76)

ここで， R は外部導体の半径， r は中心導体の半径である．マイクロ波電力を P とすると，

P=

V2
2Z

26

(2.77)

となり ，（2.76）式，（2.77）式より V を求め，（2.74）式，（2.75）式よりアンテナに与える電場を求める．ここで， Ex を r = ∆x 2 の電場で， Ey を r = ∆y 2 の電場で代表させると，

E x (IFED , JFED, KFED) = −E x ( IFED − 1, JFED, KFED)=−

V (n∆t )∆x  ∆x 
ln  
2  a (2.78)E y (IFED , JFED, KFED) = − E y (IFED , JFED − 1, KFED)=−(

V (n∆t )∆y  ∆y 
ln  
2  a (2.79))

ここで，式中の IFED, JFED, KFED は電場を励振する座標を与える．

２．３．５

FDTD コードの妥当性の検証

ここでは本研究で使用した FDTD コードの妥当性を検証する為にエネルギー保存の確認を行った．
プラズマが無い場合の電磁気的エネルギーの保存は，

∂ 1
1
2
2
(
E
×
H
)
⋅
d
A
=
ε
E
+
µ
H
0
0
∫
∫
dV
A

∂t V  22



(2.80)

で与えられる．
（2.80）式の左辺の E x H はポインティングベクトルと呼ばれ単位時間に A に垂直な方向へ流れる電磁波のエネルギー を示している．右辺は系の内部に蓄えられるエネルギーの時間変化を示している．検証は Fig. 2-7 に示す体系で行った．真空領域の端の破線は吸収境界条件を示している．
計算条件としては，メッシュサイズが X，Y，Z 方向とも 1 x 10-3 m，タイムステップが 1 x 10-12 sec，マイクロ波周波数が 2.45 GHz である．アンテナ前方の 8 セル x 8 セル x 8 セル （Calculation area）の領域について（2.80）式を適用した．Fig. 2-8 に計算結果 を示す．両者は時間的にもエネルギーの振幅においても良い一致を示している事が分かる．これにより本研究で使用した FDTD コードの妥当性が確かめられた．

27

Fig. 2-4 Arrangement of electric field and magnetic field in FDTD method28

Fig. 2-5 Supply of microwave by using coaxial line method

Fig. 2-6 Model of antenna

29

Fig. 2-7 Calculation system for checking the validation of the FDTD code

Fig. 2-8 Comparison of time variation of internal energy with energy flux
from results of the FDTD code

30

２．４

PIC 法と FDTD 法のカップリング

本研究では，２．１で述べたプラズマ解析の PIC 法と２．２で述べたマイクロ波解析の FDTD 法を組み合わせる事 ２−１７，１８）によりプラズマの挙動とマイクロ波の伝播を同時に解析する．Fig. 2-9 にカップリングコードのフローチャートを示す．FDTD コードで計算されたマイクロ波の電磁場 は PIC コ
ードに渡され，その電磁場 の効果が粒子の運動に影響する．また，電子の挙動により計算される 電流密度が FDTD コードに渡され，電磁波に粒子の効果が加えられる．このように PIC コード，FDTD コードの間で物理的に相互作用を及ぼすデータをやり 取りすることで ，プラズマおよびマイクロ波の挙動を解析する．なお，タイムステップについては PIC 法では（2.17）式により決定され，FDTD 法では（2.68）により 決定され計算の安定化条件が異なる．必ずしも両者のタイムステップを一致させる必要は無いので，本研究では PIC 法の計算負荷を軽減するために FDTD コードを連続（２〜２０回）で計算させた後に PIC コードを１回計算する事とした．カップリングコード は 3 次元 X − Y − Z 座標系を採用する事により，3 次元的構造を持つアンテナの解析を可能とした．
データの出力方法について，特に電子エネルギーの空間分布 Ee については エネルギー分布関数の平均から算出するのではなく，PIC 法のアルゴリズムを応用した方式を採用した．エネルギー分布関数から平均を算出する場合は出力する領域ごとに，全粒子に対して位置に関しての判定を行う必要がある．この方式で Ee を出力する領域が多い場合には計算速度の低下につながる ．これに 対して PIC法を利用した方式の場合の Ee の計算は次式より求まり，Ee (m ) = ∑ Eei Ai
i

∑A

i

(2.81)iとなる．ここで Eei は i 番目の電子のエネルギー，ｍは格子点の位置，Ai は(2.21)式で計算される格子点 m への重みを表している． この方式では全粒子に対して１回だけ，どのセルに含まれるかを判定す
れば全ての領域の Ee が得られるという利点がある．

２．４．１

カップリングコードの妥当性の検証

ここではカップリングコード の妥当性の検証を行うために，分散関係式との比較を行った．プラズ
マ中でのマイクロ波の振舞いはプラズマの密度，磁場成分そしてマイクロ波の周波数によって複雑に変化する．しかも各パラメータ ーにより，マイクロ波の伝播領域，非伝播領域が存在し，その境界ではカットオフや共鳴といった現象が起こる．ここでは，カップリングコード を用いて， B = 0 でのマ
イクロ波の伝播についてシミュレーション を行い，その結果と分散関係式との比較，考察を行う．また，プラズマに一様な磁場を与え，その強度を変化させた場合の電子エネルギーの変化を調べた．

２．４．１．１

B = 0 における円筒プラズマの分散関係式

本計算 では円筒プラズマを取り扱うが，この場合の分散関係式としては，文献２−１９）に詳しく

31

述べられている．
変数 S , D , P , R , L を，マイクロ波周波数 ω ，プラズマ周波数 ω pe および サイクロトロン周波数 ω c
を用いて，

S =1−

ωp2

(2.82)

ω 2 − ωc2

ω p ωcD=ω(ω 2 − ωc2 )
2

P =1−

(2.83)

ω pe2(2.84)

ω2(ω pe ω )R = S + D = 1−
1 + (ωc ω )(2.85)(ω pe ω )L = S − D = 1−
1 − (ω c ω )(2.86)

22と表すと，プラズマ波動の分散関係式は

[

]

(

)(

)

PN||4 + (S + P )N ⊥2 − 2PS N ||2 + N ⊥2 − P SN ⊥2 − LR = 0(2.87)

で与えられる．ここで N || , N ⊥ は磁場に平行および垂直な方向の屈折率を表し，

N\\ =
2

c 2k ||

2(2.88)

ω2
2
c2 k⊥
2
N⊥ =
ω2

(2.89)である．
一様密度の円筒プラズマを考え， k z ,θ を z 方向の波数，位相とするとき，電磁場をE = E( r )e i( k z z + mθ −ω t )B = B( r )e i( k z z +mθ −ωt )(2.90)(2.91)

で与え，波動方程式を解くと， Er , E? , Br , B? は E z および Bz で与えられ，次式を満たす ２−２０）．

(∇ + k )(∇ + k )P = 0
2
⊥

2
⊥1

2
⊥

2
⊥2

z

(2.92)

ここで，∇ 2⊥ =

2⊥1

また， k ,k1 ∂  ∂  m2
r  −r ∂r  ∂r  r 2
Pz = iωBz2は以下のコールドプラズマ分散関係式の解である．
⊥232(2.93)(2.94)

k ⊥4 − p0 k ⊥2 + q 0 = 0(2.95)

ここで，p0 ,q 0 は(2.96)~(2.97)式を用いて表され，

P
1
F+
S
P

(2.96)

P
ω 4D2  F − 2 Sc 

(2.97)

ω 2S− k z22c(2.98)

p0 = F +

q0 =

F=
である．

円筒導波管における電場の境界条件は，導波管の半径を a とすると，Eθ ( r = a ) = 0 ，Ez ( r = a ) = 0
で与えられるので，この境界条件をもとに(2.92)式を解くと，

Dm∂J
J m ( k ⊥ a) + ( N ⊥2 − S ) m ( k ⊥ a) = 0r∂rが得られる．ここで， J m ( k ⊥ a ) は m 次のベッセル関数であり ， k ⊥ =(2.99)N ⊥ω
である ． B = 0 を (2.83)c∂J m( k ⊥ a) = 0 が得られる ．前者はカットオ∂rフを与える式なので，後者について考え，本研究で取り扱うマイクロ波のモード数 m = ±1 の解であるk ⊥ a = 1.83 が得られる ．これを式に代入すると D = 0 となり ，解として N ⊥ = S および2

N ⊥2 + N||2 = P(2.100)に代入することで， B = 0 における円筒プラズマでの分散関係式が得られる．

２．４．２

カップリングコードを用いた計算

はじめに，カップリング コードを用いて円筒導波管内 におけるプラズマとマイクロ波の分散関係について計算した結果について述べる．計算体系を Fig. 2-10 に示す．マイクロ波は導波管左端で反射，右端で吸収境界条件を用いる．計算条件を Table 2-1 に示す．磁場は計算領域すべてで B = 0 である．
計算ではプラズマ密度を 0, 1.0×1016 , 3.0×1016 m-3 と変化させマイクロ波の伝播解析を行なった．マイク
ロ波の波数は X = 0.023 m, Y = 0.023 m（計算体系中央 ）における Ey 成分をプラズマ領域（Z = 0.26 ~ 0.36m）内で時間平均して算出した．理論値には，２．４．１．１で述べた（2.100）式の分散関係式 を用いる．カップリングコードによる計算結果および理論値を Fig. 2-11 に示す．図からわかるように計算結果が B = 0 における分散関係の理論値と良い一致を示すことが確認できる．
また，プラズマ密度が 1.0×1017 m-3 の場合の計算も行った．Fig. 2-12 にマイクロ波の電場強度分布を示す．この図から導波管内に励起された波が， z = 0 .26 m 以降のプラズマの領域で大きく減衰してい33

ることが分かる．カップリングコードは円筒導波管内のプラズマ中の マイクロ波の伝播，遮断を良く再現している．
次に，磁場強度を変化させた場合の計算結果について述べる．プラズマ領域に一様な磁場を Bx を与え，磁場強度を ω c ω = 0.99（ Bx = 0.15 T）および ω c ω = 0.33 ( Bx = 0.05 T)と変化させ電子エネルギーの上昇を比較した．Fig. 2-13 は ω c ω = 0.99 の場合，Fig. 2-14 は ω c ω = 0.33 の場合の計算結果を示している． ω c ω = 0.99 の条件は共鳴の条件であり， ω c ω = 0.33 の条件の結果と比べ電子エネルギーの上昇が大きい事が分かる．この結果からカップリングコードは共鳴現象をよく再現している事が分かる．
円筒導波管 内のプラズマの分散関係式との比較，及び磁場強度を変化させた 場合の計算結果からカップリングコードの妥当性が確認された．

34

Fig. 2-9 Flow chart of coupling code35

Fig. 2-10 Calculation system for for checking the validation of the coupling code

Table 2-1 Calculation condition36

Fig. 2-11 Comparison of numerical results with analytical results for B = 0

Fig. 2-12 Distribution of microwave electric field intensity for B = 0 and n = 1. 0 × 1017 m-3

37

Fig. 2-13 Electron energy distribution of ? c/?=0.99

Fig. 2-14 Electron energy distribution of ? c/?=0.33

38

参考文献
２−１）

内藤裕志， 粒子シミュレーションの基礎, プラズマ・核融合学会誌, 74, 470 (1998)

２−２）

田中實，山本良一， 計算物理と計算化学 ，海文堂 (1998)

２−３）

石黒静児， 静電粒子シミュレーション

２−４）

丹慶勝市，奥村晴彦，佐藤俊朗，小林誠，“Numerical Recipes in C （日本語版） ，, プラズマ・核融合学会誌 , 74, 591 (1998)

技術評論社 (1993)
２−５）

M. Hayashi, “Determination of Electron-Xenon Total Excitation Cross-Sections, From Threshold to
100 eV, From Experimental Values of Townsend’s α” J. Phys. D: Appl. Phys., 16, 581 (1983)

２−６）

M.S. Dababneh, W. E. Kauppila, J.P. Downing, F. Laperriere, V. Pol, J. H. Smart and T. S. Stein,“Measurements of total scattering cross sections for low-energy positrons and electrons colliding
with krypton and xenon” Phys. Rev. A, 22, 1872 (1980)

２−７）

J. J. Szabo, Jr., M. Martinez-Sanchez and O. Batishchev, “Numerical Modeling of the Near-Anode
Region in a TAL Thruster”, AIAA 00-3653 (2000)

２−８）

南部健一， コンピュータアナリシスシリーズ 7 原子･分子モデルを用いる
数値シミュレーション

２−９）

日本機会学会編，コロナ社 (1996)

C. B. Opal, W. K. Peterson and E. C. Beaty, “Measurements of Secondary-Electron Spectra
Produced by Electron Impact Ionization of a Number of Simple Gases”, The Journal of Chemical
Physics, 55, 4100 (1971)２−１０） M. Srendra, D. B. Graves and G. M. Jellum, “Self-Consistent Model of a Direct-Current Glow
Discharge : Treatment of Fast Electrons”, Physical Review A, 41,1112 (1990)２−１１） V. Vahedi and M. Surendra, “A Monte Carlo Collision Model for the Particle -in-cell Method:
Applications to Argon and Oxygen Discharges”, Computer Physics Communications 87, 179(1995)
２−１２） 電気学会・マイクロ波プラズマ調査専門委員会， マイクロ波プラズマの技術 ，オーム社(2003)
２−１３） 橋本修，安部琢美， FDTD 時間領域差分法入門 ，森北出版 (1996)２−１４） 宇野亨 ， FDTD 法による電磁界およびアンテナ解析 ，コロナ社 (1998)
２−１５） K. S. Lee, “Numerical Solution of Initial Boundary Value Problems Involving Maxwell’s Equationsin Isotropic Media”, IEEE Tarns. Antennas Propagation, 14, 302 (1966)２−１６） G. Mur, “Total-field Absorbing Boundary Conditions for the Time-Domain Electromagnetic Field
Equations”, IEEE Trans. Electromagnetic Compat., 40, 100 (1998)
２−１７） V. P. Gropinath and T. A. Grotjohn, “Three-dimensional Electromagnetic PIC Model of a Compact
ECR Plasma Source”, IEEE Trans. Plasma Sci., 23, 602 (1995)
２−１８） T. A. Grotjohn, “Modeling the Electron Heating in a Compact Electron Cyclotron Resonance Ion
Source”, Rev. Sci. Instrum. 67, 921 (1995)
２−１９） R. C. Cross and D. Miljak, “Effects of electron mass on Alfven waves in a cylindrical plasma,”
Plasma Phys. Control. Fusion, 35, 235 (1993)39

２−２０） T. H. Stix, “Waves in Plasmas”, Springer-Verlag, New York (1992)

40

第３章

10cm 級 マ イ ク ロ 波 放 電 型 イ オ ン ス ラ ス タ（MUSES-C model）に関する解析

３．１

緒言

この章では宇宙科学研究所（現宇宙研究開発機構）によって開発され，MUSES -C 計画で使用されているマイクロ波放電型イオンスラスタについてカップリングコード を使用して解析した結果について述べる．マイクロ波放電型イオンスラスタに関する数値解析は平山 ３−１）や船木 ３−２）によって行われているが，平山の場合はマイクロ波を減衰しない波として扱う仮定を導入している．船木の場合は FDTD
法を用いているが，誘電テンソ ルを導入する事によってプラズマのモデル化を行っている．これら２つの研究ではマイクロ波と電子の同時解析は行われていない事が課題であった．本研究で使用するカップリングコードではマイクロ波と電子の同時解析が可能である．本章ではグリッド部におけるマイクロ波の境界条件として吸収条件と反射条件の比較を行った ３−３)．そして EEDF（電子エネルギー分布関数）と電子の軌道の計算結果からマイクロ波放電型イオンスラスタにおけるプラズマ生成の特性について議論する ３−４)．

３．２ 計算体系
Fig. 3-1(a)に計算体系を，Fig. 3-1(b)に放電室内 の磁場強度分布を示す．マイクロ波はモノポールアン
テナにより供給される．マイクロ波の周波数は 4.25GHz である．アンテナの長さはアンテナからのマ
イクロ波の放射を最適化する為に 4.25GHz のマイクロ波の波長（0.071 m）の 1/4 程度になるように設計されている．その為に円筒導波管内に挿入されているアンテナの長さは 0.018 mである．円筒導波管は直径 0.046 mで長さは 0.26 m である．Z = 0.26 m以降が放電室となり，直径が 0.045 m から 0.12 m
へと広がるホーン型の形状になっている．放電室内にはリング状の Sm-Co 磁石が２列配置されており，
ECR に必要な磁場が生成される ．
マイクロ波周波数 4.25GHz に対応する ECR 領域の磁場の強さは 0.15
T である．磁場は文献３−５）の方式より算出した．
本計算では領域ごとに計算コードの取り扱いが異なり，放電室の入り口（Z = 0.26 m）までは FDTD
コード によりマイクロ波のみを 計算し，放電室内ではカップリングコードにより マイクロ波とプラズ
マが同時解析される．メッシュサイズについては X，Y 方向については 2 x 10-3 m であるが，Z 方向に
ついては円筒導波管内では 6 x 10-3 m とした．アンテナ付近ではマイクロ波の変化が大きく，体系全体に与える影響が大きいのでメッシュサイズをより細かく設定する必要があり 2 x 10-3 m とした．放電室
内も円筒導波管からホーン型への形状の変化が大きい為にメッシュサイズを 2 x 10-3 m とした．メッシ41ュ数は X 方向に 65，Y 方向に 75，Z 方向に 78 である．

３．３ 計算条件
計算条件を Table 3-1 に示す．FDTD コードに対するタイムステップは第２章で述べた Courant の安定化の条件より 1.0 x 10-12 sec とする． PIC コードに対しては 2.0 x 10-11 sec とし，FDTD コードにより２０回連続でマイクロ波を計算したのちに PIC コードによりプラズマを１回計算する事とした．推進剤としてはキセノンを想定し，船木らの実験 ３−６）から放電室内の圧力は 6.65 x 10 -3 Pa(5 x10 -5 Torr)と報告されているので，中性粒子密度は 1.604 x 1018 m-3 とした．初期の電子密度は 1.0 x 1016 m-3 とし，放電室内に一様分布とした．計算では 4,000,000 個の粒子が取り扱われ，それぞれの粒子は 900,000 個を代表する超粒子である．電子の初期エネルギーは 0 eVとした．
粒子に関する境界条件としては磁石，放電室壁に達した電子は吸収され速度と電荷共に 0 とする，
すなわち計算の対象から除かれるとして取り扱った．マイクロ波に対する境界条件として円筒導波管
内では反射境界条件とし，Z= 0.3m のグリッド 部には吸収境界条件と反射境界条件の二つの場合で計算を行った．

３．４ 計算結果と考察
３．４．１

グリッド（Z = 0.3 m）の境界条件を変化させた場合の計算結果の比較

はじめにプラズマが無い場合のマイクロ波の伝播についての解析を行った．Fig. 3-2(a)，(b)に t = 1.25
nsec における ZY 断面（X = 0.064 m）でのマイクロ波の電場強度分布を示す．アンテナから放射されたマイクロ波は t = 1.25 nsec に Z = 0.3 mの境界に達する．Fig. 3-2 (b)の反射境界条件では吸収境界条件の場合と比較して，反射したマイクロ波が入射波と重なり電場強度が Z = 0.295 m付近で上昇している．
Fig. 3-3 (a)，(b)には t = 3 nsec における電場強度分布を表している．3 nsec は 4.25 GHz のマイクロ 波周
波数では約１３周期分に相当し，定常状態に達していると考えられる．Fig. 3-3 (a) の吸収境界条件の計算結果ではマイクロ波は進行波として伝播している事が分かる．反射境界条件を課した Fig. 3-3 (b)の場合，計算領域はすべて 金属壁で囲まれていることになる．一般に金属壁 で囲まれた領域内に投入された 電磁波は反射を繰り返し，電磁波の波長と金属壁間の距離が最適な場合は定在波が形成され高い電場強度が得られる共振器（Resonance cavity）として機能する事が知られている ３−７）．本計算の場
合は金属壁 で囲まれた領域が円筒型 からホーン 型へと変化するような複雑な形状の為に共振器としては機能しない事が分かる．
次に電子を放電室に配置した場合の結果について述べる．Fig. 3-4(a)，(b)に吸収境界条件の場合の t =100 nsec 後の電子エネルギー Ee の空間分布を示す．Fig. 3-4(a)から電子エネルギーの上昇が ECR 領域
のみで起こることが確認できる．反射境界条件の計算結果である Fig. 3-5(a)，(b)の場合でも電子エネル

42ギーの上昇は ECR 領域のみでしか起こらない．両者の結果を ZY 断面（X = 0.064 m）における電子エ
ネルギーの空間分布の時間毎の平均値 eave によって比較する．Z-Y 平面における電子エネルギーの平均値 eave は次式により算出した．
ny

nz

j

k

eave = ∑ ∑ E e (j, k)/ ny ⋅ n z(3.1)

ここで ny，n z は Fig. 3-5，6 の Y と Z 方向のグリッドの数を表している．また， Ee (j,k)は YZ 平面における電子エネルギー平均値である ．この eave を時間ごとにプロットしたものが Fig. 3-6 である．この図からも 反射境界条件の場合が電子のエネルギーの上昇が大きい事がわかる．両者の計算結果を比較す
ると，Ee は吸収境界条件においては 1 1 eVほどであるが ，反射境界条件では 1 6 eVと高い事が分かる．
また，吸収境界条件では最も Ee が高い場所は Z = 0.28 m付近であるが，反射境界条件では Z = 0.29 mとグリッド側へ移動している．これは,グリッドによって反射されたマイクロ波の電場強度の強い位置
が，吸収境界条件の場合に比べ前方に移動した事によるものと 考えられ る．この結果はグリッド側でプラズマが生成される可能性を示しており，引き出されるまでのイオン の損失の低減による 性能の向上につながるものと考えられる．

３．４．２

電子の挙動解析

次に計算時間を 800 nsec まで延長した場合の電子の挙動について述べる．ここでは特に EEDF（電子エネルギー 分布関数）と電子の軌道，電子エネルギーの時間変化に着目した．直流放電型イオンスラ
スタでは Brophy モデル ３−８）と呼ばれる放電プラズマのモデルが確立している．本研究の対象である
マイクロ波放電型イオンスラスタでは放電プラズマのモデル は確立しておらず，現状では Brophy モデ
ルを基本とした モデル化が行われている．直流放電型イオンスラスタにおいてモデル化が可能な理由
は，ホロカソードから抽出され電離に寄与する電子の電子電流と新たに生成された電子の電子電流とを明確に区別できる事による．マイクロ波放電型イオンスラスタでは電離に寄与する電子と生成される電子を明確に区別できないためにモデル化が困難である．そこで重要と考えられるのが EEDF である．EEDF を求める事により，電離に寄与する電子，すなわちキセノンの場合では 12.13 e V 以上のエネルギーを持つ電子の存在率が分かればマイクロ波放電型イオンスラスタの放電プラズマのモデル確
立につながると考えられる．実験により直接 EEDF を計測する事（例えば Druyvesteyn 法による）も試みられているが，ECR 領域内にプローブを挿入して計測を行う場合は，プローブ挿入がプラズマ生成に大きな影響を与える事が問題となる．本研究 では粒子法を採用した事の利点を生かし，マイクロ波
放電型イオンスラスタにおける EEDF を計算により求めた．また，マイクロ波放電型イオンスラスタ
における電子の軌道に関しても 多粒子の軌道を追跡する事により放電室内での電子の挙動の解析を行った．

43

３．４．２．１

EEDFの解析結果

まず,Fig. 3-7 に反射境界条件の下での 800 nsec での電子エネルギー Ee の空間分布を示す．Fig 3-5（100 nsec）の場合と比べて ECR 領域での電子エネルギーは 20 eVを超えており，12.13 eV以上のエネ
ルギーを持つ電子は多数存在していると予想される．次に，EEDF の解析対象とした領域を Fig. 3-8 に示す．本計算では ECR 領域１と ECR 領域から Z+方向に 0.02 m離れた領域２との比較を行った．それ
ぞれの領域は円筒型で約 90,000 個の電子が含まれている．Fig. 3-9 に EEDF の計算結果を示す． 12.13
e V 以上のエネルギーを持つ電子の存在率は領域１は 10.9 %，領域２で 2 .4%となり，領域１では領域２の約４倍の確率で 12.13 eV以上のエネルギーを持つ電子が存在する．また，分布の形状に関しては，
分布のピークが高エネルギー側に移動せずに，高エネルギー部分が伸びる（テールが形成される ）傾向にある．これは本計算で用いた圧力が 6.65 x 10 -3 Pa と低い事により衝突頻度が低く，衝突による電子のエネルギー の損失が少ないことが原因である．結果として衝突によりエネルギーを失わなかった
電子が高エネルギー電子として存在し EEDF には高エネルギーのテールが現れた．西山が行った高エネルギー電子のエネルギー分布関数のモデル計算 ３−９）においても圧力が低い場合は高エネルギーのテールが形成されると報告されており，EEDF は圧力に対して大きく影響されると考えられる．

３．４．２．２

電子の軌道と電子エネルギーの時間変化

次に電子の軌道や電子エネルギー などの個々の電子の挙動ついて 述べる．初めに，放電室内での電子の挙動の傾向を把握する目的で多粒子の軌道を追跡した．５０個の電子の軌道を Fig. 3-10 に示す．
Fig. 3-10(a)から電子は Z 方向に沿って移動していることが分かる．Fig. 3-10(b)は放電室正面から見た場
合を示している．電子は磁石上方の限られた領域で運動している事が分かり，導波管がある 放電室中央にはまったく 移動していない 事も分かる．これは 磁気ミラーの影響であると考えられる．磁気ミラーの効果に関しては直流放電型 でも見られ，ホロカソードから抽出された電子の陽極への流入を防いでいる ．また，中央ほど磁場は弱くなるために 電子は放電室中央にも移動し易く，直流放電型では放電室全体でのプラズマ生成が可能となる．一方，マイクロ波放電型における磁石配置の意味は磁気ミ
ラーの効果よりも，まず ECR の条件を満たす事にある．本計算の場合では必要とされる磁場強度は 0.15
T であり ，この磁場強度は先の直流放電型にくらべれば非常に高い．その分，磁気ミラーの効果も高
まり電子の閉じ込めは向上するが，電子の運動は磁石近傍に制限され，電子が磁石上方に局在化 することになっている．この事がマイクロ波放電型イオンスラスタにおいて密度が均一でない事の要因であると考えられる．
また，Fig. 3-10 (b)に示す電子から５個の電子の軌道を詳しく見ると Fig. 3-11 に示されるように磁石列間を往復しながら左回りに移動している．この運動は Fig. 3-12 に示す∇B×B ドリフトにより説明される．この運動により高エネルギーの電子は壁面に衝突すること無く，電離衝突を行い，ECR により再び加速され，電離衝突を行うというプロセスを繰り返す．この効果は予測されていたが，カップリ
ングコードは∇B×B ドリフトの様子を良く再現している．國中によれば ３−１０），このような電子の軌

44

道は無限軌道と呼ばれ現在のマイクロ波放電型イオンスラスタの性能を達成するには不可欠 な効果である．
次に，５０個の電子のうち２つの電子について軌道とその電子が感じる磁場強度|Be|と電子のエネル
ギーEe の時間変化ついての詳しく見ることにする．Fig. 3-13 に２つの電子の軌道（Z-Y 平面）を示す．
X = 0.064 m の断面における放電室内の磁力線形状と ECR 領域も共に示してある．Fig 3-14, 15 は電子
１，２が感じる磁場|B|と電子のエネルギーEe の時間変化を表している．横軸はマイクロ波の周期で規格化した時間であり，t = 1 はマイクロ波の１周期に相当する．電子１，２ともに電子エネルギーの上
昇は|Be|が 0.15 Tの場合のみ起こることが分かる．電子１に関しては初期位置において|Be |が高く，|Be|が減少しながら t = 20〜25 において |Be| = 0.15 T の領域を通過するが，エネルギーの上昇は見られない．
その後 t = 70 以降に ECR 領域に長時間留まり，Ee が 40 eV程度上昇する．電子２では初期位置での|Be|は低いものの徐々に|Be|が上昇し，t = 150 において ECR 領域に到達する．そして停留時間が電子１の場合と比べてわずかにもかかわらず急激にエネルギーが上昇する．これらの結果は，電子が ECR 領域に到達する事のみが電子のエネルギー上昇の条件ではない事を意味している．最も理想的な ECR は電子とマイクロ波が同じ周波数で回転し，電子とマイクロ波の位相差が-90°になった場合に実現されるが，電子とマイクロ波の位相差が+90°の場合は減速される事になる．電子が ECR 領域に到達した場
合に必ずしも最も理想的な ECR の条件が満たされるとは限らず，エネルギーの上昇幅が停留時間に比
例しない事が起こると考えられる．
Fig. 3-16 は 800 nsec 間に起きた電離衝突の発生位置を示している．Fig. 3-16(a)は Z-Y 平面に対して，Fig. 3-16(b)は X-Y 平面に対する投影図である．Fig. 3-10 (b)で見たように，電子の軌道が磁石上方に制限される為に電離衝突にも同様の傾向が見られる．この結果はマイクロ波放電型イオンスラスタでは放電室中央 でのプラズマ生成は期待できない事を示唆している．この原因はマイクロ波の導入に用いられる 円筒導波管の使用にある ．導波管を使用する事による利点は，アンテナから放射されたマイク
ロ波を低損失で遠方に（ここでは放電室）に伝播させる 事が可能となり アンテナがプラズマに直に接しない事である．しかし，導波管が設置された場所は磁石の配置が不可能 であり，結果として Fig. 3-16で示される ようにプラズマの生成が出来ない領域となる ．従って，放電室全体でプラズマ生成を行っている直流放電型に比べ推力密度が低くなる．

３．５ 第３章のまとめ
本章では MUSES-C 計画で使用されているマイクロ波放電型イオンスラスタについてカップリング
コードを使用して解析し，以下の結果を得た．
１）グリッド部の境界条件を吸収条件と反射条件とした 場合，反射条件の方が電子エネルギーの上昇が高い．また，電子エネルギーの上昇の位置がグリッド側へ移動していた．
２）EEDF の解析結果からは，ECR 領域では非 ECR 領域に比べて高エネルギー電子の存在確率が約45４倍高く，放電室圧力が低い場合には高エネルギーのテール が形成される．
３）電子の軌道を解析した結果から電子は磁石が作る強いミラー磁場に閉じ込められ，電子の移動範
囲は磁石上方に制限される．この結果，電離衝突も磁石上方の ECR 領域でのみ起こり，放電室中
央での電離衝突は起きていない．

46

47

(a)

(b) External magnetic field configuration

Fig. 3-1 (a)Calculation system of microwave ion thruster (MUSES-C )(b)

Table 3-1 Calculation condition

48

(a)

(b)

Fig. 3-2 Electric field intensity at 1.25 nsec(a) Absorption boundary, (b)Reflection boundary

(a)

(b)

Fig. 3-3 Electric field intensity at 3 nsec(a) Absorption boundary, (b)Reflection boundary

49

(a)

(b)

Fig. 3-4 Spatial distribution of electron energy for absorption boundary at 100 nsec(a) color map, (b)contour

(a)

(b)

Fig. 3-5 Spatial distribution of electron energy for reflection boundary at 100 nsec(a) color map, (b)contour

50

Fig. 3-6 Averaged electron energy on X = 0.064 m

Fig. 3-7 Distribution of electron energy at 800 nsec

51

Fig. 3-8 Region used for calculating EEDF

Fig. 3-9 EEDF at 800 nsec

52

53

(a)

(a) Z-Y plane, (b) X-Y plane

Fig. 3-10 Trajectory of 50 electrons

(b)

(a)

(b)

(c)

(d)

Fig. 3-11 Electron trajectory (a) ~ (d)Fig. 3-12 Schematic of∇ B×B drift between magnet lows54Fig. 3-13 Trajectory of electron in discharge chamber,magnetic force line and ECR layer

55

Fig. 3-14 Time variation of electron energy Ee and magnetic field strength |Be|
(Particle 1)

Fig. 3-15 Time variation of electron energy Ee and magnetic field strength |Be|
(Particle 2)

56

57

(a)

during 800 nsec (a) Z-Y plane, (b) X-Y plane

Fig. 3-16 Ionization collision sites

(b)

参考文献
３−１）

平山美晴，中北正人， マイクロ波放電型イオンスラスタにおける電子サイクロトロン共鳴のシミュレーション

３−２）, 日本航空宇宙学会論文集 , 47, 268 (1999)

船木一幸， ECR プラズマ源内のマイクロ波伝播解析, 平成１４年度宇宙輸送シンポジウ

ム (2000)
３−３）

H. Masui, T. Tanoue, H. Nakashima and I. Funaki, ”Numerical Analysis of Microwave Ion Engine
using Electromagnetic Particle Code”, Thin Solid Films (掲載決定 )３−４）

H. Masui, N. Ya mamoto, H. Nakashima and I. Funaki, “Analysis of Plasma Behavior in M icrowave
Discharge Ion Engine”, Proceedings of Asian Joint Conference on Propulsion and Power,
AJCPP2005-22091 (2005)

３−５）

W. J. Harrold, “Calculation of Equipotentials and Flux Lines in Axially Symmetrical Permanent
Magnet Assemblies by Computer”, IEEE Transactions on Magnetics, 8, 23 (1972)

３−６）

I. Funaki, H. Kuninaka and K. Toki, “Plasma Characterization of a 10-cm Diameter Microwave
Discharge Ion Thruster”, J. Propulsion and Power, 20, 718 (2004)

３−７）

岡田文明， マイクロ波工学

３−８）

J. R. Brophy, P. J. Wilbur, “Simple Performance Model for Ring and Line Cusp Ion Thrusters ”,

基礎と応用, 学献社 (1997)

AIAA Journal, 23 1731, (1985)
３−９）

西山和孝，

電子サイクロトロン共鳴型イオンスラスタの放電プラズマ工学研究科博士論文
３−１０）國中均，

東京大学大学院(1997)

無電極マイクロ波イオンスラスタの研究・開発(1998)

58, 日本航空宇宙 学会誌, 46, 174

第４章

４．１

マイクロ波放電型中和器に関する解析

緒言

この章では宇宙科学研究所（現宇宙研究開発機構）によって開発され，MUSES -C 計画で使用されている マイクロ波放電型中和器についての解析結果を述べる．MUSES-C 計画においては，第３章で述べたマイクロ波放電型イオンスラスタと共にマイクロ波放電型中和器が搭載されている．イオンエ
ンジンでは噴出したイオンビームに対して同量の電子を噴出し，衛星本体の電気的中性を保つ中和器とよばれる電子源が必要である．中和器が故障し，衛星からイオンだけが放出されることになると，衛星本体の電位が下がりイオンビームの抽出が不可能となる．直流放電型では中和器にもホロカソードが使用されており，イオンスラスタ本体と同様に電極の損耗や電源構成の複雑化を故障要因として持っている．ヨーロッパで盛んに研究されている高周波放電型イオンスラスタでも中和器にはホロカソードが使用されておりホロカソードが持つ問題を解消できていないのが現状である．
MUSES -C 計画において使用されるマイクロ波放電型中和器は，マイクロ波投入のための L 型アン
テナとミラー磁場形成のための磁気回路から構成され，ECR（ Electron Cyclotron Resonance）放電を利用しプラズマを生成する．電子の抽出は中和器 とイオンビームとの電位差により行われる．マイクロ波放電型中和器の電流値，耐久性は十分に達成されているが，その内部のプラズマ生成機構は未だ完
全には解明されてはいない．特に電子が電離衝突に至るまでのエネルギー上昇のプロセスの把握が重要であると考えられている ．実験 ４−１）では，カットオフ密度以上のプラズマが生成されていることが示唆されているが，中和器のサイズの制約上，実験による検証は非常に困難であり，数値計算が有効であると考えられる．
マイクロ波放電型イオンスラスタに関する数値計算の報告例を第３章で紹介した．一方，マイクロ波放電型中和器に関しての報告は文献４−２）よって行われているが，この研究以降はマイクロ波放電型中和器についての数値解析は積極的には行われてこなかった．これはマイクロ波放電型中和器がECR という共鳴現象とアンテナを利用した放電が組み合わさった複雑な体系であるという 点にある．
また，ECR，アンテナ，ミラー磁場の相互作用を扱った体系はプラズマ工学の観点からも報告例が少ない．
本章では，マイクロ波放電型中和器において，磁場配位とアンテナ位置により電子エネルギーの上
昇と電離衝突発生位置がどう変化するかを，数値計算により明らかにする事を目的とする．本コードを用いてマイクロ波放電型中和器を対象とし，磁場配位とアンテナ位置を変化させた場合のプラズマ生成のメカニズムの解明は，今後のアンテナを利用したマイクロ波放電型プラズマ源の開発 ４−３，４）に有益であると考えられる．本章では，最初に計算体系，計算条件について述べ，次にマイクロ波放59

電型中和器において磁場配位，アンテナ位置を変化させた場合の電子エネルギーの上昇のシミュレー
ション結果について述べる．また，電子エネルギー分布関数（EEDF）と電離衝突発生位置についても合わせて議論する．

４．２ 計算体系
計算体系を Fig. 4-1 に示す．今回は磁場配位とアンテナ位置を変化させ，以下に示す，
Case 1 アンテナ と ECR 領域が接近したモデル
Case 2 放電室内 に ECR 領域が無いモデル
Case 3 アンテナ と ECR 領域を離したモデルの３つの体系での比較を行った．Case 1 は「MUSES-C」に用いられている 中和器をモデルとし，外部磁場は放電室周囲に離散的に設置されている 6 個のサマリウムコバルト磁石により形成される．マイクロ波はバックヨークと L型アンテナによって形成される 同軸ケーブル中を TEM 00 モードで伝送される．Case1 において放電室内に挿入されているアンテナの長さはフロントヨーク出口からアンテナ先端までが 0.008 mである．Case 2 は ECR領域が放電室内に存在しないように磁石を３個にしたモデル，
Case 3 ではアンテナ先端を Case 1のモデルから 0.003 m 離した． それぞれの体系はメッシュサイズ
が X，Y，Z 方向ともに 5 x 10-4 m であり，メッシュ数は X 方向に 60， Y 方向に 60，Z 方向に 44 である．
Fig. 4-2 (a)に Case 1，3 で用いた外部磁場形状を示す．Case 1 では ECR 領域（0.15 T）はバックヨークの前方のアンテナの近傍に存在している．Fig. 4-2 (b)に示す様に Case 2 では 3 個の磁石を用いた
為に非対称の分布であり，ECR 領域は存在せずアンテナ近傍での磁場強度は 0.04 T 程度となった．

４．３ 計算条件
ここで用いた計算条件を Table 4-1 に示す． タイムステップについては FDTD 法が 5.0x10-13 sec であり，PIC 法では 1.0x10-11 sec である．カップリングコード上では FDTD 法を 20 回計算した後に PIC 法
よる計算を１回行う．マイクロ波周波数は，4.25 GHz，マイクロ波電力は 5W である ．アンテナの根元に与えた電場は 42,000V/m である．中性粒子としてはキセノンを採用し，中性粒子数密度は 3.21x1020
m-3（圧力は 1.33 Pa を仮定）とした．電子の初期配置は放電室内に一様分布とし，速度分布としては 1
e V のピークを持つ Maxwell 分布とした．放電室内の初期プラズマ（電子）の密度は 1.0x1017 m-3 とし放電室内で一様の分布とした．計算で取り扱った超粒子数は 2,846,000 個で，超粒子１個は 125,000 個の粒子を代表している．
マイクロ波に対する境界条件については，ヨーク と金属壁で完全反射条件 ，オリフィスで吸収条件とした．アンテナ部はアンテナに平行な電場成分を 0 とすることでモデル化を行った．粒子に対する60

境界条件として 新たにシースを模擬した粒子境界条件を導入する．第２章で本研究では電子のみを取り扱うとしたが ，実際の放電室内には，電子，イオン，中性粒子が存在し，電子はイオンに比べ速度が速く放電室壁面にはイオンより多くの電子が到達する．よって放電室壁面 には電子を跳ね返すようにシースが形成され，電子が放電室内に閉じ込められる．この影響を考慮するため，壁面において高エネルギー電子がシースに入射したとき，シースポテンシャル Vs を超えたエネルギーであれば壁面損失となり，エネルギーが低ければ反射される粒子境界条件を導入した ４−５）．

∗ E ≤ Vs [eV ] ⇒ 反射
∗ E > Vs [eV ] ⇒ 損失
本計算では，シースポテンシャル Vs は実験によって報告された 値 20 e V を用いた ４−６）．
粒子の入射エネルギーを E[eV ]，入射速度を v e ，反射後のエネルギーを E′[eV ]，反射後の速度を v ′e
とすれば，

(

) 2q E(4.1)() 2q E ′(4.2)

v e2 = vex2 + vey2 + v ez2 =v ′e 2 = v′ex2 + v ′ey2 + vez′2 =mmと表わされる．ここで，入射速度の各成分が持つエネルギーの割合は，反射後も保存されるとし ，反射後の速度の各成分はそれぞれ符号が逆転すると仮定すれば ，v′ex2 v ex2=
v′e 2 v e2

v′ey2
v′e 2

=

(4.3)

v ey2

(4.4)

ve2

v′ez2 v ez2=
v′e 2 ve2

(4.5)

となり，反射後の速度の各成分は

v ′ex = −v ex

v ′e
ve

(4.6)

v ′ey = −vey

v′e
ve

(4.7)

v ′ez = −v ez

ve′
ve

(4.8)

と表わされる．ここで，高エネルギーの粒子は多くのエネルギーを失い，低エネルギー粒子は少ない
エネルギーを失うものとすれば，反射後の粒子のエネルギーは，

61


E
E ′ =  1 −  E
 Vs 

[eV]

(4.9)

E′
E
= 1−
E
Vs

(4.10)と仮定される．よって，（4.1）式，（4.2）式より，

v′e
=
ve

であるから，シースによって反射された粒子の速度成分は以下の式で表わされる．

v ′ex = −v ex 1 −E
Vs(4.11)

v ′ey = −vey 1 −E
Vs(4.12)

v ′ez = −v ez 1 −E
Vs(4.13)

４．４ 計算結果と考察
４．４．１

電場強度の空間分布

Fig. 4-3 に電場強度|E|の空間分布を示す．Fig. 4-3(a)は Case 1, 2 での t = 40.0 nsec，Fig. 4-3 (b)は Case
3 での t = 42.5 nsec における電場強度分布である ．両者ともに共通しているのはアンテナ近傍での電場強度が高い事である．本計算ではプラズマを含んでいるにもかかわらず，Case 1 と２では電場強度には違いが見られなかった．Case 1 と２ではアンテナ形状が同じでありマイクロ波の放射は変わらない事に起因しているが，プラズマがマイクロ波に与える影響がアンテナの高い電場強度に比べ低い事も要因の一つとして考えられる． Case 1 では ECR 領域での|E|は 9000 V/m以上であるに対して，Case 3
では ECR 領域での|E|は 1000~3000 V/mとなった．

４．４．２

電子エネルギーの空間分布

Fig. 4-4 に各体系 での t = 200 nsec における電子エネルギー Ee の空間分布を示す．Fig. 4-4 (a)は Case
1 での電子エネルギーの上昇を示している． Ee の上昇はアンテナ近傍と ECR 領域で見られ， ECR
領域内では 3 e V 以上となった．本計算では放電室内の圧力が高いために衝突によってエネルギーが低く抑えられたものと考えられる．Case 2（Fig. 4-4 (b)）の結果でも Ee の上昇は見られるものの最も高いアンテナ近傍の領域でも 2.5 e V と Case 1 より低い．Case 3（Fig. 4-4 (c)）では Ee の上昇はアンテナ近傍に限られ，ECR 領域でのエネルギーの上昇は見られない．この理由として Fig. 4-3(b) から分かる62ように，ECR 領域内でのマイクロ波の電場強度が Case 3 では小さい事が挙げられる．

４．４．３

電子エネルギー分布関数（EEDF）

次に各体系における EEDF（200 nsec 間の平均）について述べる．本計算では Fig. 4-5 に示す領域Ⅰ〜Ⅲの 0.025 m（5 セル）四方の立方体ごとに EEDF を出力した．電離エネルギーEi (=12.13 eV)は以下の図中に示されている．Fig. 4-6(a) に示される Case 1 での EEDF を見ると ECR が存在する領域Ⅰ，アンテナが挿入されている領域Ⅱでは Ei より高いエネルギーを持った電子の存在比率が高いことがわかる．Case 2 (Fig. 4-6 (b))ではアンテナが挿入されている領域Ⅱにおいて最も高い．Fig. 4-6 (c)に示す Case
3 での EEDF では，アンテナが挿入されている領域Ⅲで Ei を超える電子の存在比率は高い．一方,ECRの存在する領域Ⅰでは Case 1 の場合と比べ Ei を越える高エネルギー電子がほとんど存在しない．

４．４．４

電子の軌道とエネルギーの時間変化

次に Case 1 での電子の軌道を Fig. 4-7 に示す．Fig. 4-8(a)は電子エネルギーEe と磁場強度|Be|の時間
履歴，Fig. 4-8(b)は電子エネルギーEe と電場強度|Ef|の時間履歴を示している．横軸は時間を示している，ここでは 100 nsec（4.25GHz のマイクロ波の周期で 425 周期分）について示している．アンテナ近傍の初期位置(X = 0.017 m, Y = 0.019 m, Z = 0.013 m)から電子は磁力線に沿ってフロントヨークの方へ移動し，フロントヨーク近傍で反射され ECR 領域を横切り，バックヨーク近傍で反射して，以降これを繰り返し，100 nsec 間でヨーク間を２往復している．これを Fig. 4-8(a)と照らし合わせると，この電子の場合は ECR 領域を横切る t = 30 nsec において 30 eV程度の急激な Ee の上昇が見られ，その直後に電離衝突を起こし Ee が減少するが，Fig. 4-8 (b)を見ると t = 65 nsec でアンテナ近傍を通過し，強い電場によって再び 10 eV程度の Ee の上昇が見られる．

４．４．５

電離衝突発生位置

Fig. 4-9, 10,11 に t = 200 nsec までに起きた電離衝突の位置を示す．Fig. 4-9 に示す Case 1 の場合では電離衝突の発生は ECR 領域，アンテナ近傍に限らず，アンテナ先端（Z = 0.012 m）以降にも発生している．Case 2（Fig. 4-10）ではアンテナ先端で電離衝突は起こっているが衝突数は少ない．Case 3（Fig.
4-11）ではアンテナ近傍で電離衝突は起きているが ECR 領域では電離衝突は見られない．Case 1，3に共通して Z=0.015 m を境に急激に電離衝突数が減少しているが，これは Fig. 4-7 に示す磁力線の形状から Z=0.015 m 以降では磁力線が Z 軸に対してほぼ垂直になる為に Z 方向への電子の移動が妨げられた結果であると考えられる.

４．４．６

考察

初めにマイクロ波放電型中和器のプラズマ生成の機構について考える．Fig. 4-9 に示した様に Case 1においては電離衝突発生の位置はアンテナ近傍と ECR 領域のみでなくアンテナ前方に拡大している．

63

Fig. 4-7 に示した電子の軌道からも分かるように，フロントヨークとバックヨークの間に存在する電子はミラー磁場の効果によって閉じ込められ，アンテナ近傍の強い電場と ECR によって効率良くエネル
ギーが供給される．その後，磁力線に沿って移動し，その途中で衝突を起こすために電離衝突発生の位置がアンテナ前方にまで広がると考えられる．中和器でのプラズマの維持はミラー磁場に閉じ込められた電子の ECR およびアンテナ近傍の電場によるエネルギーの上昇，衝突によるエネルギー損失，
再びエネルギーの上昇というサイクルの繰り返しによると考えられる．
次に ECR 領域の有無について考えると Fig. 4-10 から ECR 領域が存在しない場合でもアンテナ近傍の強い電場のみで電離は可能な事が分かる．Fig. 4-8 からアンテナ近傍の強い電場でも 10 e V 程度のエ
ネルギー上昇は可能であり，この事が Case 2 で電離衝突が発生したことを裏付けている．この点は第３章で述べたマイクロ波放電型イオンスラスタにおいては電子が ECR に達しない限りエネルギーが上昇しない点と異なる．しかし，アンテナ近傍の強い電場での Ee の上昇は ECR 領域での Ee の上昇に比べ小さい．この理由としては電子がアンテナ近傍に存在した場合に常に電場が磁力線に直交し，電子を加速する方向に強いとは限らない為である ．磁場を用いないマイクロ波放電型プラズマ源は文献４−７）で紹介されているが，Ee の上昇と電離衝突の頻度の点から比較するとプラズマの生成の効率
は ECR 領域が存在する方が良いと言える．
アンテナ位置についてはマイクロ波の電場は Case 1，3 ともに放電室内の ECR 領域に到達している
事が Fig. 4-3 から分かる．しかし ，電場強度が，ECR 領域内で Case 1 においては 9000 V/m以上，Case
3 においては 3000 V/m以下と大きく異なる事で電子のエネルギー上昇および電離衝突の発生でも大きく異なる．この原因を考えてみる．アンテナ工学的にはマイクロ波の波長をλとすれば，アンテナ長がλ/4 の場合に良く放射される．本計算で使用した 4.25GHz のマイクロ波のλ/4 は約 0.017 m であり，
ここで採用したアンテナ長（放電室入り口からアンテナ先端まで）はλ/8（= 0.008 m）に近い．従って，λ/8 アンテナはマイクロ波を遠方まで到達させる能力ではλ/4 の場合に比べて劣り，アンテナとECR 領域の距離が離れた場合の電子エネルギー上昇は期待できない．しかし，アンテナが放射に最適な長さ，形状を持たない場合でも，ECR 領域とアンテナとの距離を接近させる事でアンテナ近傍の強い電場によって電子は加速され，電子エネルギーは上昇する．また,これを実用化の点から考えると中和器は小型軽量であることが望ましく，結果としてλ/8 アンテナの使用は中和器の小型化に貢献している．

４．５ 第４章のまとめ
本章では MUSES-C 計画で使用されているマイクロ波放電型中和器についてカップリングコードを使用して解析し，以下の結果を得た．
１）マイクロ波放電型中和器 のプラズマ 生成については，ヨーク間に存在する電子のエネルギーがECR とアンテナ近傍の強い電場によって上昇し，その電子による電離衝突によってプラズマが生

64

成，維持されている．
２）MUSES -C 計画で使用されている中和器の L 型アンテナはλ/8 に近い長さである．よって，マイク
ロ波の放射の効率は低いと考えられ，ECR 領域とアンテナの距離が離れた場合は電子のエネルギ
ー上昇が望めない．しかし，アンテナと ECR 領域を接近させる事によりプラズマ生成を可能としており，結果としてマイクロ波放電型中和器の小型化，軽量化に寄与している．

65

(a)

(b)

(c)

Fig. 4-1 Calculation geometry(a) Case 1, (b) Case 2, (c) Case 366

(a) Cases 1 and 3 (6 magnets)(b) Case 2 (3 magnets)

Fig. 4-2 (a) Magnetic field configuration(a) Case 1 and Case 3, (b) Case 2

67

Table 4-1 Calculation condition

68

(a) Cases 1 and 2 (t = 40.0 nsec)(b) Case 3 (t = 42.5 nsec)

Fig. 4-3

Distribution of electric field strength(a) Cases 1 and 2 at 40.0 nsec and (b) Case 3 at 42.5 nsec

69

(a) Case 1(b) Case 2(c) Case 3

Fig. 4-4 Distribution of electron energy Ee (t = 200 nsec)(a) Case 1, (b) Case 2 and (c) Case 3

70

Fig. 4-5 Definition of calculation region for EEDF

71

(a) Case 1(b) Case 2(c) Case 3

Fig. 4-6 EEDF averaged over 200 nsec(a) Case 1, (b) Case 2 and (c) Case 372

Fig. 4-7 Trajectory of electron and line of magnetic force at X = 0.017 m

73

(a)

(b)

Fig. 4-8 Time variations of the electron energy Ee , the magnetic field strength |Be|and the electric field strength |Ef|(a) Ee and |Be|, (b) Ee and |E f|

74

75

(a)

(a) Z-Y plane, (b) X-Y plane

Case 1during 200 nsec

Fig. 4-9 Ionization collision sites for

(b)

76

(a)

(a) Z-Y plane, (b) X-Y plane

Case 2 during 200 nsec

Fig. 4-10 Ionization collision sites for

(b)

s

77

(a)

(a) Z-Y plane, (b) X-Y plane

Case 3 during 200 nsec

Fig. 4-11 Ionization collision sites for

(b)

参考文献
４−１）

I. Funaki and H. Kuninaka, ‘’Overdense Plasma Production in a Low-power Microwave Discharge
Electron Source”, Jpn. J. Appl. Phys. 40, 2495 (2001)

４−２）

吉永睦生，佐鳥新，船木一幸，小野寺範義，石川芳男，
マイクロ波中和器のプラズマ解析

４−３）, 平成１０年度宇宙輸送 シンポジウム (1999)

T. Miyamoto, K. Mii, T. Nishijima, H. Ijiri, Y. Takao and H. Nakashima, “Development of a New
Microwave Discharge Type Ion Engine”, Vacuum, 73, 391 (2004)

４−４）

Y. Takao, H. Masui, T. Miyamoto, H. Kataharada, H. Ijiri and H. Nakashima, “Development of
Small-Scale Microwave Discharge Ion Thruster”, Vacuum, 73, 449 (2004)

４−５）

J. C. Adam, A. Heron and G. Label, ”Study of Stationary Plasma Thrusters Using Two-dimensional
Fully Kinetic Simulations”, Phys. Plasmas, 11, 295 (2004)

４−６）

小野寺範義，竹ヶ原春樹，西山和孝，船木一幸，國中均，
マイクロ波放電型中和器の電子放出機構 ，日本航空宇宙学会論文集, 49, 27 (2001)４−７） C. Fusellier, L. Wartski, J. Aubert, C. Schwebel, Ph. Coste and A, Chabrier, “A Microwave Plasma
Cathode Electron Gun for Ion Beam Neutralization”, Rev. Sci. Instrum., 69, 1153 (1998)

78

第５章

アンテナを利用したマイクロ波放電型
プラズマ源に関する解析

５．１

緒言

この章ではアンテナを利用したマイクロ波放電型プラズマ源に対する解析について述べる．第３章
で述べたように，MUSES-C 計画において用いられているマイクロ波放電型イオン スラスタではプラズマの生成が不均一である．この原因は導波管 の使用によって放電室中央にプラズマの生成が出来ない領域が生ずる事に起因していた．そこで，導波管を用いない方式がいくつか考案されている．その一つとしてマイクロ波放電型中和器の特性を利用した方式がある．第４章で述べたようにマイクロ波
放電型中和器は高密度プラズマが生成されており，この特性を利用する事で導波管を用いずに，放電
室内に均一プラズマを生成できると考えられる．具体的には複数本のアンテナを放電室内に挿入する事によって均一プラズマの生成を行う方式である．
これは マルチモノポール方式 ５−１）と呼ばれ Fig. 5-1に示すような構造を持つ．複数本のアンテナによりプラズマを生成する事に加えて，複数本のアンテ
ナへのマイクロ波の分配にはコネクタ等の結線によるものではなく，アンテナからアンテナへとマイ
クロ波の中空での伝播を利用して分配する方式を採用している．これにより，マルチモノポールアン
テナシステム には従来のマイクロ波放電型イオンスラスタに必要であった，DC ブロッ クと呼ばれるマイクロ波電源への直流成分の印加を防ぐ機能を内蔵している．このマルチモノポール方式を利用して，マイクロ波放電型イオンスラスタの性能向上，特に推力密度の向上が期待されている．
また，マイクロ波放電型中和器を高密度イオン 源と見なし，小型マイクロ波放電型イオンスラスタ
としての応用も研究されている．この背景には，現在の宇宙開発における衛星の大型化，研究開発費の高騰化という問題がある．これらの問題が積極的な宇宙開発の妨げとなっている．また，そのほとんどが国家プロジェクトの一環として運営されてきた事も問題であり，民間企業や大学にとって宇宙開発は参画し難い分野であった．そこで注目されたのが 10~100 kg の小型衛星の開発 ５−２）である．衛星を小型化 することで民間企業や大学でも参入できるような 宇宙開発プロジェクト ５−３，４）の規模になる．また，衛星を小型化し，複数個の衛星を打ち上げるという事も考えられており衛星の打ち上げや運用の中で起こりうる多くのリスクを分散化 するという点からも注目されている．こうした背景のもと小型マイクロ波放電型イオンスラスタは小型衛星用推進機として期待されており，現在その性能は徐々に実用レベルに近づきつつある．小型マイクロ波放電型イオンエンジンについての実験や性能については付録にて詳しく述べる．
これらの新型，改良型イオンスラスタに共通している事は，アンテナを放電室に直接挿入することによるプラズマ生成方式を採用しているという点である．このアンテナを利用する方式の場合，アン79

テナを放電室内に複数本に配置する事で均一化 が望める事以外にカットオフの制限を受けないという利点がある．第４章でマイクロ波放電型中和器の実験結果からカットオフ密度を超えるプラズマの生成が確認されている事を述べた．マイクロ波放電型中和器に限らず，小型の高密度プラズマ源は多く開発されているが ５−５），これらの小型で永久磁石を用いたプラズマ源の場合，磁場が急激に変化する．
その為，第２章で述べたような外部磁場 B とマイクロ波の波数 k の方向の厳密な区別が難しく，B と
k の方向を分類する第２章のような方法での解析は難しい．従って高密度である理由は，分散関係の観点からは明らかにされていないが，マイクロ波がカットオフ密度以上のプラズマ中にでもスキンデ
プス長は進入する事が可能である為とされている．プラズマでのスキンデプス長 d は次式であらわされる ５−６）．

()

1

δ = c / ωp − ω 2 2
2

(5.1)

ここで c は光速，? p はプラズマ角周波数，? はマイクロ波角周波数である．スキンデプスの効果は，金属に電磁波を放射した場合によく見られ，電磁波のほとんどは金属表面で反射されるが，スキンデ
プス長の分だけは金属中に進入可能である．これをプラズマに当てはめると，プラズマはマイクロ波
周波数がカットオフ周波数以上ならばマイクロ波を吸収，透過するように振る舞い，カットオフ周波数以下ならば金属のようにマイクロ波を反射するようになる．しかし，金属と同様にマイクロ波のほ
とんどは反射されるが，スキンデプス長の分はプラズマ中に伝播する事になる．大きなプラズマ源である場合はこの効果の影響は少ないが，小型プラズマ源である場合はスキンデプスの効果が大きくなる．(5.1)式を用いて中和器の場合の d を計算すると，? = 4.25 x 109 x 2p rad/sec，プラズマ密度 n = 1 x1018 m-3 ，? p =5.6 x 1010 rad/sec の条件では 0.005 mとなる．放電室の半径が 0.009 mの中和器では d の影響が大きい事が分かる．この利点を利用すれば，マイクロ波周波数に依存せずにイオンスラスタに必要な高密度プラズマを得る事が可能となり，マイクロ波電源の電源効率が良い周波数帯でのプラズマ生成が可能となる．4.25 GHz のマイクロ 波電源では電源効率が約 60 %と低く，現在情報通信の分野において研究が進んでいる効率の良い周波数帯の電源が利用可能となれば，イオンエンジンの大幅な性能向上が可能となる．
前章でのマイクロ波放電型中和器の解析においては，放電室とヨークが一体となっており磁場形状の任意の変更は困難であり，またアンテナ形状と磁場形状との関係の変更も困難であった．今後のアンテナを利用した新型および改良型イオンスラスタの開発を考えると，より詳細な解析，そして設計指針を得る必要があると考えられる．特に前述のマルチモノポールシステムを利用した 10 cm級イオンスラスタについては，MUSES-C 計画のマイクロ波放電型イオンスラスタと同等の性能に到達できていない．そこで本章ではカップリングコード を利用して，アンテナを利用したマイクロ波放電型プ
ラズマ源についての解析を行った．

80

Fig. 5-1 Schematic of multi-monopole antenna system81

５．２ 一様磁場中での挙動解析
始めに，一様磁場中においた電子の挙動解析を行った．ここでは前章までのように多粒子の計算ではなく ，３個の粒子だけを 追跡するものとし，アンテナと電子の距離，アンテナ挿入方向と磁場の方向の関係（平行，垂直），磁場強度を変化させ計算を行った．計算体系を Fig. 5-2 に示す．アンテナの形状はモノポールで，長さは 0.015mである．このアンテナはマイクロ波周波数 2.45GHz の波長λに対して 1/8 の長さを持つ．電子の初期配置はそれぞれ，電子１がアンテナ先端，電子２はアンテナ先端から Z 方向に 0.005 m 離れた位置，電子３は電子１から Y 方向に 0.025 m 離れた位置である．計算条件
としてはマイクロ波周波数が 2.45 GHz，アンテナの根元に与えた電場は 10000 V/mである．タイムス
テップは 5 x 10-13 sec である．最初にアンテナから放射されたマイクロ波の電場について約１周期間の強度分布を Fig. 5- 3 (a) ~ (g)に示す．いずれの時刻においても Z = 0.025 m 以降の特にアンテナ正面で電場強度が低く，かつ急激に下がっている．また，アンテナ先端の上下方向（±Y 方向）には強い電場が保たれている．これはモノポールアンテナ の指向性 に起因しており，アンテナ正面へ電磁波 の放射が低い事を示している．Fig. 5-4(a) ~ (g)には Fig 5-3 と同時刻の電場ベクトルを示している．アンテナ近傍の電場ベクトルを確認し易くする為にアンテナ上の電場ベクトルは除いて表示してある．この図からアンテナ近傍では Ez 成分が少ない事が分かる．

５．２．１

アンテナと磁場が平行な場合

まず，アンテナと磁場が平行な場合の計算について述べる．この計算において磁場は+Z 方向に一様で,その強度は+ 0.0875 T（2.45 GHz に対する ECR を満たす磁場強度）である．Fig. 5-5 は電子１〜３の軌道を示している．電子１ではラーマ半径が大きくなっている事が分かるが ，電子２，３ではラーマ
半径がほとんど変化していない．これらの電子のエネルギーEe と電子が感じる電場強度|Ef|の時間変化
を Fig. 5-6 に示す．Fig. 5-6 において横軸はマイクロ波周波数で規格化した時間である．Fig. 5- 6 (a)に示す電子１に関しては常に電場強度が高く，時間と共に電子のエネルギーが急激に増加する．電子２(Fig. 5- 6 (b))では電子１の場合に比べ電場強度が低く，エネルギーの上昇が見られない．Fig. 5- 6 (c)に示す電子３では僅かであるが電子のエネルギーの上昇が確認できる．
電子２では t = 15 において電子３と同じ程度の電場強度があるにもかかわらずエネルギー の上昇は見られない．このことは，電子２と電子３に位置における電場方向成分が Fig. 5-4 で示したように異なる事からきており，電子２の位置では，
Z 方向磁場に垂直で ECR による電子の加速に寄与する電場成分 Ey が低い．
これらの結果から電子エネルギーの上昇はアンテナと ECR 領域との距離に強く依存する事が分かる．
これは実験でも報告されており ５−７），アンテナを設置する場合にはアンテナと ECR 領域を一致させる必要がある．

82

５．２．２アンテナと磁場が垂直な場合

次にアンテナと磁場が垂直な場合の計算について述べる．この計算において磁場は+Y 方向に一様で,その強度は+ 0.0875 T である．Fig. 5-7 は電子１〜３の軌道を示している．３つの粒子の中では電子
１（アンテナ先端）において最もラーマ半径が大きい．しかし，Fig. 5-8 (a)に示す Ee の時間変化を見ると，アンテナと磁場が平行な場合に比べ，電子のエネルギーの上昇は低い．これは ，ECR において電子を加速する電場成分の減少により説明される．いま，磁場は+Y 方向に存在するので ECR において電子を加速できる電場は Ex，Ez である．Ey の電場は磁場と平行なため ECR の加速には寄与しない．
Fig. 5-4 から分かるように，
アンテナ近傍では電子を加速する電場成分の一つである Ez が存在しない．
電子は Ex のみで加速されるために，磁場とアンテナが平行な場合と比べてエネルギーの上昇が低い．
電子２（Fig. 5-8 (b)），３（Fig. 5-8 (c)）について Ee と|Ef|の時間変化を見ると，電子２に対してアンテナと磁場が平行な場合（Fig. 5-6 (b)）と異なりエネルギーの上昇が高い事が分かる．この事から，ア
ンテナ正面では Y 軸に対して垂直な電場成分が低いけれども存在しており，この成分が ECR により電子を加速した為と考えられる．

５．２．３

マイクロ波電場強度依存性

５．２．１の計算ではアンテナからの距離に電子のエネルギーが強く依存する事が確認された．ここでの計算ではマイクロ波電場強度を上昇させた場合のアンテナから離れた位置（電子２，３）での電子エネルギー上昇について調べる．磁場の向きは５．２．１と同じく+Z 方向，磁場強度は 0.0875 Tである．アンテナに供給する電場 Eant を 20000 V/m，40000 V/mと増加させた場合の計算を行った．Eant= 20000 V/mの場合の計算結果を Fig. 5-9 に示す．Fig. 5-9 (a) に示す電子２のエネルギーの時間変化
からは，電子のエネルギーに変動があるものの 10000 V/mの条件（Fig. 5-6 (b)）と同様に低く，1 eV程度である．電子３（Fig. 5-9 (b)）ではエネルギーが 3eVまで上昇する．40000 V/mの場合では，電子２のエネルギーの時間変化は Fig. 5-10 (a)に示すように 20000 V/mの場合と同様にエネルギーの変動が見られる．その幅は 20000 V/mよりも大きいが，結果として電子２のエネルギーはほとんど上昇せず 1 eV程度にとどまっている ．電子３（Fig. 5-10 (b)）の場合はエネルギーの上昇が大きく，最大で10 eV ほどになる．この結果はアンテナに与える電場が高ければ電子３の位置でもキセノンの電離エ
ネルギー以上のエネルギーを持つ電子を生成可能であることを示している．しかし，マイクロ波の電力を上昇させ，アンテナからの距離が離れた電子３のような位置でプラズマ生成を行う場合は，カットオフ密度の問題がある．これは前述したように，マイクロ波放電型中和器でカットオフ密度を超えるプラズマが生成可能な理由は，マイクロ波がプラズマ中へスキンデプスの長さだけ進入できる事によるものである．アンテナと ECR の距離をスキンデプス長以上離した場合は，カットオフ密度に達したプラズマでマイクロ波は反射され，遠くに到達しないと考えられる．従って，高密度プラズマを得るためには，アンテナを ECR 領域に接近させる方式の方が効果的であると考えられる．

83

５．２．４

磁場強度依存性

次に,磁場強度依存性について結果を述べる．一般にマイクロ波放電の場合は ECR が起きる条件として，マイクロ波周波数とサイクロトロン周波数の一致が重要であり，これにより効率の良いプラズ
マ生成が可能となる．しかし，近年，放電室半径が数 mm（短ギャップ長）となるようなマイクロプ
ラズマにおいては 2nd ECRと呼ばれる通常の ECRに必要な磁場の 1/2の強さの方が高密度プラズマが得られると報告されている ５−８）．ECR 条件の場合，電子のエネルギーの上昇が急激に起こるために，
ラーマ半径も急激に大きくなる．放電室半径が小さい場合は，放電室壁への電子の損失が大きくなり，
結果として ECR 条件では高密度プラズマが得られない．2nd ECR 条件の場合，電子エネルギーの上昇
は ECR 条件に比べれば低いものの，それに伴ってラーマ半径の増加も抑えられ，電子の閉じ込めが向上する事で放電室壁への電子の損失が減少し，高密度プラズマが生成可能となる． ここでの計算では磁場の方向を５．２．１と同じく Z 方向のみとし，磁場強度 0.04375 T(0.0875 T x 1/2)と 0.175 T(0.0875
T x 2)について 計算を行った．初期位置は，いずれの条件においても Fig 5-1 の電子１の地点（アンテナ先端）である．Fig. 5-11 は ECR 条件にある 場合の電子の軌道図（X-Y 平面）を示している．この場
合は電子の軌道は徐々に大きくなり，ほぼ真円の軌道を描く．次に，Fig 5-12 に磁場強度 0.04375 T（2nd
ECR）の場合の計算結果を示す．Fig 5-12(a)に示す軌道図では，ECR 条件の場合（Fig 5-11）と比べ電子の軌道は楕円を描く．また，ラーマ半径の最大値で比較すると ECR 条件の場合よりも大きい事が分かる．これは，この電子が-Z 方向（アンテナの根元），つまり電場が高い方向へ動いている為である．
エネルギーの時間変化（Fig 5-12(b)）からは ECR 条件（Fig 5-11(b)）の場合より緩やかにエネルギーが上昇している事が分かり，2nd ECR でも高いエネルギーが得られる事が確かめられた．最後に Fig
5-13 に示す磁場強度 0.175T の場合では，ラーマ半径の増加は見られないが，ラーマ半径が小さい事によって，常にアンテナ近傍に存在している為に電子が感じる電場強度は先の２条件の場合より高く，
10 eV程度電子エネルギーの上昇が見られる．しかし，マイクロ波からエネルギーを得る効率が ECRの条件の場合に比べ低い為にエネルギーの上昇は低い．

84

Fig. 5-2 Calculation geometry for antenna system85

(a) 6 nsec(b) 6 .05 nsec(c) 6 .10 nsec(d) 6 .15 nsec(e) 6 .20 nsec(f) 6 .25 nsec(g) 6 .30 nsec

Fig 5-3 Time evolution of electric field strength(a) 6.0 nsec ~ (g) 6.30 nsec

86

(a) 6 nsec(b) 6.05 nsec(c) 6.10 nsec(d) 6.15 nsec(e) 6.20 nsec(f) 6.25 nsec(g) 6.30 nsec

Fig 5-4 Time evolution of electric field vector(a) 6.0 nsec ~ (g) 6.30 nsec

87

Fig. 5-5 Trajectory of electron at Bz = 0.0875 T(a) Particle 1(b) Particle 2(c) Particle 3

Fig. 5-6 Time variation of electron energy Ee and electric field strength |Ef| at Bz = 0.0875 T(a) particle 1, (b) particle 2, (c) particle 388

Fig. 5-7 Trajectory of electron at By = 0.0875 T(a) Particle 1(b) Particle 2(c) Particle 3

Fig. 5-8 Time variation of electron energy Ee and electric field strength |Ef| at By = 0.0875 T(a) particle 1, (b) particle 2, (c) particle 389

(a) Particle 2(b) Particle 3

Fig. 5-9 Time variation of electron energy and electric field strength at Bz = 0.0875 T and Eant = 20000 V/m(a) particle 2, (b) particle 3(a) Particle 2(b) Particle 3

Fig. 5-10 Time variation of electron energy and electric field strength at Bz = 0.0875 T and Eant = 40000 V/m(a) particle 2, (b) particle 3

90

(a)

(b)

Fig. 5-11 Trajectory of electron 1 at Bz = 0.0875 T(a) trajectory and (b) time variation of electron energy and electric field strength

(a)

(b)

Fig. 5-12 Trajectory of electron 1 at Bz = 0.04375 T(a) trajectory and (b) time variation of electron energy and electric field strength

(a)

(b)

Fig. 5-13 Trajectory of electron 1 at Bz = 0.175 T(a) trajectory and (b) time variation of electron energy and electric field strength91

５．３ 永久磁石の磁場中での解析
５．２では少数の電子の挙動に着目して計算を行い，アンテナを利用したマイクロ波放電 におけるいくつかの設計指針を得た．それは１）アンテナと ECR を接近させる，２）アンテナに対して磁場は平行に生成する，という事である．磁場強度依存性の結果から 2nd ECR の効果も利用できると考えられるが，本計算が対象とする放電室半径 2 cm ~ 10cmのイオンエンジンでは ECR の条件の方が効果的であると考えられ，上述の２条件が設計指針となる．５．２では一様磁場という理想化された 条件であったので，ここでは永久磁石が作る磁場中 での計算を行う．また，多粒子を計算する事で電子エネルギーの空間分布，EEDF，電離衝突位置を算出し，設計指針となるデータの取得を行う．先の２条件を満たす体系として Fig. 5- 14 で示す磁石とアンテナの配置（Case 1）が考えられる．マイクロ波周
波数は 2.45 GHz，アンテナはモノポールで，アンテナ長は５．２の計算と同じく 0.015 m（?/8）とした．磁石は Sm-Co マグネットを考え，厚さが 0.01 mのヨークを使用する．アンテナは磁場に対して平行に挿入されており条件１を満たす．磁場形状と磁場強度は Fig. 5-15 のように 形成されており，ECR
領域は Y = 0.02 mの位置に存在している．条件２を満たすために，アンテナは ECR 領域と接するように Y = 0.02 m の位置に設置した．ここでは Fig. 5-16 に示す体系（Case 2）の結果との比較を行った．
Case2 は Case 1 の磁石とヨークを 90 度左に回転したモデルとなる．Fig. 5-17 に磁場形状と磁場強度を示す．アンテナは磁石とヨーク間から挿入されている．アンテナ長は 0.02 mであり，Fig. 5-16 からアンテナ先端が ECR 領域に接している事が分かる．計算条件を Table 5-1 に示す．両体系ともプラズマ領
域は 0.03 m x 0.03 m x 0.03 m であり，計算で取り扱った超粒子数 は 4,000,000 個であり，超粒子 1 個は125,000 個の実粒子を代表している． マイクロ波に対する境界条件はプラズマ領域の端で吸収境界条
件，粒子に対しては第４章で述べたシースを模擬する粒子境界条件 を適用した．

５．３．１

電子エネルギーの空間分布

Fig. 5-18 に Case 1 場合の 100 nsec における電子エネルギー Ee の空間分布を示す．Fig. 5-18 (a)は X =0.015 m の断面，Fig. 5-18 (b)は Z = 0.02 m の断面における結果を示している．Fig. 5-18 (a)から，アンテ
ナ先端で電子エネルギーが高い事が分かる． ECR に沿った領域で Ee の上昇が見られ，Z = 0.015mの±Y 方向に Ee の高い領域が広がっている．また，アンテナ先端より先の Z = 0.02 m 以降の ECR 領
域では Ee の上昇は見られない．アンテナ正面でのエネルギー上昇が見られない原因は，５．２で述べたようにモノポールアンテナのアンテナ正面へのマイクロ波の放射が低い特性によるものである．Fig.
5-19 は Case 2 の Ee についての結果を示している．Case 2 の場合はアンテナ先端で Ee の上昇が高く，ECR 領域に沿った領域での Ee の上昇は低い． ECR 領域に沿った領域は５．２．２で述べたような磁場とアンテナが垂直の領域であり，かつ，アンテナから離れている為に電子の加速が行われない．ア
ンテナから離れた ECR 領域に存在しているエネルギーの高い電子は，アンテナ先端で加速され，磁力線に沿って移動してきた電子であると考えられる．

92

５．３．２電子エネルギー分布関数（EEDF）

次に EEDF の解析結果について述べる．Case 1 では Fig. 5-20 に示す領域１，２，３，４において EEDFを出力した．なお，出力した EEDF は 100 nsec 間で平均化されている．また，図中の右上に示される
数値はそれぞれの領域に含まれる 12 eV 以上のエネルギーを持つ電子の割合を示している．Fig. 5-21に Case 1 の EEDF を示す．アンテナに近い領域２において高エネルギーの電子の存在確率が高く，その割合は 12.3 %である．領域１，３においても高エネルギーの電子の存在確率は高く，５．２．１の電子エネルギーの空間分布（Fig. 5-18 (b)）と同様の傾向が見られる．５．１の結果からアンテナと ECR領域が接するような位置でのみ 電子エネルギーが上昇する事が分かっており ，領域１，３に存在する
高エネルギーの電子は領域２から移動したものと考えられる．領域４においても 高エネルギーの電子が存在しているが，領域１，２，３に比べると高エネルギーの電子の存在確率は低い．
次に Fig. 5-22 に Case 2 における EEDF の出力領域を示す．Fig. 5-23 に Case 2 における EEDF の出力結果を示す．領域１，２にはほぼ同じ割合で高エネルギーの電子が存在している．領域３では高エネ
ルギーの電子は存在しておらず，初期に与えた Maxwell 分布のままであった．これは５．１で述べたモノポールアンテナ正面での電場強度が低い事に起因している．また，領域１，２と３の EEDF の形状を比較すると，領域１，２においては高エネルギー領域においてテールが形成され，領域３の Maxwell
分布からのずれが大きい事が分かる．
これらの解析結果から，Case 1 と Case 2 のアンテナ近傍での EEDF を比較すると，Case 1 の方が高エネルギー電子の存在確率は高い．

５．３．３ 電離衝突発生位置
Fig. 5-24 は Case 1 について ，Fig. 5-25 は Case 2 について 100 nsec までに起きた電離衝突の発生位置を示している． Case 1,2 共にアンテナ近傍で電離衝突が発生しており，Fig. 5-18，19 に示したような
電子エネルギーの高い領域に電離衝突の発生が集中している．
電離衝突の発生回数は Case 1 では 1263，
Case 2 では 319 である．よって電離衝突の発生が多い Case 1 の方が，プラズマの生成には適しているといえる．

５．４ 第５章のまとめ
アンテナを用いたマイクロ波放電型プラズマ源の解析を行い，以下の結論を得た．
１） ECR 条件を満足する一様磁場（0.0875 T）中のモノポールアンテナの解析から，電子のエネルギ
ーはアンテナの近傍でのみ上昇し，アンテナと磁場の向きが平行である場合が電子のエネルギーの上昇が大きい事が分かった ．アンテナに与える電場を大きくする事によって，離れた距離にある
ECR 領域でも電子のエネルギーが上昇する．磁場強度を変更した結果からは 2nd ECR の効果も確認でき，この効果はより小型のプラズマ源への応用が可能と考えられる．

93２） 永久磁石での解析では，アンテナが磁場に対して平行に挿入されている Case 1 と磁場に対して垂直に挿入されている Case 2 において計算を行い，比較を行った．電子エネルギーの空間分布とEEDF から，Case 1 の方がキセノンの電離電圧を超えるエネルギーをもった電子が多く存在していることが分かった．

94

Table 5-1

Calculation condition95

Fig. 5-14 Calculation geometry for Case 1

(a)

(b)

Fig. 5-15 Magnetic field configuration and ECR layer for Case 1(a) Z-Y plane at X = 0.015 m and (b) X-Y plane at Z = 0.015 m

96

Fig. 5-16 Calculation geometry for Case 2

(a)

(b)

Fig. 5-17 Magnetic field configuration and for Case 2(a) Z-Y plane at X = 0.015 m and (b) X-Y plane at Z = 0.03 m

97

(a)

(b)

Fig. 5-18 Spatial distribution of electron energy for Case 1 at 100 nsec(a) X = 0.015 m, (b) Z = 0.015 m

98

(a)

(b)

Fig. 5-19 Spatial distribution of electron energy for Case 2 at 100 nsec(a) X = 0.015 m, (b) Z = 0.03 m

99

Fig. 5-20 Calculation region of EEDF for Case 1

Fig. 5-21 EEDF for Case 1100

Fig. 5-22 Calculation region of EEDF for Case 2

Fig. 5-23 EEDF for Case 2

101

Fig. 5-24 Ionaization site for Case 1 during 100 nsec

Fig. 5-25 Ionaization site for Case 2 during 100 nsec

102

参考文献
５−１） T. Miyamoto, K. Mii, T. Nishijima, H. Ijiri, Y. Takao and H. Nakashima, “Development of a New
Microwave Discharge Type Ion Engine”, Vacuum, 73, 391 (2004)５−２） 竹ヶ原春樹 ， 電気推進を中心とした 衛星・探査機の推進システムの現状と将来展望, 日

本航空宇宙学会誌, 49, 138 (2001)５−３） N. Kumagai, K. Kawahara, K. Sato, K. Tamura, H. Takegahara, M. Igarashi, M. Sugiki and H.
Hashimoto, “Research and Development Status of Low Power Pulsed Plasma Thruster System forMicro-Lab Sat II”, IEPC-03-0202 (2003)５−４） S-W. Kim, S. Satori, Y. Itoh, H. Okamoto, T. M. Sugiki, Y. Aoki, A. Nagata and T. Kizaki“Development of Microwave Engine”, IEPC-03-251 (2003)５−５） K. Miyake, T. Mikami, H. Tahara, K. Ogata, Y. Matsubara, S. Nogawa and J. Ishikawa,“Improvement of the Gas Efficiency in a Microwave Plasma Cathode Type Oxygen Ion Source”, Rev.
Sci. Instrum., 65, 1313 (1994)５−６） J. Ishikawa, Y. Takeiri and T. Takagi, “Axial Magnetic Field Extraction-type Microwave Ion Sourcewith Permanent Magnet”, Rev. Sci. Instrum., 55 449 (1984)５−７） N. Yamamoto, H. Kataharada, T. Chikaoka, H. Masui, H. Nakashima and Y. Takao“Effect of Antenna Configuration on Thrust Performance in a Miniature Microwave Discharge Ion
Engine”, IEPC-2005-036 (2005)５−８） E. Morisaki and H. Fujiyama, “Development of Transported Coaxial ECR plasma Source for Inner
Coating of Metallic Tubes”, Surface Coating and Technology, 98, 834 (1998)

103

第６章

結論

本論文 では，マイクロ波放電型宇宙推進用プラズマ源における電子とマイクロ波の相互作用を明らかにする事を目的とし，PIC 法と FDTD 法を組み合わせたカップリングコードにより数値解析を行った．具体的には，１）10 cm級マイクロ波放電型イオンスラスタ（MUSES -C model），２）マイクロ波
放電型中和器，３）アンテナを用いたマイクロ波放電型プラズマ源を対象として 計算を行い，以下の結論を得た．

１）10cm 級マイクロ波放電型イオンスラスタ（ MUSES -C model）に関する解析
MUSES-C 計画で使用されている 10 cm級マイクロ 波放電型イオンスラスタについての数値解析を行った．
はじめに，グリッド部の境界条件を吸収条件と反射条件とした場合での比較を行い，反射条件の方が電子エネルギーの上昇が高く，また，電子エネルギー の上昇の位置がグリッド側へ移動する事が分かった．反射境界条件における EEDF の解析結果からは，ECR 領域では非 ECR 領域に比べて高エネ
ルギー 電子の存在確率が約４倍高く，放電室圧力が低い場合には高エネルギーのテールが形成される事が分かった．電子の軌道を解析した結果から，電子は磁石が作る強いミラー磁場に閉じ込められ，電子の移動範囲は磁石上方に制限される．この結果，電離衝突も磁石上方の ECR 領域でのみ起こり，放電室中央 での電離衝突は起きていない事が分かった．この不均一 なプラズマ生成が，直流放電型イオンスラスタに比べ推力密度が低い原因であると考えられる．

２）マイクロ波放電型中和器に関する解析
MUSES -C 計画で使用されているマイクロ波放電型中和器についてカップリングコード を使用して解析を行った．ここでは磁場形状やアンテナ形状を変化させ，電子エネルギーの上昇と電離衝突発生位置を比較する事により，マイクロ波放電型中和器のプラズマ生成の特性についての知見を得た．
プラズマ生成については，電子の軌道と電子が感じる磁場・電場強度を解析することにより，ヨーク間に存在する電子のエネルギーが ECR とアンテナ近傍の強い電場によって上昇し，その電子による電離衝突によってプラズマが生成，維持されている事が分かった．アンテナ長がλ/8 の場合は，マイクロ波の放射に優れるλ/4 の長さのアンテナに比べ，マイクロ波の放射の効率は低い．従って，アンテナと ECR 領域の距離が離れた場合は，電子のエネルギー上昇が望めない．しかし，アンテナと ECR
領域を接近させる事により プラズマの生成が可能となり ，結果としてマイクロ波放電型中和器の小型化，軽量化に寄与している．

104３）アンテナを利用したマイクロ波放電型プラズマ源に関する解析
ここでは，小型マイクロ波放電型イオンスラスタやマルチモノポールアンテナシステム を利用した
イオンスラスタ の設計指針を得る目的で，アンテナを用いたマイクロ波放電型プラズマ源の解析を行った．
はじめに，ECR の条件を満たす一様磁場中のモノポールアンテナの解析を行い，電子のエネルギーはアンテナの近傍でのみ上昇し，アンテナと磁場の向きが平行である場合の方がその上昇が大きいという結果を得た．アンテナに与える電場を大きくする事によって，離れた距離にある ECR 領域でも電子のエネルギーが上昇する事が確認されたが，アンテナを ECR 領域に接近させスキンデプスの効果を利用する方が高密度プラズマの生成には有効であると考えられる．磁場強度 を変化させた結果からは2nd ECR の効果も確認でき，この効果はより小型のプラズマ源への応用が可能と考えられる．永久磁石での解析では，アンテナが磁場に対して平行に挿入されている Case 1 と，磁場に対して垂直に挿入されている Case 2 について計算した．電子エネルギーの上昇，EEDF，電離衝突 からプラズマ生成メカ
ニズムの比較を行った．電子エネルギーの空間分布と EEDF から Case 1 の方がキセノンの電離電圧を超えるエネルギーをもった電子が多く存在していることが分かった．電離衝突発生数も Case 1 の方が多い．これらの結果はアンテナを用いたイオンスラスタの設計の際の指針を与える．

また，今後の課題としては以下の 4 項目が考えられる．
１）アンテナを用いたマイクロ波放電型プラズマ源のプラズマ密度の把握
マイクロ波放電型イオンスラスタでは要求される推力（MUSES-C 計画では 8 mN）からマイクロ波の周波数が決定される．この理由は，ECR 放電ではプラズマ密度がマイクロ波の周波数によって決定される事に起因している．MUSES -C 計画ではマイクロ波周波数 4.25 GHz を採用する事で推力 8mN が達成されている．4.25 GHz のマイクロ 波電源の効率は第５章で述べたように約 60 %程度であり，システム全体の効率向上を考えた場合は改善すべき 点である．アンテナを用いたイオンスラス
タではスキンデプス の効果により，カットオフを超えるプラズマ密度を得る事が可能であり，必要な推力を得るために要求される プラズマ密度は，より低いマイクロ波の周波数でも得られる事になる．この点を利用し，マイクロ波電源の効率の良い周波数帯でプラズマ 生成を行えば，マイクロ波
放電型 イオンエンジンの飛躍的 な推進効率の向上が望める．今後はアンテナを用いたマイクロ波放電型プラズマ源の放電室形状やアンテナ位置によってプラズマ密度がどの程度まで大きくできるかを解析する必要がある．

２）カップリングコード の Full-PIC での計算本研究 の目的は電子とマイクロ波の相互作用に重点を置いた解析にあった．対象となるプラズマ
源全領域での計算を行う目的から，PIC 法でメッシュサイズはデバイ長程度に設定するという条件を課さなかった．この条件を課さない事により，PIC 法にはイオンを含めず，静電場を考慮しないこと105とした ．この問題点 として，金属壁 でのシースを模擬する境界条件を必要とする 事と定常状態に至る計算が不可能 である事が挙げられる．１）で述べたような プラズマ密度の上限を決定するような計算では定常状態に達する必要があり，定常状態に至るには イオンの運動によって生ずる電荷分離を考慮する必要がある．本研究 で開発したカップリングコードではイオンを考慮する事も可能である． Full-PIC で計算を行う場合には，計算を可能とする計算体系の選択と大容量のメモリを使用する為に計算の並列化が必要である．

３）様々なアンテナ形状での解析
アンテナは，その形状により指向性 や励起可能なモード が異なり，その事がアンテナを利用したプラズマ源においてはプラズマの生成に大きな影響を与える．本論文で紹介した小型マイクロ波放電型イオンスラスタ ではスター やディスク，一方マルチモノポールアンテナシステム ではループア
ンテナ が採用されている．今後は，これらの形状のアンテナに対しても 計算を行い，実験と比較を行う事で最適なアンテナ形状についての指針が得られるものと考えられる．

４）アンテナでの電子の再結合による損失の評価
FDTD 法ではアンテナに平行な電場成分をゼロにすることにより（すなわち金属として）アンテ
ナをモデル化し，一方 PIC 法においてはアンテナの存在は考慮されていない．本計算の結果からはアンテナを ECR 領域に接近させる事が重要であると述べた．従って，電子のアンテナ表面での再結合による損失を評価する事も重要であると考えられる．今後，PIC 法においてもアンテナを金属として取り扱う（すなわち表面での電子再結合を考慮する）モデルの導入が必要となる．

106

謝辞

本研究を行うにあたり，御指導頂いた九州大学大学院総合理工学研究院，中島秀紀教授に深く感謝いたします ．実験と計算の両方行えた事は大きな経験になりました ．また，貴重な御意見を頂きました，九州大学大学応用力学研究所，中村一男教授，九州大学大学院工学研究院，麻生茂教授に感謝いたします．研究と教育に対する姿勢を学ばせて 頂きました，九州大学大学院総合理工学研究院，渡辺幸信助教授に深く感謝いたします．本研究で使用した FDTD コードをご提供頂きました，宇宙航空研
究開発機構 ，船木一幸助教授に深く感謝いたします ．大分工業高等専門学校 ，鷹尾良行助教授に深く感謝いたします ．研究にとどまらず 様々な事を議論できたおかげで知見が広がりました．梶村好宏講
師に深く感謝いたします．１年半だけでしたが，常に高い目的意識をもって研究に取り組む姿を見て，自分の甘さを痛感する次第です．いつも励まして頂いた山本直嗣助手に深く感謝いたします．残り少ない時間で御恩をわずかでもお返しできればと 思います ．研究に集中できる 環境を整えてくださった
岡安悦子秘書に深く感謝いたします ．中島研究室の宮本尚使氏に深く感謝いたします．宮本さんと実験をしていた頃を懐かしく思います ．現在宇宙航空研究開発機構にて研究中 の桑野寛久氏に深く感謝いたします ．研究の場を移してからは連絡が減りましたが，きっといい 成果を出してくれると期待しています．現在スウェーデンに留学中の林真照氏に深く感謝いたします ．計算機 の事でいろいろとお世話になりました．修士２年の篠原大介氏，児玉章裕氏 ，和泉敬介氏，谷所正彦氏，新屋敷圭祐氏，修士１年の小川信司氏，加来大介氏 ，笠岳幸氏，近岡貴行氏，野田賢治氏，宗貞伸孝氏，学部４年の西島康太氏 ，川渕亮氏，三好誠氏，近藤慎哉氏，金川隆保氏に感謝いたします．後輩の皆さんのおかげで年の差を気にする事無く，研究室生活を送れました．ありがとうございました．そして ，私の研究に関わって頂いた卒業された先輩，同期，後輩の皆さんに 感謝しています．
最後に私を支えてくれた父，母，弟，祖父，祖母に感謝しています．これから恩返しが出来るようにがんばります．

107

付録 A 小型マイクロ波放電型イオンスラスタに
関する実験

A．１

緒言

本論文の第５章においてマイクロ波放電型中和器から電子ではなくイオンを加速放出することで ，
小型のイオンスラスタとして作動させる取り組みがあることを述べた．この背景には，第４章で述べたように，マイクロ波放電型中和器が小型低電力にもか かわらず高密度 プラズマが得られる 特性を持つ事，および宇宙開発において衛星の小型化が積極的に多方面で行われている事がある．小型衛星のミッション では，１）従来の大型衛星に比べ大幅なコストの低下が見込まれる．２）複数の衛星を打ち上げる事でリスクの分散を行える．３）複数の衛星を利用したネットワーク運用を可能とするという利点があり，国家プロジェクトとしてだけではなく，民間企業や大学のプロジェクト も立ち上げられている．ここでは小型衛星への搭載を目指した小型マイクロ波放電型イオンスラスタの実験 A−１，２），
特にアンテナ形状が小型マイクロ波放電型イオンスラスタの性能に与える影響 A−３）について述べる ．

A．２

実験装置

本研究の実験装置の概要図 を Fig. A-1 に，写真を Fig. A-2 に示す．実験装置は主に真空排気系，マ
イクロ波伝送系，イオンスラスタ本体から構成されている．以下，それぞれの系統について説明する．

A．２．１

真空排気系

真空容器は SUS 製で，内径 60 cm，長さ 100 cmの円筒型である．真空容器には大小合わせて 11 箇所の穴が開いており，前後の穴にはイオンスラスタ設置用のフランジ及び正面観察用ガラス 窓フランジが設置されている．
容器側面には真空排気用 (Turbo molecular pump，Cryo pu mp) のフランジが 2 つ，
排気バルブ (Mechanical booster pump) が１つ設置されている．またプローブ差込用フランジが 1 つと
側面観察用 ガラス窓フランジが１つ設置され，予備の穴が３つある ．尚，全実験 を通してこの真空容器は電気的にアースされ基準電位となっている．
使用した真空排気装置を以下に示す．
1.

Rotary pump

310 l/ min : 2 台

2.

Mechanical booster pump

108 m3 /h

3.

Turbo molecular pu mp

300 l/sec

108

4.

Cryo pu mp2000 l/sec

Rotary pump と Mechanical booster pump は粗排気用として使用し，1 Torr まで真空を引き，その後 Turbo
molecular pump によって 10-4 Torr 台まで引く．
最後に Cryo pump を使用して 10-6 Torr 台まで真空を引く．
到達圧力は 7×10-7 Torr，推進剤 Xe ガス 0.2 sccm流入時は 1.2 x 10-5 Torr であった．
真空度の測定には Mechanical booster pump 作動時までの低真空領域 (大気圧〜10-3 Torr) ではピラニ
真空計を使用し，それ以上の真空度領域では電離真空計を使用した．

A．２．２

マイクロ波伝送系

本研究で使用したマイクロ波の周波数は 2.45 GHz である．一般に 2.45 GHz の発振に用いられるマグ
ネトロン管タイプのマイクロ波電源 の場合，低電力領域 では発振が安定しないという特性を持つ．そのため ，本研究 では固体発振器 により発振されたマイクロ波をアンプにより 増幅する形式のものを使用した．
マイクロ波伝送系のフローチャートを Fig. A -3 に示す．マイクロ波は発振器から発振された後，アッ
テネタ −（減衰器）を経て，アンプによって増幅される ．その後アイソレータ，方向性結合器を経て特性インピーダンス 50 O の N 型同軸ケーブルより出力される．真空容器内へは N 型真空用コネクタを使用し，真空容器内で N 型ケーブルからセミリジッドケーブルへと変換され，SMA コネクタを経てスラスタ部へ導入される．投入電力は，パワーモニタにより 計測された 入射波と反射波の差より算出される．反射波はダミーロードにて熱として外部へ放出される．
また，アンテナとの整合はマイクロ波電源から真空容器までの N 型ケーブルの長さと真空容器からスラスタ部までのセミリジッドケーブルの長さを変化させて 行っている．

A．２．３

スラスタ本体スラスタ本体は ISAS/JAXA による小惑星探査機「HAYABUSA」に搭載されている中和器 (電子源) を参考に製作されている．Fig. A-4 にスラスタ本体の断面図を，Fig. A-5 にスラスタの写真を示す．スラ
スタは大きく放電室部，磁気回路部，アンテナ部，グリッドから構成される．放電室 は Al 製のベース
によって固定されている．
放電室部は高さ 12 mm/内径 21 mm の A l 材側壁部 と内径 35 mm/外径 42 mm の Al 材磁石フォルダー部からなる．マイクロ波周波数 2.45 GHz に対応した ECR 磁場強度 0.0875 Tを実現するため，側壁と磁石フォルダーの隙間にサマリウムコバルト (Sm-Co) 磁石が 4 個充填されている．
Sm-Co 磁石はすべて同型の直方体で寸法 5×10×12 mm，長手方向磁化，その表面磁場強度は約 3000
Gauss であった．これを リング状に均等に並べ，カスプ磁場を形成する．希土類磁石はキュリー温度が200 ℃と低温であり，放電室内で直接プラズマに接触すると，大きな熱流束を受けることになるので，磁石の劣化を防止するためにも放電室側壁部が必要となる．放電室前方と後方のヨーク，永久磁石で109磁気回路を構成させるためにヨークには SUSB (磁化軟鉄 ) を用いている．放電室内磁場形状を Fig. A -6に示す．ECR 領域はアンテナが直接接するような場所に形成されているのがわかる．
Mo 材アンテナは SMA コネクタ（female）によって固定されており，マイクロ波電力を放電室内部に伝え放電を起こすものである．アンテナとしてはスター，ディスク，クロス（Fig. A -7）の３種類を使用した．イオンスラスタとマイクロ波電源は Fig. A -8 に示す DC ブロック により絶縁されており，アンテナに高電圧が印加された場合でも，マイクロ波電源には直流高電圧は印加されない．
Table A-1 に 2 枚グリッドシステムの各種パラメータを，Fig. A-9 にグリッドの写真を示す．グリッド
材質は Mo である．グリッド間にはマイカ板が挟まれており，グリッド間の電気絶縁及びグリッド支持を担っている．またアルミナねじによってスクリーングリッド，マイカシート ，アクセルグリッドを支持している．

A．３

実験結果と考察

イオンエンジンの性能を表すパラメータである，イオン生成コスト ε C ，推進剤利用効率η u ，推力
F，比推力 Isp，推進効率η t は以下のように定義される．

Pi − Pr
Ib

εC =

Ib
e
m&
mi

ηu =

2miVb
e
F
I sp =m& g

F = Ib

ηt =

ηu
ε
1+ C
Vb

(A.1)

(A.2)

(A.3)
(A.4)

(A.5)

ここで Pi はマイクロ波の入射電力，Pr は反射電力である．Ib はイオンスラスタから引き出されたイオ& は推進剤流量，Vb はビーム 電圧である．g はンビーム電流である．e は電荷，mi はイオンの質量， m
重力加速度である．イオンビーム電流の計測は Fig. A -10 に示す回路で行った．イオンビーム電流はス
クリーングリッドの電源に流れる電流値からアクセルグリッド の電源に流れる電流を引いたものとした．
Fig. A-11 は３種類のアンテナに対する入射電力とイオンビーム 電流の関係を示している．入射電力& =0.010 mg/sec（Fig. A -11(a)）においてイオンビーム電流はディスクに対して 5.4 mA，スター
10 W， mでは 5.0 mA，クロスでは 4.5 mA であり，ディスクアンテナにおいて最も高い電流値が得られた．一方

110

& =0.029 mg/sec の場合（Fig. A -11(b)）はイオンビーム電流はディスクで 11.3 mA，で，入射電力 14 W，m
スターでは 13.1 mA ，クロス では 12.6 mA であった．

m& =0.029 mg/sec においては（Fig. A-11(c)），ス

ターにおいて最も高いイオンビーム 電流値が計測され，ディスクが最も低い値となった．性能が逆転した理由としては，マイクロ波とプラズマとのカップリングとおよびアンテナ表面での再結合が影響している．ディスクアンテナは ECR 領域に広く接しており，マイクロ波電力が効率良く吸収される．
しかし ，ディスクアンテナ は３種類のアンテナの中で最も大きい表面積 を持つ為にイオンと電子の表
面再結合による 損失も大きい．しかし，推進剤流量 が低く，プラズマ密度も低い場合は，ディスクア
ンテナ とプラズマとの結合が良く，アンテナ表面での再結合による 損失よりもプラズマの生成が勝る
ので，他の 2 つのアンテナよりも高いイオンビーム電流が得られたと考えられる．
Fig. A-12 はスターアンテナでのイオンスラスタの性能を示している．最も良い推進性能としては電力 10 W，推進剤流量 0.02mg/sec において推進効 率 0.56，推進剤利用効率 0.65，イオン生成コスト 1,000
W/A，推力 0.56 mN，比推力 2.5×103 sec が達成された．実用化の為には，イオンスラスタの性能改善と小型中和器の開発が必要である．加えて，軽量で効率の良い電源を使用する事も重要である．

A．４

付録 A のまとめ

小型マイクロ波放電型イオンエンジンの推進性能に対するアンテナ形状の影響を調査した．イオン
ビーム 電流が，ディスク，スター，クロスの３種類のアンテナ形状に対して計測され，最適なアンテ
ナ形状は推進剤流量に依存する事が分かった．この原因は，マイクロ波とプラズマの結合効率の向上とアンテナ表面での再結合 の減少を同時に達成できない事によると考えられる．結果として ，ディスクタイプは低推進剤流量において性能が良く，スターアンテナは高推進剤流量において性能が良かった．最も良い性能として推進効率 0.56，推進剤利用効率 0.65，イオン生成コスト 1,000 W/A，推力 0.56
mN，比推力 2.5×103 sec がスターアンテナにおいて達成された．

111

Fig. A -1

Overview of experimental setup

Microwave generator

Vacuum chamber

DC power supply

Fig. A -2

Experimental setup112

： Forward and Backward Microwave： Forward Microwave
Du mmyL oad

O scillator

A ttenu ato r

Micr owave
AMP

Cir cula to r

Vacuum Chamber

Directi onal
Couple r

FWD

Fig. A -3

REF

Flowchart of microwave transmission system

Fig. A -4

Cross section of thruster

113

E ngineDC blockFig. A -5

Fig. A-6

Photo of thruster

Magnetic configuration in discharge chamber

114

(a)

(b)

(c)

Fig. A-7 Photos of antennas(a) star (b) disc (c) crossFig. A -8 Schematic of DC block

Table A -1

Parameter of grid system

115

Fig. A -9

Photo of grid (Screen, Accele.)

Fig. A-10 Schematic of electric circuit.

116

(a)

(b)

(c)

Fig. A-11 Extracted ion beam for different levels of mass flow rate(a)m& =0.010 mg/sec (b) m& =0.020mg/sec (c) m& =0.029 mg/sec

117

Fig. A-12 Engine performance of star antenna for three levels of mass flow rate

118

参考文献
A−１） 増井博一， 小型マイクロ波イオンエンジンに関する研究 ，九州大学大学院修士論文(2003)A−２） Y. Takao, H. Masui, T. Miyamoto, H. Kataharada, H. Ijiri and H. Nakashima, “Development of
Small-scale Microwave Discharge Ion Thruster”, Vacuum, 73, 449 (2004)A−３） N. Yamamoto, H. Masui, H. Kataharada and H. Nakashima, “Antenna Configuration Effects on Thrust
Performance of Miniature Microwave Discharge Ion Engine”, Journal of Propulsion and Power (to be
published)

119