J. Plasma Fusion Res. Vol.97, No.10 (2021)５
６８‐５７１

小特集

イオン源の数値シミュレーション技術の進展と応用
Progress and Application of the Numerical Approaches for Ion Sources

１．はじめに
1. Introduction
柴 田 崇 統，吉 田 雅 史１）＊
SHIBATA Takanori and YOSHIDA Masafumi１）＊
高エネルギー加速器研究機構，１）山口大学（原稿受付：２
０
２
０年６月２
０日）

１．
１ イオン源の用途・構造と仕組み

数枚取り付けられており，電極間に（一般的に 10−100 kV）

１．
１．
１ イオン源の用途の静電圧を印加することで目的のイオンをビームとして引イオン源は，核融合プラズマ加熱，大型加速器，がん治

き出すことができる．図１には，一般的なイオン源の模式療など医療機器，探査衛星のエンジン，イオン注入など半図を示した．

導体製造など，多分野における粒子源として用いられる．
核融合や医療（陽子線・ホウ素中性子捕捉治療などのがん

１．
２ 水素の負イオン（H−）

治療）分野では，核融合反応を起こすためや，陽子線・中１．
２．
１ 負イオンとは？

性子線の生成のために，水素ガスを導入したイオン源が用水素の場合を例にとると，原子核（陽子）と外殻電子１いられることが多く，正イオン（陽子）または負イオン

つが電気的に釣り合うことで，原子の状態を形成する．こビームが引き出される．加速器の場合，陽子加速器であれのとき電子は，水素の基底準位のエネルギーを占め，そのば，同様に水素イオン源が用いられるが，研究対象となるエネルギー利得は 13.6 eV である（図２）．この原子に新た反応過程によって多様なガス種が使用される．探査衛星のな電子が付着することで，負に帯電したイオン（負イオはエンジンなどでは，プラズマ点灯が容易であることから，最外殻の電子状態１つに空きがあり，そこに電子が付されアルゴンやキセノンが主に用いられる．イオン注入では，ることで全体として比較的に負イオンになりやすい．一半導体に不純物として含まれるホウ素やアルミニウム，リ方，水素やアルカリ金属（I 族）元素のように，同じ周期のンなどのイオンビームが生成される．
１．
１．
２ イオン源の構造・ビーム生成の仕組み
イオン源の構造は用途によっても多様だが，一般的には金属製のチャンバ，ガス導入系，プラズマ点灯のためのパワー入力系（フィラメント，高周波アンテナなど），プラズ
マ閉じ込め用の磁石，イオンビーム引出用の孔の空いた電極によって構成される．
チャンバ内部にガスを導入し，放電開始電圧（アーク電圧や高周波の電磁場）を印加することでプラズマが生成される．プラズマは，カスプ磁石やコイルなどによる磁場によってチャンバ内に閉じ込められる．正イオン源などの場
合，プラズマ中には正イオンと電子のみが主として存在す
図１

る．イオン源チャンバの一端には，孔が空けられた電極が複イオン源の構造とイオンビーム引出の概念図．

High Energy Accelerator Research Organization, Tsukuba, IBARAKI 305-0801, Japan＊現在の所属：宇部工業高等専門学校

corresponding author’s e-mail: takanori.shibata@kek.jp

５６８!2021 The Japan Society of Plasma
Science and Nuclear Fusion Research

Special Topic Article

1. Introduction

T. Shibata and M. Yoshida

水素負イオンの生成過程には，（１）体積生成過程［３］，および（２）表面生成過程［４］の２通りの方法がある．
反応過程の詳細については第２章で述べるが，（１）H−体
積生成では，まず振動励起準位の高い水素分子を生成す
る．これは，比較的エネルギーの高い電子と基底状態の分子の衝突による励起によって積極的に起こる（振動励起過程）．体積生成の２番目の（解離性付着）過程では，振動励起分子 H2（!!）にエネルギーの低い電子が付着すること
で，H−が生成される．上記の過程を促進するため，負イオ
ン源には高温と低温のプラズマ領域を形成する．図３には，体積生成型負イオン源における H−生成機構の模式図
を示した．
図２

一方，図４に示す（２）H−の表面生成過程でも，この２温負イオンの構造とエネルギー準位の概念図．

度 領 域 が 活 用 さ れ る．表 面 生 成 過 程 で は，セ シ ウ ム
中で原子核の価数が小さいものでも，最外殻に別の電子が（Cs）などをチャンバ内に導入し，イオン源引出部の電極

１つ入ることで，エネルギー的により安定した状態を作る（モリブデンなど）に蒸着させることで，金属表面の仕事関ことができる．新たに１つ電子を付した際のエネルギー利数を低下させる．この表面に，水素原子や正イオンが入射得を，電子親和力と呼び，これが正の元素は負イオンになすることで，金属中の自由電子が遷移し，H−が生成されり得る．

る．表面生成H−の親粒子である原子や正イオンは，ドライ

１．
２．
２ なぜ水素負イオンを用いるのか？

バ領域の高温・高密度プラズマで積極的に生成される．一核融合の中性粒子ビーム入射（Neutral Beam Injection;
NBI）加熱では，数 100 keV−1 MeV まで加速したイオンビームを中性化してコアプラズマ中に入射する．ビーム中性化は，中性ガスを充填したセル内にビームを通し，イオンの運動量を中性粒子に移行させるが，MeVクラスのビー
ムの場合，正イオンによる中性化効率は０に近い．負イオンの場合は，中性化効率が６０％近いため，大強度 NBI では
H−が採用されている［１］．
大型陽子加速器や医療用加速器などでは，線形加速器（リニアック）で直線的に加速したイオンビームを，シンク
ロトロンなど円形加速器内で周回する陽子に加えること
で，高輝度のビームバンチを形成する．このとき，周回す
る陽子ビームに，新たに陽子を加えようとするとクーロン反発によって大部分が加速空洞内部に損失してしまう．損
失を防ぐため，最近の大強度加速器では H−のビームをリ図３体積生成型イオン源における H−生成機構．

図４

表面生成型イオン源における H−生成機構．

ニアックで加速し，周回する陽子ビームと合流させてから，ストリッパフォイルなどで電子を剥ぎ取る方式を取る．
このように，H−生成の利点は，特に高エネルギーの物理分野で顕著にみられる．これらの領域では，ビーム大強度化が目的とされることが多く，負イオン源内部で如何に大量のH−を作り，効率よく引き出すかが，装置全体の性能を左右する［２］．
１．
２．
３ 負イオン源の構造
水素の負イオン源は，正イオン源（陽子源）と異なり，水
素原子や正イオンに対して電子付着を行う機構が必要である．これを実現するため，負イオン源にはフィルター磁場と呼ばれるイオン源チャンバを横切る磁場を印加する．この磁場によりイオン源プラズマは２つの領域に分割される．パワー入力系近傍では，高温・高密度のドライバ領域
プラズマが生成され，それがフィルター磁場を横切って引出孔近傍に拡散することで，低温・低密度の引出領域プラ
ズマが形成される．
５６９

Journal of Plasma and Fusion Research Vol.97, No.10 October 2021

方，引出孔付近で生成された負イオンは，電子温度が高い

より，イオン源引出孔近傍で生成した表面生成 H−の電荷と電子衝突によって破壊されるため，引出領域では，プラ分布を考慮した２次元ポテンシャル場を輸送過程に取り入

ズマが低温である必要がある（表面生成型イオン源の場合れられるようになった［９，１０］．

も，体積生成過程は起こるため，引き出されるビーム電流

２０１０年以降，核融合用イオン源では，ラングミュアプの数割程度は体積生成由来の H−である）．

ローブや分光計測などのイオン源実験による結果との再現

体積生成，表面生成型の負イオン源いずれにしても，ド性を得るため，電子の輸送モデルに実体系情報を取り入れ

ライバ領域でプラズマ生成を担うのは主として電子と正イ（後にKEIO-MARC; Kinetic model of Electrons in IOn source

オンである．そのため，負イオン源のプラズマが特別不安with Multi-cusp ARC discharge）コードと命名），また中性

定ということは無い．しかし，ビーム引出孔近傍の電位分

粒子輸送モデルはダイバータ解析に倣い，原子・分子励起−

布には，電子，正イオン，H による空間電荷分布が影響す準位を計算する衝突輻射モデルと統合された［１１‐１４］．まるため，ビーム品質には H−輸送過程の制御が重要となる．

た，H−輸送モデルでは３次元 PIC-MCC モデルが開発され，プラズマ中のビーム引出境界にあたる等電位面（プラ１．
３ 負イオン源のシミュレーション

ズマメニスカス）の形状によって，引き出されるH−の発散

１．
３．
１ シミュレーション立上げの背景（ハロ）成分が形成される機構を解明した［１５，１６］．

冒頭でも述べた通り，水素の負イオンは，核融合 NBI

核融合用イオン源シミュレーションから得られた知見の MeV クラスビームや，大強度陽子加速器など，高エネル

は，他分野のイオン源へと応用された．欧州原子核研究機

ギーの物理領域で重要な役割を果たす．これらの分野にお構（CERN）との共同研究では，大型ハドロン衝突型加速

ける大型装置では，ビーム強度の増加や，エミッタンス低器（LHC）の線形加速器 LINAC4 用に開発が進められる高減（ビームの高収束性）などのビーム性能向上がしばしば周波放電型（RF-ICP）イオン源に対し，チャンバ内の時間目的とされている．上記のために，イオン源内部の物理機

変化する電磁場分布の影響を計算可能とした EM（Electro-

構を定量的に理解し，既存の装置を如何に改良すればよい

Magnetic）PIC-MCC モ デ ル が 確 立 さ れ た
［１７，
１８］．ま

か指針を得る上で，イオン源のシミュレーションは有効なた，前述の KEIO-MARC コードは，医療応用として，ホウ手段と言える．本小特集では，核融合分野とも関連が深い

素中性子捕捉治療や陽電子断層撮影に用いられる体積生成

水素の負イオン源を例として，負イオン源内部の物理機構型負イオン源のプラズマ分布解析に用いられ，H−ビームと対応したシミュレーションモデルを解説していく．

出力増加に向けた設計指針検討に活用されている［１９，２０］．

負イオン源のシミュレーションは，１９９０年代に核融合

上記，近年のシミュレーション応用の成果は，文献［２１］にNBI 装置として，仏 Ecole Polytechnique で試験が行われてまとめられている．

いた Camembert III イオン源
［５］における実験結果を解析

１．
３．
３ イオン源物理過程との対応するためにHatayamaらによって立ち上げられた［６，７］．こ

図５には，イオン源プラズマの主な物理機構・領域と，のとき，核融合ダイバータプラズマの挙動の解明に活用さ対応するシミュレーションモデルを示した．第２章では，れていた中性粒子輸送モデルや Monte-Carlo 法による H−

イオン源全体に関わる原子・分子過程のモデリングを紹介

輸送モデルが初めてイオン源分野に導入された．イオン源

する．ここでは，プラズマの生成に関連する反応過程のほ−か，H−体積生成過程や分光計測による結果を理解するた

内部の原子・分子および体積生成された H の輸送過程を−解析することで，効率の良いH ビーム引出のためには，引めの励起原子・分子に対する解析も説明する．

出孔付近のプラズマ温度が 0.5−1 eV である必要性を数値

第３章では，ドライバ領域と呼ばれるアーク電圧や高周的に明らかにした．

波などを入力してプラズマを生成する領域のシミュレー

１．
３．
２ イオン源シミュレーションの実践

ションを紹介する．この領域では，放電用のフィラメント

上記結果を元に，シミュレーションの活用によって装置や高周波アンテナ，カスプ磁場等の影響で，プラズマのエ

設計指針を得るため，２０００−２０１０年ごろまでは，主としてネルギー分布関数は平衡状態のMaxwell分布から逸脱する

核融合NBI用の水素負イオン源を対象とした種々の計算モ

ため，個々の粒子に着目した解析が行われる．この結果を，デルが開発され，あるいはダイバータプラズマ解析から導第２章の原子・分子モデルと組み合わせることにより，−入された．H の輸送モデルでは，荷電粒子の空間電荷効果

H−体積生成に必要な振動励起分子や，H−表面生成の親粒によるポテンシャル場の影響を計算するための Particle-In-子である水素原子・正イオンの粒子束が計算可能である．

Cell（PIC）法と，衝突を含めた粒子の輸送過程を扱うため

表面生成型負イオン源では，イオン源引出部の電極に Csの Monte-Carlo Collision（MCC）法を統合させた PIC-MCCなどを蒸着させることで，金属内部の自由電子が外部に遷モデルが開発された［７］．また，２
００６年ごろから，大型イ

移し易い状態を作り出す．H−の生成効率は，Cs が電極表

オン源のプラズマを対象に，ドライバ領域における電子の面に 0.5−0.6 原子分蒸着し，かつ Cs と電極間に水素原子が輸送モデリングが，日本原子力研究開発機構（JAEA），核

局在化することで促進される．上記，金属表面からプラズ

融合科学研究所（NIFS）協力のもと開発され，上述の中性

マ中の原子・正イオンへの電子遷移に関するシミュレー

粒子輸送，H−輸送モデルと統合された［８］．H−輸送モデ

ションは，第４章にて解説される．

ルは，従来の PIC-MCC モデルにスケーリング則を適用す生成された H−は，イオン源静電加速部に印加された高ることで，計算領域・粒子数が飛躍的に増加した．これに電圧により，ビームとして引き出される．ビームの初期発
５７０

Special Topic Article

1. Introduction

図５

T. Shibata and M. Yoshida

負イオン源内の領域と，本特集におけるモデリング解説の対応．

［４］Yu.I. Belchenko et al., Investiya of USSR Academy of Science Ser. Fiz. 37, 2573 (1973).
［５］M. Bacal et al., Rev. Sci. Instrum. 67, 1138 (1996).
［６］A. Hatayama et al., Rev. Sci. Instrum. 73, 910 (2002).
［７］T. Sakurabayashi et al., J.Appl. Phys. 95, 3937 (2004).
［８］N. Takado et al., J. Appl. Phys. 103, 053302 (2008).
［９］A. Hatayama, Rev. Sci. Instrum. 79, 02B901 (2008).
［１０］S. Kuppel et al., J. Appl. Phys. 109, 013305 (2011).
［１１］K. Miyamoto et al., J. Appl. Phys. 93, 845 (2003).
［１２］T. Shibata et al., J. Appl. Phys. 114, 143301 (2013).
［１３］T. Yamamoto et al., Rev. Sci. Instrum. 85, 02B125 (2014).
［１４］S. Fujita et al., AIP Conf. Proc. 2011, 050016 (2018).
［１５］K. Miyamoto et al., AIP Conf. Proc. 1515, 360 (2013).
［１６］S. Nishioka et al., J. Appl. Phys. 119, 023302 (2016).
［１７］M. Ohta et al., Rev. Sci. Instrum. 85, 02B113 (2014).
［１８］K. Nishida et al., J. Appl. Phys. 119, 233302 (2016).
［１９］M. Onai et al., Rev. Sci. Instrum. 87, 02B127 (2016).
［２０］S. Yamada et al., AIP Conf. Proc. 2052, 969994 (2018).
［２１］A. Hatayama et al., New. J. Phys. 20, 065001 (2018).

散角，指向性が電極設計と適合していないと，ビームハロと呼ばれる周辺部のイオンが電極と衝突し，絶縁破壊や電極自体の損傷を引き起こす．ビーム発散角は，イオン源内部におけるプラズマメニスカスの形状によって決まる．この形状は，電極間の印加電圧と制限電流と関係を記述したChild-Langmuir 則，およびプラズマ中の空間電荷効果によって決まる．第５章では，引出孔周辺部のモデリングにより，電位分布とプラズマ輸送，H−生成過程を同時に計算する手法を紹介する．

参考文献［１］C.F. Barnett, Atomic Data for Controlled Fusion Research,
ORNL-6086 v1 (1990).
［２］日本原子力研究所東海研究所大強度陽子加速器施設開発センター 大強度陽子加速器プロジェクトチーム，
JAERI-Tech 2003-044（p.1-788），日本原子力研究所．
［３］M. Bacal and G.W. Hamilton, Phys. Rev. Lett. 42, 1538
(1979).

５７１

J. Plasma Fusion Res. Vol.97, No.10 (2021)５
７２‐５７６

小特集 イオン源の数値シミュレーション技術の進展と応用

２．イオン源プラズマの原子・分子過程
2. Atomic and Molecular Processes in Ion Source Plasmas
柴 田 崇 統，澤 田 圭 司１）
SHIBATA Takanori and SAWADA Keiji１）高エネルギー加速器研究機構，１）信州大学（原稿受付：２
０
２
１年６月２
０日）

核融合プラズマ加熱や大型陽子加速器・ホウ素中性子捕捉治療機器などに応用される水素負イオン源では，プラズマパラメータが核融合ダイバータプラズマと近く，電離や解離などの原子・分子過程が類似していることから，核融合周辺プラズマモデリングで培われた水素原子・分子の衝突輻射モデルや中性粒子輸送コードが大いに活用されている．本章で解説するイオン源の中性粒子シミュレーションモデルでは，負イオン源に特徴的な負イオンの体積生成過程・表面生成過程や，それらに影響を与える電子・イオンの電磁場中での輸送・非平衡
エネルギー分布の影響を取り入れるモデル拡張が近年（といっても１５年くらい）の間に進められてきた．また，イ
オン源分光計測結果を再現するモデル開発を行うことで，シミュレーションからイオン源設計指針を与えることが可能となった．
Keywords:
ion source, negative ion, atomic and molecular processes, collisional radiative model, simulation

２．
１ 負イオンの体積生成と表面生成
表１

負イオン源プラズマ中の主要な原子・分子過程を表１に−示す．負イオン源におけるH 生成方法は，体積生成と表面

負イオン源における主な原子・分子過程（添え字の＊は励
起準位を表す）．

生成の２種類に分けられる．体積生成による負イオンのみ
をビームとして引き出すイオン源を体積生成型負イオン源とよび，一方で表面生成によって生成された負イオンを引き出すものを表面生成型イオン源とよぶ．
２．
１．
１ 体積生成型負イオン源
体積生成型負イオン源では，主に，水素分子への電子の解離性付着（表１の過程（３））により負イオンが生成される［１］．振動量子数が !"!"の高振動励起分子では小さな振動量子数の分子より解離性付着が数桁も起こりやすい．高振 動 励 起 分 子 は，電 子 衝 突 に よ り 電 子 基 底 状 態!
!
）か ら 生 成 さ れ た 励 起 分 子 H2（#!!'orH2（%!!&!
）に戻るときに生成さ$!"'）が光を放出して H2（%!!&れる（図１）．このためプラズマ生成部では励起分子の生成が可能となるように 10 eV 以上の高い電子温度が実現するように設計される．
電子の解離性付着は電子のエネルギーが小さいほど起こりやすい．振動励起分子 H2（!"!"）は電磁場に捉われる
ことなく拡散して低電子温度の引出領域に到達し，解離性付着により H−が生成される．
２．
１．
２ 表面生成型負イオン源
表面生成型負イオン源の多くにはセシウム（Cs）オーブ
ンが取り付けられており，高温蒸気として導入された Cs
がイオン源電極（モリブデンなど）に蒸着すると電極金属
High Energy Accelerator Research Organization, Tsukuba, IBARAKI 305-0801, Japancorresponding author’s e-mail: takanori.shibata@kek.jp

５７２!2021 The Japan Society of Plasma
Science and Nuclear Fusion Research

Special Topic Article

2. Atomic and Molecular Processes in Ion Source Plasmas

T. Shibata and K. Sawada

ることが示唆され［６］，水素分子の解離や電離の理解のため，電子励起状態のみを区別したモデルが開発された［７］．
その後，非接触ダイバータプラズマにおける分子活性化再結合の理解のため，さらに振動励起状態を区別したモデルが開発された［８］．分子のモデルの整備はその後 Wünderlich，Fantz らによっても進められている［９］．
体積生成過程による負イオン生成では，分子について電子・振動状態のみを区別した扱いがなされてきたが，解離性付着による負イオン生成の速度係数は回転状態により数桁も変化することが知られている．電子・振動・回転状態を区別するモデルは，最近，澤田らにより開発されている
［４］
．
２．
２．
２ 中性粒子輸送コード
中性粒子輸送コードでは，粒子源から放出された原子や
図１

分子を，表１の反応を考慮しながら追跡する．分子放出の

水素分子のポテンシャル曲線，および H−体積生成過程における振動励起分子生成の模式図．

位置・速度・角度や反応の生起などに乱数を用いて多数個の粒子の追跡を行い，用意した空間メッシュや速度空間表面の仕事関数が低下する．ここに入射した水素原子や正メッシュを用いての追跡粒子の運動の情報を集計する．こイオンが金属中から自由電子を受け取って負イオンが生成

れにより各種の粒子密度や速度分布を算出する．

される［２，
３］（表１の過程（１６））．電極金属表面の過程につ従来，負イオンプラズマでは，基底状態原子と電子基底

いては，本小特集第４章にて詳細な解説がなされている．

状態の振動状態 "!!の分子を追跡粒子としていたが，畑山 ら は 核 融 合 用 負 イ オ ン 源 Camembert III を 対 象 と し

２．
２ 原子・分子過程のモデリングて，振動励起分子も追跡粒子とした中性粒子輸送コードの負イオン生成量の評価のためには，体積生成型では振動開発を行い，負イオンの生成を調べた［１０］．プラズマ中の励起分子密度の空間分布，表面生成型では水素原子・各種

原子や分子が感じる電子温度・密度は，それらの空間的な正イオンの密度の空間分布とそれら粒子の速度分布が重要

運動により刻々と変化する．追跡粒子の各位置での励起状である．表１ではこの観点から反応を整理している．これ

態を経由する生成・消滅は衝突輻射モデルで計算される．

らの反応の起こりやすさを表す反応断面積の情報は文献追跡粒子とプラズマ粒子の弾性散乱については，現在コー［４］にまとめられている．ここではこれらの諸量を算出すドの整備が進められている．

る際に使われる衝突輻射モデルと中性粒子輸送コードを簡

イオン源用途では，柴田らにより，反応レートの計算に単に説明する．

荷電粒子の計算モデル（第３章）から得られたアーク放電

２．
２．
１ 水素原子と水素分子の衝突輻射モデル型イオン源の非平衡電子速度分布関数を代入するモデル改

例えば表１の過程（７）では基底原子からの直接的な電離

良やイオン源壁表面における粒子反射係数を考慮したモデに加えて励起原子からの電離もある．この計算では各励起ル改良がなされている［１１‐１３］．計算機コードの応用例に状態の原子密度を求める必要がある．分子過程でも同様に興味のある読者は参照されたい．

励起状態を考慮する必要がある．
よる励起原子密度の時間変化を記述した連立微分方程式か２．
３ 高周波放電型・表面生成型イオン源 CERN
LINAC4 における分光計測とシミュレーション

らなる．これを解くと各励起状態の原子密度が算出されこれまで負イオン生成を念頭に置いて原子・分子のモデる．また得られた励起原子密度を用いて，励起状態を経由

ルを示したが，これらのモデルは原子・分子の発光強度のする電離や再結合の大きさや，励起原子からの発光強度が計算にも利用される．ここではその応用について紹介す計算される．連立微分方程式を解く際，通常は基底原子密る．

水素原子衝突輻射モデルは，電子衝突や光の自然放出に度，電子温度・密度を与え，これらの時間変化に対する励イオン源のプラズマ挙動解明・装置設計への応用にとし起原子密度の緩和時間が短いことから，これらのパラメーて欧州原子核研究機構（CERN）における大型加速器 LHCタの元での定常状態を仮定して励起原子密度を求める（準の線形加速器（LINAC4）用イオン源の解析を紹介する．

定常近似）．

LINAC4 用イオン源は，高周波放電型（Radio Frequency;
RF）表面生成型イオン源である．

核融合分野では，水素プラズマ中の Balmer 線の発光強

図２に，LINAC４用に開発が進められている表面生成型度の計算のため，澤田，藤本らによって分子解離に伴う励−＊−起原子生成（H2＋e →H＋H ＋e ；添字の＊は励起準位

RF 負イオン源の模式図を示す．イオン源チャンバ周囲にを表す）がモデルに追加された［５］．核融合プラズマで計取り付けられた RF アンテナコイル（図 で は 6 turn）に測された水素原子 Balmer 線発光強度をこのモデルを用い!#" !2 MHz，振幅 !!!100−200 V の交流電圧を印加するて解析したところプラズマ周辺部に中性水素分子が存在すことで，高周波プラズマを生成する．また，アンテナコイ
５７３

Journal of Plasma and Fusion Research Vol.97, No.10 October 2021

フィルタは Balmer-!光波長 650−660 nm のみを透過するものであり，PMT 入射窓に取り付けた．PMT はオフラインで，複数のレーザー光，タングステン（W）ランプ（発光強度が既知の光源）を用いた波長校正および絶対感度校正を行った．絶対感度校正では，２次電子電流値に対する入射光子数の依存性グラフを取り，実際のイオン源プラズマ
診断時に PMT に入射する Balmer-!光子の総数を測定可能とした．この計測ポートに対し，Balmer 発光線強度を計測したところ，図４に示すような，大小の非対称な強度ピー
クが RF の１周期の間に現れた．図３に示したように RF
図２

プラズマでは１周期の間に ICP-RF 電場 !"が２回最大値を

CERN LINAC4 用 RF イオン源，および外部 RF アンテナに印加電流の時間変化と特徴的な位相 A，B．

取ることから，プラズマ密度の変化に伴い，Balmer 発光線も同じ強度の２つのピークを取りそうだが，軸上の分光ルに流れる交流電流の１周期にわたる時間変化と後述の議ポートからの計測ではそのようになっていない．

論で発光線強度ピークに対応する位相（A，B）を図２下部

２．
３．
１ RF プラズマシミュレーションモデル

に示す．従来，LINAC4 イオン源では，RF パワー投入時のRF イオン源のプラズマ解析では，以下の３つのモデル

プラズマ点灯に関する不安定性があり，連続運転を行う際を同時に扱う必要がある；の稼働率が８割程度に留められる問題があった．この原因（モデル A）ICP-RF 電磁場や外部電磁場の形成過程

をシミュレーションから探ることが慶應大グループに依頼（モデル B）電子・正イオンの輸送過程された．

（モデル C）原子・分子の輸送過程

RF イオン源のプラズマ生成機構は，大きく２段階に分

図５はシミュレーションモデルの模式図である．上記，

けられる．図３左側に示す「E モード（容量結合プラズ（モデル A），（モデル B）では，第３章で解説されるイオ

マ）」状態では，アンテナコイルに印加した交流電流 "#'$%ン源ドライバ領域の荷電粒子輸送モデルが確立しており，によって，誘導結合型（Inductively Coupled Plasma; ICP）

LINAC4 で は，EM PIC-MCC（Electro-Magnetic Particle-‐RF 電 磁 場 が 形 成 さ れ る．こ の と き，プ ラ ズ マ 電 流

In-Cell with Monte-Carlo Collisions）法に基づいた計算モデ

"
(%
"
)
&"
!
" は，ICP-RF 電場 !" による荷電粒子の加速を受け，ルが適用された．詳細な説明は第３章および参考文献［１４］イオン源周（"）方向に形成される．イオン化が十分起こり，プラズマ密度が成長すると，"
(%
"
)
&"
!
#
'$
%の状態に移" !"行する．すると，ICP-RF 電磁場の形成は， "
#
'$
%ではなく
"
(%
"
)
&"
!
(%
"
)
&"
!
" によって担われる．一方， "" 自身は ICP-RF
電場による加速を受けるために共鳴状態が形成され，プラ
ズマ密度は飛躍的に成長する．これを「H モード（誘導結合プラズマ）」と呼ぶ．
図２に示すイオン源チャンバ #軸上の分光計測ポートには，光ファイバ・光学フィルタを介して，光電子増倍管（PMT）が 取 り 付 け ら れ，そ の 出 力 を 測 定 す る こ と でBalmer-!光に対する強度とその時間変化を調べた．光学

図４

図３

RF プラズマ E モードおよび H モードにおけるコイル電流，
プラズマ電流および電磁場の振動方向．

LINAC4 イ オ ン 源 フ ォ ト マ ル 分 光 計 測 か ら 得 ら れ た
Balmer-!線強度の時間構造．

図５

５７４

RF プラズマシミュレーションモデルの模式図．

Special Topic Article

2. Atomic and Molecular Processes in Ion Source Plasmas

T. Shibata and K. Sawadaに預けるが，EM PIC-MCC モデルでは，荷電粒子１つ１つ

ドから得られた原子密度と分子密度を用いて，電子衝突にの運動方程式を，衝突過程（電離などの非弾性衝突やクーよる基底原子からの励起原子生成および分子解離に伴う励ロン衝突）を含めて逐次計算を行うと同時に，プラズマ電

起原子生成を考慮した水素原子衝突輻射モデルにより計算

流"
'%
$
(
&$によってイオン源各位置に形成される ICP-RF 電

される．図７−９より，実験の Balmer-"線発光強度の増加

磁場を，時間変化を考慮した Maxwell 方程式より計算すは電子の平均エネルギーの増加によることが明らかになっる．その結果，イオン源の各位置における電子と正イオン

た．発光強度が変化しているのは，電子温度（平均エネルの平均エネルギーと密度の分布を求めることができる．

ギー）とそれに伴う電離による電子密度の変化による．原

（モデ ル C）の 原 子・分 子 の 輸 送 過 程 は，前 述 の と お子・分子密度はこの間一定とみなせる．図１０は，ビュー

り，原子・分子の衝突輻射モデルと中性粒子輸送コードにポートから観測される Balmer-"線強度の RF２周期間にわより扱われる．EM PIC-MCC モデルから得られた電子・正イオンの平均エネルギーと密度の空間分布をこれらに代入することで，原子・分子密度の空間分布を求めることができる．この情報を再度 EM PIC-MCC モデルに代入すること
で逐次計算を行う．原子・分子モデルから，さらにイオン源各位置におけるBalmer-"線発光強度を算出し，分光計測
結果の解明を試みた．
分光計測で得られた Balmer 強度のピーク構造の解明のため，軸方向（図２の $方向）への輸送機構の可能性を検討した．高周波アンテナコイルの端子間に印加する電圧をシミュレーション体系に取り入れ，それによってイオン源
内部に形成される容量性電場（!&）（図６）をシミュレー
ションに取り入れた［１４］．その結果，RF プラズマ中の電子が $軸方向加速され，図３に示した #方向のプラズマ電

図７

RF プラズマシミュレーションから予測された位相 A，B
における電子平均エネルギーの空間分布．図の R ＝ 0 はイオ
ン源中心軸，Z ＝ 0 は天板位置に対応する．

図８

位相 A，B における RF イオン源内部の Balmer-"線光子生成レートの空間分布．

図９

位相 A，B における RF イオン源内部の Balmer-"線光子生成レートの空間分布．

流が損なわれるため，プラズマ不安定性によって点火効率が低くなることが示唆された．図７は，イオン源に印加した RF １周期中のある位相 A とその逆相（位相 B）における，イオン源内部の電子平均エネルギーの空間分布（計算
領域は，図２における #!!，$!!領域）である．
２．
３．
２ 水素原子密度・Balmer-"光子生成量の空間分布シ
ミュレーション
上述の中性粒子輸送コードにより，原子・分子密度の空間分布を算出した．さらに水素原子衝突輻射モデルを用いて主量子数 %!#の励起原子密度を求め，Balmer-"（%!#から %!"への遷移）発光強度を算出した［１５］．図８に，図７の位相 A，B における基底原子密度の空間分布を示す．
また，図９には，イオン源各位置における Balmer-"発光強度分布（単位時間・単位体積当たりの光子放出）を示す．
RF プラズマ中の Balmer 線発光強度は，中性粒子輸送コー

図６

RF アンテナコイルへの印加電圧によってイオン源チャン
バ内部に生じる容量性電場（Ez）の解析結果．

５７５

Journal of Plasma and Fusion Research Vol.97, No.10 October 2021

たる時間変化の計算結果である．位相 A，B において発光
線強度ピークが現れる時間構造が実験とシミュレーション
でおおよそ一致したことから，２．
３．
１節末部に述べた，容量性電場による !軸方向へのプラズマ加速による RF プラズ
マ H モード形成が損なわれる機構の妥当性が確 認 された．本結果を踏まえ，LINAC4 イオン源では，ガス流量調整等によるプラズマ点火効率の増強や，コイル巻き数・構造の見直しによる容量性電場緩和対策が進められた
［１６］．

２．
４ おわりに負イオン源における負イオン体積生成量および表面生成量を解析する上で，原子の衝突輻射モデルおよび電子状
態・振動状態分子を扱う衝突輻射モデル，および中性粒子
輸送コードの活用は極めて重要である．このような中性粒

図１０ シミュレーションによる RF サイクル２周期の Balmer-!
線強度の計算結果．

子に対するモデルの信頼性を得るためには，実験による測定結果と計算結果を直接比較し，仮定した物理モデルを確

ence Ser. Fiz. 37, 2573 (1973).
［３］M. Bacal and M. Wada, AIP Conf. Proc. 1515, 41 (2013).
［４］K. Sawada and M. Goto, Atoms 4, 29 (2016).
［５］T. Fujimoto et al., Nucl. Fusion 29, 1519 (1989).
［６］T. Fujimoto et al., J. Appl. Phys. 66, 2315 (1989).
［７］K. Sawada and T. Fujimoto, J. Appl. Phys. 78, 2913 (1995).
［８］K. Sawada and T. Fujimoto, Contrib. Plasma Phys. 42, 603(2002).
［９］D. Wünderlich and U. Fantz, Report IPP 10/18 (2001).
［１０］A. Hatayama et al., Rev. Sci. Instrum. 73, 910 (2002).
［１１］T. Shibata et al., J. Appl. Phys. 114, 143301 (2013).
［１２］T. Shibata et al., Plasma Fusion Res. 9, 1401011 (2014).
［１３］M. Onai et al., AIP Conf. Proc. 1869, 030043 (2017).
［１４］S. Mattei et al., Rev. Sci. Instrum. 85, 02B115 (2014).
［１５］T. Shibata et al., AIP Conf. Proc. 1655, 020008 (2015).
［１６］J. Lettry et al., AIP Conf. Proc. 1655, 030005 (2015).

認する必要がある．本章で紹介した例では，近年核融合や加速器分野で多く活用される RF 型負イオン源に対して分光計測を行うとともに，衝突輻射モデル・中性粒子輸送
コード（本章）＋荷電粒子輸送モデル（第３章）による計算を行った．測定および計算結果の比較から，これらによる計算が実験結果をおおむね再現することを確認した．
このような成果は，最近の計算機性能向上によるところもあるが，近年の中性粒子と荷電粒子モデリング双方における１つ１つのモデリング改善の積み重ねが大きく貢献し
ている．

参考文献
［１］M. Bacal and G.W. Hamilton, Phys. Rev. Lett. 42, 1538
(1979).
［２］Yu.I. Belchenko et al., Investiya of USSR Academy of Sci-

５７６

J. Plasma Fusion Res. Vol.97, No.10 (2021)５
７７‐５８０

小特集 イオン源の数値シミュレーション技術の進展と応用

３．ドライバ領域におけるモデリング
3. Modeling in Driver Region
星野一生
HOSHINO Kazuo
慶應義塾大学理工学部（原稿受付：２
０
２
１年６月２
８日）

負イオン源内部では，アーク放電や高周波放電によりプラズマは高エネルギーテールを持った非平衡状態にある．さらにカスプ磁場，フィルター磁場による，複雑な磁場形状における輸送を考える必要がある．そのため，負イオン源内のプラズマ解析には運動論的シミュレーションが有効である．本章では，このような運動論的シ
ミュレーションコードとして，アーク放電型マルチカスプ負イオン源を対象とした KEIO-MARC コード，RF-ICP
負イオン源を対象とした EM-PIC-MCC コードについて，モデル概要を解説し，それらを使った解析例について紹介する．
Keywords:
negative ion source, kinetic simulation, KEIO-MARC, ion source with multi-cusp Arc-discharge, EM-PIC-MCC, RF-ICP

３．
１ はじめに析，電子エネルギー分布関数（EEDF: Electron Energy Dis-

第１章で述べられたように，負イオン源装置内はフィル

tribution Function）解析に用いる．

ター磁場によりドライバ領域と引出領域に分けられてい粒子軌道は以下の運動方程式により追跡する．

る．ドライバ領域では，アーク放電や高周波放電によりプ#!#"!#"!%! " #"$!! ! "!
# !$##"#

ラズマが生成される．生成されたプラズマは，Maxwell
分布から逸脱した非平衡状態となり，高エネルギーテール（１）を持つ．この高エネルギーテールは，水素負イオン生成に

ここで，!，#，"，"，!は，それぞれ，質量，粒子速度，

必要な振動励起分子の生成等において重要な役割を果た

電荷，電場，磁場である．右辺の最終項は弾性・非弾性衝

す．また，フィルター磁場やカスプ磁場が存在する場合，

突による運動量変化を表す．装置内のほとんどの空間ではプラズマの輸送過程や壁への損失過程は複雑になる．

準中性条件が成り立つとし，通常の解析では電場について

このような系では，個々の荷電粒子による運動論的な振は無視する．実磁場配位のデータを与えることで，カスプる舞いが重要となる．そこで，本章では，主にドライバ領

磁場によるミラー効果も含め軌道追跡が行われる．ローレ

域におけるプラズマを解析するために開発された運動論的ンツ力による旋回運動は数値計算上誤差を生みやすいが，シミュレーションコードとして，アーク放電型イオン源

Buneman-Boris 法［６］を用いることで精度良く粒子軌道をを対象とした KEIO-MARC コード［１‐３］，高周 波 放電型

追跡する．カスプ磁場を越えて壁境界に到達した粒子につ（RF-ICP）イ オ ン 源 を 対 象 と し た EM-PIC-MCC コ ー ド

いては，シース電位より高いエネルギーの粒子は損失し，［４，５］について，その概略とシミュレーション例についてエネルギーの低い粒子は反射される．
衝突過程はモンテカルロ法を用いて模擬する．電子‐電紹介する．

子間および電子‐H＋間のクーロン衝突は，２体衝突モデル

３．
２ アーク放電型イオン源のモデリング［７］により模擬する．電子とその他の水素粒子（H，H2，

３．
２．
１ シミュレーションモデル

H+2 ， H+3 ）との衝突については，およそ５４０種の反応を考慮

３．
２．
１．
１ 電子輸送及び電子エネルギー分布関数解析モデル

可能であり，ヌル衝突法［８］により多数の反応を計算リ

KEIO-MARC コード（Kinetic model of Electrons in IOn

ソース的に効率よく取り扱うことができる．背景粒子（H，source with Multi-cusp ARC discharge）［１‐３］は，テスト粒

H2，H＋，H+2 ，H+3 ）の密度，温度については入力データとし

子モンテカルロコードであり，３次元の実空間および速度て与える．

空間で荷電粒子の軌道追跡を行う．実形状で，カスプ磁場

フィラメントから発生した熱電子は，シース領域で即座を含む実磁場形状を考慮可能である．アーク放電型イオンにアーク電圧まで加速されるとして，１次電子の初速度を源のモデリングでは，準中性を仮定し，主に電子の軌道解与える．その後，中性粒子との衝突によりイオン化が起き，

Faculty of Science and Technology, Keio University, Yokohama, KANAGAWA 223-8522, Japanauthor’s e-mail: hoshino@appi.keio.ac.jp

５７７!2021 The Japan Society of Plasma
Science and Nuclear Fusion Research

Journal of Plasma and Fusion Research Vol.97, No.10 October 2021

発生した２次電子と合わせ，系内が定常状態になるまで計算を行う．
以上の一連の計算で得られた粒子情報を統計処理することで，EEDF の他，電子密度分布，電子温度分布，装置内の平均流速分布などを評価することができる．
３．
２．
１．
２ 振動励起分子および負イオンの評価
KEIO-MARC では電子輸送の解析のみ行うため，負イオンの生成については直接評価できない．負イオン生成量の評価には，さらに２段階の解析を行う［９］．
まず，負イオンの生成に重要となる振動励起分子密度については，レート方程式を用いて評価する．#!!""#の
図１

準位の振動励起分子に対する合計１５本のレート方程式を連アーク放電型マルチカスプ負イオン源モデル［９］．

立して解くことで，各振動励起準位の密度を得る．通常，
レート方程式で用いる反応レート係数には，Maxwell 分布で積分された平均的な係数が用いられることが多い．一方，前述のとおり，イオン源内では電子は非平衡状態であり，高エネルギーテールが原子・分子過程において重要な役割を果たす．そこで，各種過程の反応レート係数 !('$%
は，KEIO-MARC で 評 価 し た ド ラ イ バ 領 域 の EEDF，!はエネルギー）を用いて，以下により評価する．
（"!$'#!('$%!!"!$"('$%#
!$
#'#
!$
&!
'#

（２）

ここで，"('$%#!$はは各種反応過程の反応断面積，#'#!$は !に対応する電子の速度である．

図２

KEIO-MARC で得られた EEDF 及びその z 位置依存性［９］．

次に，KEIO-MARC で得られた引出領域の EEDF 及び
レート方程式により評価した振動励起分子密度を用いて，れると考えられ，イオン源内の EEDF に対して中性粒子と粒子バランスから水素負イオン密度を評価する．生成項の非弾性衝突が重要な役割を果たしていることが示唆されは，各振動励起状態の分子と電子による解離性付着反応とる．

なる．消滅項としては，電子脱離反応，相互中性化反応，結その後，$!40 mm 付近まで進むと衝突・緩和過程によ

合性脱離反応，非結合性脱離反応に加え，輸送損失項を加りピークは小さくなっていくが，高エネルギー成分は残っえている．これら，生成項と消滅項のバランスから，負イ

たままである．さらに引出方向へ進み，プラズマ電極付近

オン密度を評価する．

（$!10 mm）になると電子密度が大幅に減少し，特に高エ３．
２．
２ シミュレーション例

ネルギー成分が著しく減少する．これは，$!40 mm 付近

KEIO-MARC を使った EEDF 解析の例として，中性子補から引出領域に向けて存在するフィルター磁場の影響であ足療法（BNCT: Boron Neutron Capture Therapy）による

る．フィルター磁場領域では，高エネルギー電子は透過しガン治療で用いられる水素負イオン源（住友重機械工業株

にくく，フィルター磁場領域を横切る前に磁力線に沿って式会社）の解析結果［９］を紹介する．このイオン源は，円筒壁へと損失する．この結果，高エネルギー成分が低減する

形のマルチカスプアーク放電型のイオン源であり，その概ので，引出領域では低温電子領域が形成され，効率的な負略図を図１に示す．ビーム引出軸方向を $とし，$!0 mm

イオン体積生成を実現している．
次に，３．
２．
１．
２節の手順により評価された負イオン密度がプラズマ電極である．解析では，アーク電流を 2240 W，を図３に示す［１０］．ここでは，モデルの中で大きな仮定のアーク電圧 80 V，水素ガス圧 1.6 Pa としている．
解析の結果得られた各 z 位置における EEDF を図２に示一つである負イオンの閉じ込め時間を横軸に，縦軸に負イす．全体的に熱緩和した低温成分と，そこから逸脱した高オン密度をプロットしている．また，異なる３ケースのエネルギー成分を持ち，電子が非平衡状態であることがわ

フィルター磁場強度の結果が比較されている．フィルター

かる．特に，フィラ メ ン ト に 近 い $!120〜80 mm では，

磁場強度はAからB，Cと進むにつれ強くなっており，フィフィラメントから発生した高エネルギーの１次電子の寄与

ルター磁場強度が強くなるほど負イオン密度が減少していが大きく，１次電子の初期エネルギーである 80 eV 付近にることがわかる．これは，主にフィルター磁場により電子ピークを持つ．また，65 eV 付近にもピークが見られる．電密度が減少したことに起因している．本モデルでは，負イ子の非弾性衝突の中で最も頻度の高い分子のイオン化反応

オン密度を評価するために粒子バランスモデルを用いてい（非解離性イオン化）の閾エネルギーは約 15 eV であり，こるが，そこで仮定する閉じ込め時間に対して，負イオン密の反応により初期の80 eVに加え，65 eVにもピークが見ら度が強く依存していることがわかる．より適切な閉じ込め
５７８

Special Topic Article

3. Modeling in Driver Region

図４

K. Hoshino

RF-ICP イオン源のシミュレーションモデル例［１１］．

境界条件［１３］を用いる．半径方向（!方向）については，
Bayliss Turkel の放射境界条件［１４］を用い，!$!では軸対
図３

称により境界条件を定めている．

H−の閉じ込め時間及びフィルター磁場強度に対する H−密度依存性［１０］．

次に，RF ソレノイドアンテナが作る電磁場（"$#，!$#）
については，３次元実形状において高周波３次元電磁界解

時間については，外部からの引出電場のしみ込みや，負イ

析ソフトウェア ANSIS HFSS を用いた計算を行う［１５］．

オンを含むシース形成（イオン―イオンシース）など，複雑得られた電磁場分布は $方向に平均化する．

な物理過程が寄与する．第５章で述べられる引出領域のシプラズマ電流と RF ソレノイドアンテナが作る電磁場のミュレーション結果等に基づいて，適切なモデル化が必要

寄与を足し合わせ，軸対称２次元の電磁場分布として次にである．

説明する粒子輸送解析で用いる．
（２）粒子輸送解析モデル
上記で求めた電磁場中における荷電粒子 （e，H＋，H＋
2 ）３．
３ 高周波放電型イオン源のモデリング
３．
３．
１ シミュレーションモデル

の輸送解析を３．
２．
１節で述べた KEIO-MARC と同様の手順

高周波放電型（RF-ICP）イオン源におけるプラズマの挙で行う．電磁場中の荷電粒子の運動方程式は Buneman-

動を解析するためには，１）電磁場中における荷電粒子のBoris 法を用いて解く．背景中性粒子との衝突については

軌道解析，２）RF コイルによる時間変動する電磁場とプラヌル衝突法［８］を用い，弾性衝突，H，H2 のイオン化，解ズマ自身により発生する電磁場の解析，３）イオン化反応

離性イオン化，解離性再結合を考慮する．一方，荷電粒子等をはじめとする各種衝突過程，を自己無撞着に解析する

間のクーロン衝突については，現在のモデルでは無視して必要がある．この様な解析を行うために，３次元空間の軌

い る．こ れ は，RF-ICP イ オ ン 源 で は H２ガ ス 圧 が 高 い

道追跡モデル，軸対称を仮定した２次元空間における電磁（1−10 Pa）一方，プラズマ密度が低く（＜1016 m−3）
［１６］，

場解析，そして，各種衝突過程を扱うモンテカルロモデル

クーロン衝突よりも中性粒子との衝突が支配的となるためを 結 合 し た，2.5 次 元 EM-PIC-MCC コ ー ド（Electro-

である．

Magnetic Particle-In-Cell with Monte-Carlo method for Col-

この軌道追跡では，荷電粒子を超粒子として扱い，得ら

lision processes）が開発された［４，５］．以下，この EM-PIC-れた位置・速度から，各超粒子（イオン，電子）が作る電MCC モデルの概略について述べる．

流を空間格子点へ割り付ける．さらに $方向に平均化する（１）電磁場解析モデル

ことで，（１）電磁場解析で必要となるプラズマ電流 %'&の空

こ こ で は，図４に 示 す よ う な 欧 州 原 子 核 研 究 機 構間分布が得られる．

（CERN）の線形加速器 Linac4 で用いられる円筒形の RF-

ここで得られた %（１）の電磁場解析を行い，求'&を用いてICP イオン源を想定する．電磁場解析では，プラズマが作まった電磁場分布を用いて（２）の粒子輸送解析を行う．こる電磁場と RF ソレノイドアンテナが作る電磁場とを分けのステップを繰り返すことで，電磁場と粒子輸送の自己無て個別に解析する．

撞着な解析を行う．
３．
３．
２ シミュレーション例

まず，プラズマが作る電磁場（"'&，!'&）については，軸対称（%!
%$$!）を仮定し，２次元円筒座標系（!-#平面）でRF-ICP イオン源では，２つの運転モードがある．一つ以下の電磁場解析を行う．

は，低密度状態で RF ソレノイドコイルを横切る高い電位による容量結合プラズマ状態（Eモード）である．これに対%!'&%#"'&$!
%
"

（３）

%"'&
%#!'&$#!%
!#
!
'&""
%
"

（４）

し，プラズマ密度が高くなると，RF 電流に誘起される電場により放電が維持されるようになる（H モード）
．プラズ
マ密度の増加により，E モードから H モードへと E-H 遷移が起こる．この E-H 遷移の前後でプラズマの特性が大きく

ただし，"
!及び #
!は，それぞれ真空中の誘電率及び透磁率変わるため，その理解は，特にRFシステムとプラズマ負荷である．プラズマ電流 %'&は，粒子輸送解析により各格子点のインピーダンスマッチングに対して重要な課題となってに割り付けられた電流値を用いる．これを時間領域差分法

いる．

（Finite Difference Time Domain Method: FDTD 法）［６，１２］

このような E-H 遷移を理解するために，CERN の Linac4により解く．軸方向（#方向）の境界条件には，Mur の吸収

イオン源を対象とした EM-PIC-MCC コードを用いたシ
５７９

Journal of Plasma and Fusion Research Vol.97, No.10 October 2021

ミュレーション解析が行われた［５，１７］．従来の解析は，流ミュレーションモデルである KEIO-MARC コードおよび体モデルや等価回路モデルによる解析が主であり，運動論

EM-PIC-MCC コードについて，モデル概要とシミュレー

的シミュレーションとしては初めての試みであった．解析

ション例についてまとめ，医療用負イオン源や加速器用負の結果，E モードでは，RF コイルが作る軸方向（z 方向）のイオン源のシミュレーション例について紹介した．ここで強い電場 !"%!により，質量の軽い電子が電場と逆方向へ移紹介した以外にも，中性粒子輸送コードや分光データとの動する一方，イオンは慣性により余り動かないことがわ

比較のための衝突輻射モデルなどの開発も進められていかった．この電子の挙動により，!"%!と逆方向にプラズマる．また，核融合プラズマ加熱用の負イオン源シミュレー

が作る電場 !$%#が形成される．この !"%!と軸方向の電子の運ションを始め，様々なシミュレーション解析が行われ，実動が，E モードプラズマにおける重要な特徴となる．

験との比較によるモデル妥当性の検証や，多くの重要な知プラズマ密度が高くなり，この !$%#が !"%!と同程度の大き見が得られている．（例えば，［２０‐２３］）

さになってくると "方向の電場が弱まり，!方向の RF 電場!"!!が支配的になる．プラズマは RF コイルの巻かれた領域

参考文献に局在化するようになり，容量結合プラズマから誘導結合［１］I. Fujino et al., Rev. Sci. Instrum. 79, 02A510 (2008).
［２］R. Terasaki et al., Rev. Sci. Instrum. 81, 02A703 (2010).
［３］T. Shibata et al., J. Appl. Phys. 114, 143301 (2013).
［４］T. Hayami et al., AIP Conf. Proc. 1390, 339 (2011).
［５］S. Mattei et al., Rev. Sci. Instrum. 85, 02B115 (2014).
［６］C.K. Birdsall and A.B. Langdon, Plasma Physics via Computer Simulation (New York: McGraw-Hill, 1985).
［７］T. Takizuka and H. Abe, J. Comput. Phys. 25, 205 (1977).
［８］K. Nanbu et al., Proc. Soviet Union-Japan Symp. Comput.
Fluid Dynamics, (Computing Center of the USSR Academy of Sciences), 2, pp.126-132 (1989).
［９］S. Yamada et al., AIP Conf. Proc. 2052, 060004 (2018).
［１０］S. Yamada et al., Plasma Fusion Res. 14, 3401160 (2019).
［１１］W. Kobayashi et al., AIP Conf. Proc. 2052, 050009 (2018).
［１２］R.W. Hokney and J.W. Eastwood, Computer Simulation Using Particles (Bristol: IOP Publishing, 1988).
［１３］G. Mur, IEEE Trans. Electromagn. Comput. 23, 377 (1981).
［１４］A. Bayliss and E. Turkel, Commun. Pure Appl. Math. 33,
707 (1980).
［１５］A. Grudiev et al., Rev. Sci. Instrum. 85, 02B134 (2014).
［１６］S. Mochizuki et al., AIP Conf. Proc. 1655, 020016 (2015).
［１７］S. Mattei et al., Plasma Sources Sci. Technol. 25, 065001(2016).
［１８］K. Nishida et al., Rev. Sci. Instrum. 85, 02B117 (2014).
［１９］W. Kobayashi et al., Plasma Fusion Res. 14, 3401094 (2019).
［２０］N. Takado et al., J. Appl. Phys. 103, 053302 (2008).
［２１］T. Shibata et al., AIP Conf. Proc. 1515, 177 (2013).
［２２］T. Shibata et al., Plasma Fusion Res. 9, 1401011 (2014).
［２３］矢 口 真 子 他：第３６回 プ ラ ズ マ・核 融 合 学 会 年 会，
02P24（春日井，２
０
１
９年１
２月）
．

プラズマへの遷移（E-H 遷移）がシミュレーション上で確認された．
このように，RF-ICP イオン源のプラ ズ マ 解 析 に EMPIC-MCC モデルが有用であることが示されたが，高い計算コストが課題である．例えば，RF-ICP イオン源では，効率的に電力を注入するために，RF システムと RF-ICP のイ
ンピーダンスマッチングが重要となるが，EM-PIC-MCC
モデルは計算コストが高いため，RF システムモデルを組
み込み，パラメータサーベイによりその最適化を検討することは難しい．そこで，このような目的のために，整合回路の等価回路モデルを用いた０次元モデルも開発されている［１１，１８］．プラズマ部分については Maxwell 分布を仮定
した流体モデルであるが，計算コストが低く，パラメータサーベイに適している．この０次元モデルと EM-PIC-MCC
シミュレーションの比較が行われた［１１］．定常状態では両モデルでプラズマインピーダンスはよく一致し，プラズマ
インピーダンスのパラメータサーベイ等では非常に有用なツールとなることが示された．実際，プラズマインピーダ
ンスの各種回路パラメータ依存性の解析などが行われている［１９］．一方で，プラズマ着火直後は，プラズマは非平衡状態であり，両モデルで結果は大きく異なってしまうた
め，解析目的に合わせ，適切にモデルを選択する必要がある．

３．
４ おわりに本章では，主にドライバ領域を対象とした運動論的シ

５８０

J. Plasma Fusion Res. Vol.97, No.10 (2021)５８１‐５８５

小特集 イオン源の数値シミュレーション技術の進展と応用

４．プラズマとイオン源壁面との相互作用
4. Plasma-Wall Interactions in Ion Sources
和田元，剣 持 貴 弘１）

WADA Motoi and KENMOTSU Takahiro１）
同志社大学 理工学部，１）同志社大学 生命医科学部（原稿受付：２
０
２
１年６月１
１日）

イオン源の壁面状態はイオン源プラズマ中のイオン種比率を決定する上で重要な役割を演じる．特に金属壁の仕事関数は二次電子放出や負イオン生成量を大きく左右する．イオン源プラズマ‐壁面相互作用を模型化するのに必要となる素過程データは実験的に採取するのが難しいことが多く，表面反応シミュレーションが有効なツールとなるが実験と同様，低エネルギー領域においては数値シミュレーションを行うことも難しくなる．本章
では負イオン源プラズマと壁面の相互作用に関するシミュレーション研究について概略を纏める．
Keywords:
negative ion, surface production, particle reflection, ion induced desorption, sputtering

４．
１ はじめに計算するために用いられることになる．

プラズマ―固体壁相互作用と言われると，多くの読者はの損耗現象を思い浮かべることであろう．もともと装置自４．
２ イオン源プラズマー壁面相互作用シミュ
レーションの現状

体の容積／内壁表面積比の小さなイオン源性能は，壁面状

４．
２．
１ 負イオン生成電極の水素反射シミュレーション

核融合プラズマ閉じ込め装置のリミターやダイバーター等

態に大きな影響を受ける．正イオン源を用いてトカマクプ

Schneider ら［５］はエネルギーの揃った水素・重水素のラズマの中性粒子加熱を行っていた１９８０年代から，イオン

原子・分子イオンを金属表面に打ち込み，負イオン表面放源プラズマの数値シミュレーションは行われていた［１］．

出係数として測定した．Cs 付着量を変化させて仕事関数と

この中でイオン源壁面での水素原子再結合係数がイオン源

重水素負イオン収率の変化を測定するとともに，厚いアルから引き出されるイオン種を決定する重要な因子として扱

カリ金属表面での水素負イオン化収率を測定した．Hiskes

われていた．負イオン源については，壁面状態がさらに直［６］らは表面で水素粒子の後方散乱過程によって形成され接的にイオン源の性能を決定する因子になると考えられる H−の生成率を計算するため，Oen と Robinson［７］が作成た．これは１９７４年にBelchenkoら［２］が小型イオン源にセシした MARLOWE コードを用いることにより，反射粒子の2

ウム（Cs）を添加して 3 A/cm を超える電流密度で水素負速度分布関数を計算した．１９８０年代当時の計算機能力ではイオン（H−）が引き出せ得ることを実証し，引出電極に対時間のかかる計算であり，使用された粒子数はわずか１０６向した陰極表面での高効率負イオン生成の可能性を指摘しであった．図１に Hiskes の論文に掲載された表面反射水素たためである．

原子のエネルギー分布を示す．表面生成負イオンのエネル

当初からCsの役割は電極表面の仕事関数を低下させ，電ギー分布実測データと直接比較するには精度的に問題が極金属から水素原子への電子移行率を増大させることでああったが，水素負イオン収率をエネルギーの関数として理ろうと予想された．Kishinevskii［３］は単純な表面ポテン

論計算する上で，強力なツールとなり実験データを良く説

シャルモデルにもとづいて表面での水素原子負イオン化率

明した．その後核融合プラズマ研究のプラズマと壁相互作を計算し，表面仕事関数の低下が負イオン化率を指数関数用分野で粒子反射のみならずスパッタリングを計算するた的に増大させることを示した．Belchenko らの成功に触発

め，TRIM/SRIM［８］や ACAT［９］な ど の 二 体 衝 突 近 似され，BNL（Brookhaven National Laboratory）
［３］
や LBL

（BCA: Binary Collision Approximation）シミュレーション（Lawrence Berkeley Laboratory）［４］で表面生成型の負イ

コードが利用されるようになった．

オン源が開発されるとともに，水素負イオン生成電極表面

表面生成型負イオン源の，コンバーターと呼ばれる H−の素過程に関する研究も開始された．同時にプラズマ―固生成電極表面からの負イオン放出シミュレーションは，コ

体壁相互作用の研究のために開発されたシミュレーション

ンバータータイプの H−源が LANSCE（Los Alamos Neu-

モデルも，表面で負イオン化される水素原子の速度分布を

tron Science Center）でしか使用されなくなったこともあ

Doshisha University, Kyotanabe, KYOTO 610-0321, Japancorresponding author’s e-mail: mwada@mail.doshisha.ac.jp

５８１!2021 The Japan Society of Plasma
Science and Nuclear Fusion Research

Journal of Plasma and Fusion Research Vol.97, No.10 October 2021

図２

図１

Hiskes らが水素表面負イオン化率計算に使用した MARLOWE により求めた Na（点線）と Cs（実線）に対する反射
水素の角度分布（左）とエネルギー分布（右）．粒子数は106 個．［６］よ り 引 用．Reproduced from J.R. Hiskes and
P.J. Schneider, Phys. Rev. B 23, 949 (1981), with the permission of American Physical Society.

水素プラズマ中 200 V にバイアスされた Mo 表面から放出される負イオンのエネ ル ギ ー ス ペ ク ト ル．（a）
ACAT にRasser［１３］
らによって計算された表面負イオン生成率を組み合わせたシミュレーション結果．（b）
対象とした実験条
件での測定結果．（b）
は［１５］より引用．Reproduced from
M. Wada, LBL-Report, LBL-15800, Copyright (1983) U.S.
DOE.

ただし，ビームとして引き出された電流にはイオン衝撃脱離成分が大きく現れるものの，表面から放出される負イオンの角度分布を積分すると，負イオン生成率の高い反射成

り［１０］，表面生成型負イオン源から生成される H−イオン分の方が脱離成分に比べて大きいことが実験的に確認されのエネルギー分布・角度分布のシミュレーション結果はており［１８］，理論予測と良く一致する．

Hiskes の報告から３０年後に，漸く Kenmotsu ら［１１］によっ

４．
２．
２ プラズマグリッドでの水素負イオン表面生成て高精度化されたものが報告されることになる．シミュ

Okumura ら［１９］は，磁気フィルターを用いた体積生成

レーション粒子数は必要精度に応じて１０７から１０８にまで増型負イオン源にCsを導入して使用すると，通常の体積生成加している．負イオンのエネルギー分布を詳細に計算するの数倍の H−電流が得られるのと同時に，H−とともに引きには，特に表面負イオン化率のエネルギー依存性を明らか出される電子電流量が低減されることを，大型負イオン源にしておく必要がある．Seidl［１２］らは半理論，半実験的にを用いて実証した．プラズマグリッ ド（PG）と呼ばれ最小仕事関数付近の金属表面負イオン化率を与える簡略モる，イオン源プラズマと対向する引出孔の開けられた電極

デルを提案しており，詳細な理論解析を行った Rasser ら

温度に対応して引出 H−電流量が変化したことから，PG［１３］のモデルと共に良く用いられる．さらに Rasser らに

表面の仕事関数低下に伴い，水素負イオンは PG 表面で生

よって紹介された Blandin ら［１４］の低脱出速度近似式も，成されていると考えられた．実際に Leung ら［２０］は Ba

表面生成される一般の負イオンのエネルギー・角度分布をを添加した小型のイオン源を用いた実験で，プラズマ電極

計算するモデルに組み込まれることが多い．

電位で生成されたと考えられる H−イオン成分を観測して粒子反射模型と負イオン化率を組み合わせることによっいる．

て，表面生成 H−のエネルギー分布を計算した結果の一例

水素負イオンが壁から連続的に供給される状態になるを図２に示す．比較のために同一条件で採取した実測デーと，壁の近傍では正イオンと負イオンにより，電子の少な

タ［１５］を合わせて示す．実測データはCs導入前の高仕事関

いイオン―イオンプラズマが形成されることになる．水素数表面で生成される H−電流のエネルギースペクトルを示を吸着した低仕事関数表面に水素正イオン・原子が入射す

＋

＋
，H＋
2 ，H3 の表面反射により生成さる際の表面生成 H−のエネルギー分布，及び角度分布がれた H−の成分が分離されて，シミュレーション計算の有BCA に表面負イオン化率を組み合わせたモデルで計算さすものであり，入射 H

効性が確認できる．なお図２の結果では，表面に吸着されれている．その結果，表面に斜入射する水素正イオン・原

た水素が入射イオンに弾き出されてH−として放出される，子に対して負イオン転換率が高く，吸着水素がイオン衝撃

イオン衝撃脱離成分がグラフ上に現れるレベルをはる

脱離過程により H−として放出される量が，図３に示すよ

かに下回っている．同様の実測結果と SRIM シミュレー

うに低エネルギー部分で大きいことが示された［２１］．計算

ション結果との比較は，Cartry［１６］らによって黒鉛表面にされた表面生成 H−のエネルギー，角度分布によれば H−の対して行われている．黒鉛の場合は表面に吸着された水素

表面脱出速度の垂直方向成分は小さく，PG 表面近傍で H−がイオン衝撃脱離によって H−として放出される成分が大

密度が局所的に増加する．PG 表面からは水素放電からのきいことが実験的に観測されている．表面生成型水素負イ

光子により光電子が放出されることから［２２］，実際の PG

オン源においても仕事関数極小条件においては，表面吸着

近傍の H−対電子比は電子偏向磁場構造などにより大きな

水素のイオン衝撃により生成される H−が負イオン電流の影響を受けるものと考えられる．イオン源から引き出され大部分を占めるエネルギースペクトルが観測される［１７］．

る H−電流対電子電流比を低減するなどの性能最適化のた
５８２

Special Topic Article

4. Plasma-Wall Interactions in Ion Sources

M. Wada and T. Kenmotsu

イオンに対する PG 表面の Cs 損失率を計算した結果を示したものである．図４と比較して，Cs 損失の立ち上がるエネ
ルギーが低く，損失率が大きいことがわかる．したがって実際のイオン源運転において，真空リークなどがあると Cs
消費量が増大するものと予想される．

４．
３ 詳細モデルの研究
４．
３．
１ 低エネルギー領域への対応
プラズマ―固体壁相互作用に詳しい読者の中には図３のエネルギー範囲を見て，二体衝突近似を用いるには無理があることに気づかれた方も多いと思う．ところが，結果が
図３

10 eV の水素原子入射により放出される水素原子（破線）と水素負イオン（実線）のエネルギー分布．［２１］より引用．
Reproduced from M. Wada et al., Rev. Sci. Instrum. 85,
02B114 (2014), with the permission of American Instituteof Physics.

妥当かどうかを検証しようにもこのエネルギー範囲の表面
素過程データを実験的に採取することは困難である．そこ
で核融合装置のプラズマ―壁面相互作用の研究グループで行われている MD（Molecular Dynamics：分子動力学）計

め，今後さらに詳細な PG 表面負電荷生成モデルが構築・算を行い，数 eV の H−の表面生成についてのシミュレー

応用されるものと期待される．

ション結果を Rutgliano ら［２７］が報告している．ただし，計４．
２．
３ Cs リサイクリング計算

算に時間がかかることから，必ずしも統計的に十分な精度

２０２１年時点での ITER 実験計画では，中性粒子入射加熱までデータを蓄積できていない状況である．また，基本的用水素負イオン源に Cs を導入して運転することが計画さに入射イオン（原子）の反射としての負イオン化であり，れており，イオン源には制御性の高いオーブンからCsが供吸着水素の脱離までを計算するには，さらに計算対象の粒給され［２３］，イオン源運転の最適化が図られる．本来必要子数が増えることから，すぐには結果が得られないと予想となる Cs 量は PG 表面に半単原子層を形成するだけの量でされる．それでも Damone ら［２８］による２
０１５年のポテン

2良く，最小仕事関数状態に対して 0.5 mg/m 程度しか必要でないことになるが，実際のイオン源の運転ではCsを導入
後負イオン電流が減少するので，Cs を追加導入することによって水素負イオン電流を維持している．フィラメントカ
ソードを用いたアーク型イオン源の場合はカソードから蒸散されたカソード金属が Cs を覆い，PG 上 Cs が失われると考えられているが，同様の Cs 追加導入はカソードが無い
RF 負イオン源でも必要であり，また重水素で特に必要なCs 導入量が多くなると報告されている［２４］．
そこで水素放電に対向する PG 表面から，吸着 Cs がどのように失われるかについて，Cs 層‐水素（重水素）層‐Mo
バルクの構造を作って BCA を用いた計算を行った．Cs

図４

Mo 製 PG 上に吸着した Cs のスパッタリング率．重水素に対して Cs の損失率が大きくなることがわかる．Reproduced from M. Wada, Rev. Sci. Instrum. 89, 052103 (2018),
with the permission of American Institute of Physics.

図５

Mo 製 PG 上に吸着した酸素正イオンによる Cs のスパッタ
リング率．軽水素吸着層と比較して重水素吸着層の方が大きな Cs 損失率を示す．Reproduced from M. Wada, Rev.
Sci. Instrum. 89, 052103 (2018), with the permission ofAmerican Institute of Physics.

と水素の同時吸着表面の構造とその仕事関数については比較的古くから研究されており［２５］，報告されている最小仕事関数状態に近い表面構造：Cs が最上部に吸着し，Cs と母材金属原子の間に水素同位体の層が形成される構造に対して計算を行った．図４に結果の一例を示す［２６］．図４の水素―重水素の比較で，重水素運転時には PG 表面の Cs
が，より失われやすいことがわかる．これは質量効果によって衝突カスケードの最終段階で表面 Cs を叩き出す際の運動量が大きいためだけでなく，バルク層とCsの間にある水素同位体が運動量ベクトルの向きを変えるのに寄与しているためと考えられる．同様の同位体効果は，酸素不純物に対して行った計算結果に対しても見られた（図５）．
中性粒子加熱用の負イオン源の運転では真空度も良く，大きな量の不純物イオンがプラズマ中にあるとは予想できないものの，逆流正イオンなどによってイオン源構成要素が加熱するなどした際に，少量の不純物イオンがイオン源
プラズマ内に生成される可能性がある．図５は酸素不純物
５８３

Journal of Plasma and Fusion Research Vol.97, No.10 October 2021

シャルの提案と Cs-Mo 表面構造の計算以降，着実に計算をきた．二体衝突近似が使える入射粒子エネルギー領域に対進めており，今後信頼性の高い表面生成 H−のエネルギー

しては表面生成 H−スペクトルや角度分布として妥当な結分布と放出角度分布が得られるものと期待される．

果が得られ，実験データも得やすいことから実測結果と直４．
３．
２ 表面電子状態密度計算

接比較可能である．例としてCsを用いない負イオン源運転

表面の電子状態密度分布は，本質的に３次元であって水を可能とすると期待され，研究が推進されている C12A7

素原子への電荷移行を詳細に計算するためには表面構造の（12CaO‐7Al2O3）エレクトライド材料を水素プラズマ中に正確な情報が必要となる．Padama ら［２９］はCs-W系に水素保持して−80 V に保持した場合に表面から生成される H−

原子が吸着すると Cs の吸着エネルギーがどのように変化のうち，前方に放出されるエネルギースペクトルを計算しするかを計算し，水素吸着が全体としてのCsの吸着エネル

た結果を図６に示す．この計算においても，低エネルギー

ギーを下げることを明らかにしている．Cs の表面吸着エネ側においては表面に吸着されている水素，エレクトライド

ルギーは PG 表面からの Cs 損失率を計算する上で重要なパの籠構造にある水素が叩き出されて H−となる成分が含まラメータであり［２６］，水素負イオン源内の Cs リサイクリれるため，エネルギー分布・角度分布の精度が疑問視されングをモデル化するのに不可欠なデータである．また水素る．それでも図６の結果は，実験で報告されたスペクトルプラズマの原子・分子イオン・中性粒子照射を受けた Cs［３１］に近い形状を示している．

吸着表面がどのように自己組織化し，粒子バランスが形成

低エネルギー領域におけるエネルギー分布，および角度されるのかを知るためには，大規模な表面反応シミュレー

分布シミュレーションの研究は，今後も高精度化を目指し

ションを行うしか方法がない．階層に分けて得られた結果た進展が続くものと予想されるが，現在の計算精度でも負を数値データ化して，負イオン源プラズマ全体のシミュ

イオン源シミュレーションモデル全体の信頼性に大きな影

レーションモデルに組み込む作業が，今後必要となるもの

響を与えないのでは無いかと考えられる．そもそも PG 表と予想される．

面を離れるときの H−運動エネルギーは小さく，また Csの PG からの損失は高エネルギー粒子や不純物の寄与が大４．
４ イオン源寿命評価シミュレーション

きいので，こちらも低エネルギー側モデルの高精度化に

４．
４．
１ グリッド損耗

よって大きな誤差は生じないと予想される．ではイオン源

ITER の建設が進められる中，核融合発電商用化実現のプラズマ，シース領域，引出領域のモデルが完成しつつあための次ステージとしてデモ炉の仕様検討が進められている状況下で，プラズマ―壁面相互作用のパートが為すべきる．デモ炉のシステム検討にあたっては機器寿命の評価がことは無いのかと言われると，全くの逆で，多くのことに不可欠であり，装置壁面の損耗率などが候補材料に対して取り組む必要があると考えている．筆者らが負イオン源プ

再評価されることになるものと考えられる．イオン源の寿ラズマと固体壁面の相互作用研究を開始したころ，国際会

命と言えば，宇宙探査衛星に使用されるイオンエンジンの議の場で ITER 用 RF 負イオン源の研究開発を指揮する

グリッド損耗［３０］が連想されようが，同様の損耗過程がデFantz 氏に「表面シミュレーションでは綺麗な表面を対象

モ炉の場合には問題となり得る．特に中性粒子加熱の場合としているが，実際の負イオン源の中の Cs はもっと汚は運転期間を通じてビーム集束性を担保する必要があるのい．」と指摘を受けた．我々も基礎過程の実験で Cs と格闘

で，イオンエンジンと同等レベルでの寿命評価シミュレー

することが多く，確かにイオン源内に理想的な表面状況がションが行われる必要があろう．

実現されているとは考えていない．Yoshida ら［３２］の行う

４．
４．
２ 逆流イオンによる損傷

負イオン源内 Cs 輸送計算に利用可能な精度の基礎データ

光中性化セルなどの使用により，DEMO 炉設計では中性として供出するため，表面素過程実験と実イオン源実験結

粒子入射用負イオン源の運転ガス圧は低くなるものと予想される．しかし，それでも下流側での正イオン発生を０にすることは難しいので，逆流正イオンがイオン源の PG 対
向側端板に衝突することが予想される．1 MeV 以上の重水素であるので，逆流正イオン損傷防止用に端板に取り付けられる予定の Mo アーマータイルにインプランテーション（イオン注入）されることになるため，大きなスパッタリン
グが生じるとは予想されない．しかしながら，水素吸蔵が進んで材料が脆化することも危惧されるので，内部応力の発生に伴う剥がれ等の損傷可能性などについて，シミュレーション計算にもとづいて十分に検討しておく必要があ
ろう．
図６

４．
５ おわりに計算機能力の飛躍的進歩により，以前は計算不可能と考えられていた表面素過程がシミュレーション可能となって
５８４

C12A7 エレクトライド表面から放出される水素負イオンのエネルギースペクトルをシミュレーションで計算した結果．角度分布を考慮すると反射成分の寄与が小さくなり，観測されたスペクトル［３１］を再現する．

Special Topic Article

4. Plasma-Wall Interactions in Ion Sources

M. Wada and T. Kenmotsu［１５］M. Wada, LBL Report, LBL-15800 (1983).
［１６］G. Cartry et al., Phys. Plasma 19, 063503 (2012).
［１７］M. Wada et al., J. Appl. Phys. 67, 6334 (1990).
［１８］M. Wada et al., AIP Conf. Proc. 111, 247 (1984).
［１９］Y. Okumura et al., AIP Conf. Proc. 210, 169 (1990).
［２０］K.N. Leung et al., Phys. Rev. Lett. 62, 764 (1989).
［２１］M. Wada et al., Rev. Sci. Instrum. 85, 02B114 (2014).
［２２］J. Lulainen et al., Phys. Plasma 24, 103502 (2017).
［２３］M. De Muri et al., Fusion Eng. Des. 167, 112331 (2021).
［２４］B. Heinemann et al., New J. Phys. 19, 015001 (2017).
［２５］C.A. Papageorgpoulus and J.M. Chen, Surf. Sci. 39, 283(1973).
［２６］M. Wada, Rev. Sci. Instrum. 89, 052103 (2018).
［２７］M. Rutigliano et al., Surf. Sci. 693, 121545 (2020).
［２８］A. Damone et al., Plasma Phys. Control. Fusion 57, 035005(2015).
［２９］A.A. Pdama et al., e-J. Surf. Sci. Nanotechnol. 16, 391 (2018).
［３０］I. Funaki et al., 47th AIAA/ASME/SAE/ASEE Joint
Propulsion Conf. AIAA 2011-5734 (2011).
［３１］M. Sasao et al., Appl. Phys. Express 11, 066201 (2018).
［３２］M. Yoshida et al., Rev. Sci. Instrum. 91, 033313 (2020).

果をもとにした検証を，今後もさらに進めていく必要がある．

参考文献［１］C.F. Chan et al., J. Appl. Phys. 54, 6119 (1983).
［２］Yu. I. Belchenko et al . Nucl. Fusion 14, 113 (1974).
［３］K. Prelec, Nucl. Instum. Methods 144, 413 (1977).
［４］K.W. Ehlers and K.N. Leung, Rev. Sci. Instrum. 51, 721(1980).
［５］P.J. Schneider et al., Phys. Rev. B 23, 941, (1981).
［６］J.R. Hiskes and P.J. Schneider, Phys. Rev. B 23, 949 (1981).
［７］O.S. Oen and M.T. Robinson, J. Nucl. Mater. 76-77, 370(1978).
［８］J.F. Zeagler, Nucl. Instrum. Methods B 136-138, 141 (1998)［９］Y. Yamamura, Nucl. Instrum. Methods B 28, 17 (1987).
［１０］R. Keller et al., AIP Conf. Proc. 1097, 161 (2009).
［１１］T. Kenmotsu et al ., AIP Conf. Proc. 1390, 134 (2011).
［１２］M. Seidl et al., J. Appl. Phys. 79, 2896 (1996).
［１３］Rasser et al., Surf. Sci. 118, 697 (1982).
［１４］A. Blandin et al., J. Physique 37, 369 (1976).

５８５J. Plasma Fusion Res. Vol.97, No.10 (2021)５８６‐５９１

小特集 イオン源の数値シミュレーション技術の進展と応用

５．引出領域プラズマと負イオンビーム引出のモデリング
5. Modeling of Extraction Region Plasma and Negative Ion Beam Extraction
宮本賢治
MIYAMOTO Kenji
鳴門教育大学（原稿受付：２
０
２
１年６月１
７日）

負イオン源の引出部における電位構造や負イオンの引出・輸送過程のモデリング・解析の多くは，PIC シ
ミュレーションにモンテカルロ計算が付加された PIC-MCC（Particle in Cell with Monte Carlo Collision）を用いている．本章は PIC-MCC シミュレーションによる，プラズマメニスカス（負イオン引出界面；以下，メニスカス）からの負イオンの引出過程とビーム光学特性に関する計算結果例を紹介する．また表面生成負イオンが存在する
場合に，イオン性プラズマが形成される物理機構についての計算結果も紹介する．
Keywords:
plasma, negative ion, simulation, extraction, PIC, Monte Carlo, beam optics, Coulomb collision, plasma meniscus（１）式において，'
，#
*10
'%，$
"! はそれぞれプラズマ電極

５．
１ 引出領域プラズマのモデリング
本負イオン源の引出領域では，フィルター磁場や電子抑（PG）の引出孔の半径，引出電圧，負イオン電流密度を表

制磁場のような強い横磁場や，引出静電場の浸透，プラズ

す．また下付きの添え字 sim と real はそれぞれスケーリン

マ自身の空間電荷効果，正・負イオンと中性粒子との衝突

グしたシミュレーション上の値と実際の値である．スケー

過程など多くの物理過程が混在するプラズマを解析する必リング率を (とすると，#
'
.
#
/
+
- "(
'%$
'%$
.
*(
,となり，これを要がある．現在，負イオン源の引出領域の数値計算モデル（１）式に代入すると，#
$*10$
/
+
- が (を用いて表される．加の多くは，PIC（Particle-In-Cell）法，特に explicit PIC 法

速電圧 $())については，加速電圧／引出電圧の比が，ス［１‐３］を用いている．

ケーリングしたシミュレーション上の値と実際の値とが同

５．
１．
１ スケーリング

一にするように設定する．すなわち，

explicit PIC 法を用いる場合，安定な数値解を得るため

$*10$
$())$
$())$
$*10$
#
#
#
/
+
-$
/
+
- "#
.
*(
,$
.
*(
,

に，メッシュ幅はデバイ長以下にする制約がある
［２］．負

（２）

イオン源の典型的な電子密度はオーダーで 1017〜1018 m−3

次に磁場の値 ! については，ラーマー半径 '& についてなので，デバイ長は μm のオーダーとなる．一方，負イオン

'
'
#
#
&$
&$
/
+
- "(
.
*(
,となるので，次式が成り立つ．

源の引出・加速部のサイズは数 10 cm もあり，実サイズで!/+- "

モデリングすると莫大なメッシュ数を要することになる．
そのため計算資源の制約から，PIC 法を用いた数値シミュ)/+!
().*(, .*(,

（３）

レーションではサイズがスケーリングされている．またこ（３）式には速度 )が含まれているが，これはエネルギー保れに伴って，引出・加速部での電圧や磁場の値もスケーリ

存則ングする必要がある［４］．

!%)"
$/+-，
"&
" /+-

電圧や磁場の値のスケーリングの方法については，以下の通りである．まず電圧については，負イオンの電流密度

!%)"
$.*(,
"&
" .*(,（４）とパービアンスが，スケーリングしたシミュレーション上を用いて消去することができる．
（４）式で $/+-，$.*(,をの値と実際の値とが同一にするように設定する．パービア（１）式と連立すると
$*10$
$*10$
#
/
+
-，#
.
*(
,に置き換えてンスとは，空間電荷制限電流と引出電圧の 3/2 乗の比であ!/+- "(!!!#!.*(,

り，荷電粒子ビームの光学的な性能に対する指標の１つである［５］．すなわち次式が成り立つ．
"
"
#
"#
'
"#
'
"!
'%$
"!
'%$
/
+
-#
.
*(
,
"
#
#
"
"
$*10$
$*10$
#
#
/
+
.
*(
,

（５）

が得られる．
サイズ以外にも電子密度や真空中の誘電率をスケーリン（１）

グする方法がある［６］．いずれも同等であることが示されている．これまでに，スケーリングを用いた数値シミュ

Naruto University of Education, Naruto, TOKUSHIMA 772-8502, Japanauthor’s e-mail: kmiyamot@naruto-u.ac.jp

５８６!2021 The Japan Society of Plasma
Science and Nuclear Fusion Research

Special Topic Article

5. Modeling of Extraction Region Plasma and Negative Ion Beam Extraction

K. Miyamoto

荷電交換衝突等の衝突断面積の総和をレーション結果は実験結果と定性的に傾向が一致することが確認されており［７，８］，負イオン源の引出領域における,' #,(",)",*"+

電位構造解析や負イオンの引出の物理機構の解析へ十分に（９）適用可能である．

とする．ただし，,(，,)，,*+ はそれぞれ衝突 (，衝突 )，５．
１．
２ PIC-MCC について

+ の断面積とする．衝突断面積の総和 ,' の内，最大衝突 *となる衝突断面積を ,'!-(4とし，下の式で定義される偽のクーロン衝突や荷電粒子と中性粒子の衝突過程について

衝突断面積 ,.3,,を導入する．

は，モンテカルロ法を用いて計算する．このように PIC
シミュレーションにモンテカルロ計算が付加されたシミュ,.3,,#,'!-(4!,'レーションは PIC-MCC と呼ばれる．

（１０）

この ,.3,,を用いて，新たな衝突断面積の総和を次式で定義

５．
１．
２．
１ 衝突過程のモデリング
まず，クーロン衝突においては２体衝突モデル（BCM:する．

Binary Collision Model）［９］が用いられる．BCM では，セル*,'#,'",.3,,

中における２つの荷電粒子をランダムに選択し，運動量保（１１）存則，エネルギー保存則から衝突後の速度を求める．クー*は負イオンのエネルギーに依存しない．"$"!#の範囲,'

ロン衝突前後における速度変化は（６a）式〜（６d）式で表さ, ,,の一様乱数を発生させて，もし "! .3* ならば，衝突しな,',.3,,",(, ,,な ら ば，衝 突 (かったとする．も し .3
* $"!
*
,'
,'
, ",",, ",が発生したとする．また，もし .3,,* ( $"! .3,, *( )
,'
,'

れる．
#(4#

(5
(4(
+
.+*
/1$! (61
+
.+1
+
.$
61
(%
(%
（６a）
#!*
/1+'
!(4&

ならば，衝突 b が発生し，以降も同様である．

(5(#(5# (61+.+*/1$! 4(61
+
.+1
+
.$
(%
(%

なお前節１．
２で述べたスケーリングを適用する場合，衝)を実サイズと同じ値にな
突過程についても衝突周波数 $,

（６b）
#!*
/1+'
!(5&
#(6#!(% 1
#!*
/1+'
+
.+*
/1$!(6&（６c）$$(% #!(4"(5（６d）

るようにする必要がある．この場合，衝突断面積 ,や時間
ステップ幅 #'（)#!"#'）をスケーリングすることが挙げ
られる．
５．
１．
２．
２ 負イオンのエネルギー緩和過程のシミュレー+，$は図１に示す散乱角で，"$(#，($!#の範囲の一様

ション

乱数を用いて次式で与えられる．
$#$*(#， +#$(
0
*
2
(
.($

クーロン衝突や荷電交換衝突による表面生成負イオンのエネルギー緩和過程が負イオンビームの光学特性に及ぼす（７）

影響について検討した［１２］．負イオン源の引出領域をモデ

ただし一様乱数 ($の平均値と分散は正規分布を仮定した

ル化した．プラズマ電極表面から放出される表面生成負イ場合，クーロン衝突の理論から次式を満たす必要がある．

オンの初速度分布として，1.5 eV の初期温度の半マクス

$
$&
%&% '
,
.%
$
($)
(
#"， (
($
)
# ' $' % % #'
*)
#(
&
" &

ウェル分布を仮定した．これは，水素のフランクコンドン（８）原子の入射エネルギーを 2〜3 eV とし，プラズマ・壁相互作用計算コー ド（ACAT: Atomic Collision in Amorphous$'，%' は場の粒子の密度と質量を，%&はテスト粒子の質量

Target，第４章参照）での計算結果を外挿すると，エネルである．また #0は場の粒子とテスト粒子の換算質量であ

ギー反射係数が 0.5 程度［１３］と見積られることによる．まる．

た，正イオンと水素原子の温度はそれぞれ0.25，0.1 eVとし

次に荷電粒子と中性粒子の衝突過程については，ヌル衝た．

突法［１０，
１１］を用いる．ヌル衝突法の概要は以下の通りで

図３は負イオンのエネルギー分布の時間発展である．時ある．例えば負イオンと中性粒子の弾性衝突や励起衝突，間の経過と共に，負イオンエネルギー分布のピークが低い方へ移動しており，クーロン衝突により負イオンのエネル
ギーが緩和されることがわかる．
引出部における引出電極（Extraction Grid: EXG）出口と電子抑制電極（Electron Suppression Grid: ESG）出口での負イオンビームのエミッタンス図を，エネルギー緩和過程
が無／有についてそれぞれ 図４，図５に 示 し た（EXG，ESGについては図２を参照）．一般的にエミッタンス図は，楕円の形状をし，この面積（エミッタンス）が小さいほど，
負イオンビームの光学特性は良好である．エネルギー緩和過程がある場合の方が，エミッタンスは小さい．よって，

図１

クーロン衝突後の散乱角．

エネルギー緩和過程により，負イオンビームの光学特性が
５８７

Journal of Plasma and Fusion Research Vol.97, No.10 October 2021

改善することが分かる．

［１４］．このイオン性プラズマの理解は，負イオンビームの大電流化やメニスカス構造の理解と負イオンビームの光学

５．
２ イオン性プラズマの形成要因

特性の向上を図る上で，非常に重要である．

多量の表面生成負イオンが存在する場合，プラズマ電極

PIC-MCC シミュレーションにより，イオン性プラズマ

近傍に正・負イオンのみから成り，電子がほとんど存在しの検証を行った［１５］．負イオンの引出し方向を x 軸に，まないペアプラズマ（イオン性プラズマ）が観測されている

たフィルター磁場の方向を &軸に取った．
電子はクーロン衝突や中性粒子との衝突によって，（１２）式で表されるステップ幅で，フィルター磁場を横切って拡散するようにモデル化した．
$# #( ("%%" $$#%" " !%#$
'& %

（１２）

(%，#$
(
ただし !
，%'はそれぞれ拡散係数，時間ステップ
幅，０〜１の一様乱数である．（１２）式は拡散方程式から求められ，$&はフィルター磁場に沿って電子が壁で損失するに要する時間を，また $% は電子がフィルター磁場を横切って拡散するのに要する時間に相当する．すなわち !$"$%&は，電子の実効的な閉じ込め性能を示すパラメータで，この値が大きいほど電子の閉じ込め時間が長い．
なお，５．
１．
２．
１節で述べたような２体衝突モデルやヌル
図２

衝突法を考慮すると，自然に荷電粒子とのクーロン衝突や負イオンのエネルギー緩和過程のシミュレーションにおける引出領域のモデル概要図［１２］．

中性粒子との衝突の効果が電子軌道に反映されるが，ここでは，パラメータの !$"$% がイオン性プラズマの形成への&寄与について検討するために，（１１）式で表される「!$"$%&モデル」を用いた．
図６は（a）!$
!
!
"
$
"
$
!$，（b）!$!に 対 す る，PG
% #!
% #"
&
&

#&分布である．図６で，%#9 mm の位置が近傍での #%!"#&の値が４０以上の領域がイオン性プ
PG 入口である．#%!"
ラズマ層に相当する．すなわち，電子の閉じ込め時間が短
いと電子密度は低下するが，代わりに準中性条件を保持するように負イオン密度が高くなり，イオン性プラズマ層が

図３

負イオンのエネルギー分布の時間発展［１
２］
．

図４

引出部での負イオンビームのエミッタンス図（エネルギー緩和過程を考慮しない場合）
［１
２］
．

図５

引出部での負イオンビームのエミッタンス図（エネルギー緩和過程を考慮した場合）
．

図６

５８８

（b）#$
PG 近傍での nH− / ne 分布．（a）#$‖/ $⊥ ＝ 0.08，
‖/ $
⊥
＝ 1.0
［１
５］
．

Special Topic Article

5. Modeling of Extraction Region Plasma and Negative Ion Beam Extraction

K. Miyamoto

なる．このように，メニスカス形状は引出電圧に依存する．

形成されやすいことが言える．

一般的に荷電粒子ビームの光学特性と，Child Langmuir

イオン性プラズマの生成機構については他に，アーク
チャンバーとプラズマ電極の間に印加するバイアス電圧に則から求めた空間電荷制限電流と引出電圧の 3/2 乗の比よる電子の選択的な引き込みが原因であるとも考えられ（（１３）式参照）であるパービアンス値 "との間には相関性

る．表面波プローブを用いた PG 近傍の電子密度測定によがある．

り，体積生成負イオンのみの場合，電子密度はバイアス電

"! !
$!
#""!

圧に大きく依存するが，多量の表面生成負イオンが存在す
る場合はあまり依存しないことが示されている［１６］．この（１３）

実験結果からも，特に表面生成負イオンが存在する場合，

ただし，$はメニスカスと引出電極との距離である．パー電子の閉じ込め時間がイオン性プラズマの形成に寄与する

ビアンス値 "の引出電圧依存性を図１０に示す．負イオンことが示唆される．

ビームの発散角が負イオンビームの発散角が最小となる電圧 10 kV 付近で，パービアンス値は最小になり，これが最５．
３ 負イオンビーム光学特性の引出電圧依存性

適なパービアンス値に相当する．この条件よりも引出電圧

一般的に安定かつ効率よく負イオンを加速・輸送するたを低くすると，（１３）式の分母の ##!
%$が小さくなり，パービ

めに，発散角の小さい収束性の良い負イオンビームが要求

アンス値は増大する．一方，この条件よりも引出電圧を高される．そのため，PIC-MCC シミュレーションによる負イ

くすると，メニスカスと引出電極との距離が大きくなるの

オンビームの光学特性の研究が行われている．

で，（１３）式の分子の$!が大きくなり，パービアンス値は増大する．

負イオンビームの光学特性の引出電圧依存性について検討した［１７］．CERN 研究所の Linac４用負イオン源［１８］を最適パービアンス値においてメニスカス形状はほぼ平坦対象とした．負イオン源の断面図と数値計算領域を図７にであるが，PG 孔端部では湾曲している．メニスカス端部か示す．負イオンの引出し方向を x 軸に，フィルター磁場にら引出された負イオンは，この湾曲によりビームハロにな

平行な方向を y 軸とした．プラズマ電極 PG（ただし図７でりやすいことが，これまでの先行研究から明らかにされてはPEと表記）は，表面生成された負イオンを効率よく引きいる．図１１に負イオンビームの電流密度と発散角の y 方向出すためにテーパー形状となっている．また，ここでは体の空間分布を示す．図１１（b）において，負イオンビームの積生成負イオンのみを考慮した．

コアの部分（−1.5 mm!%!0.8 mm）は発散角が小さいに引出電圧を変えた場合の負イオンビームの発散角とメニもかかわらず，ハロの部分は（%!−1.5 mm，0.8 mm!%）

スカス形状を図８，図９にそれぞれに示す．引出電圧がは，発散角は大きい．PG 孔端部でのメニスカスの曲率は，10.7 kV において発散角は最小になる．これは実験結果と図７に示す，引出部側テーパーの開口角度に依存する．す一致する．この場合，メニスカス形状はほぼ平坦である．

なわち，この開口角度が大きい程，PG 孔端部でのメニスカ

引出電圧が低い場合（7 kV）には，メニスカス形状は引出

スの曲率は大きくなる．したがって，負イオンビームハロ側へ突出した凸形状に対し，引出電圧が高い場合（14 kV）を抑制するには，なるべくこの開口角度を小さくする必要には，メニスカス形状はプラズマ側へ入り込んだ凹形状に

図７

図８

図９

Linac４用負イオンの断面図と数値計算領域［１７］．

引出電圧を変えた場合の負イオンビーム発散角［１７］．

引出電圧を変えた場合のメニスカス形状［１７］．

図１０ パービアンス値 P の引出電圧依存性［１
７］
．

５８９

Journal of Plasma and Fusion Research Vol.97, No.10 October 2021

図１
２ 電子抑制磁場に垂直な平面での（a）
電子密度と（b）
負イオン
密度の空間分布［１９］．

図１３ 電子抑制磁場に垂直な平面でのプラズマ中へ浸透する電場の空間分布［１９］．

プラズマメニスカス付近における負イオン密度の空間的な非対称については，エミッタンス測定結果を基に，測定

図１
１ 負イオンビームの電流密度と発散角の y 方向の空間分布［１７］．

位置（加速部出口から約 1 m 下流）から負イオンプラズマ
側へ逆方向に負イオンビームの軌道計算した結果［２０］ともがあると考えられる．

傾向が一致している．

５．
４ 表面生成負イオンが存在する場合のメニス
カスと負イオンビームの光学特性

５．
５ おわりに引出領域に表面生成負イオンが多量に存在するイオン性

前章で体積生成による負イオンビーム軌道は，メニスカ

プラズマにおいて，PIC-MCC シミュレーションは負イオ

ス形状に左右されることを示した．これは表面生成負イオ

ンプラズマ中の電位構造やメニスカスを自己矛盾なく求めンの場合についても同様である．

ることができるため，負イオンの引出・輸送過程の解析に

表面生成負イオンが存在する場合の電子抑制磁場に垂直非常に有効である．特にメニスカスについては，従来の正な平面での電子密度と負イオン密度の空間分布を図１２に示

イオンと電子から成るプラズマにおける Child-Langmuir

す［１８］．負 イ オ ン 引 出 用 の 電 界 と 電 子 抑制 磁 場 に よ る則が必ずしも成立するとは限らないため，これに基づいて"!!ドリフトの影響で，PG 近傍ではフィルター磁場と垂負イオン源でのメニスカスを求めることは正しくない．

直な平面において，電子の空間分布が非対称性となる．こ負イオンビームの光学特性はメニスカス形状に依存すれに伴って図１２（b）に示すプラズマメニスカスの形状，まる．したがって引出領域に加速部を追加したモデリングにたはプラズマ中へ浸透する電場（図１３参照）は非対称になより負イオンビームの軌道計算を行うことで，PIC-MCC

る．したがって，負イオン密度の空間分布も非対称になる．

シミュレーションは負イオンビームの光学特性の解明への
５９０

Special Topic Article

5. Modeling of Extraction Region Plasma and Negative Ion Beam Extraction

K. Miyamoto［１０］V. Vahedi, M. Surendra, Comput. Phys. Commun. 87, 179(1995).
［１１］小越澄雄：プラズマ工学（電気書院，2010，p.131）．
［１２］K. Miyamoto et al., AIP Conference Proc. 1869, 020002(2017).
［１３］M. Wada et al., AIP Conference Proc. 1515, 59 (2013).
［１４］K. Tsumori et al., Rev. Sci. Instrum. 83, 02B127 (2012).
［１５］S. Nishioka et al., J. Appl. Phys. 119, 023302 (2016).
［１６］M. Kisaki et al., Rev. Sci. Instrum. 85, 02B131 (2014).
［１７］M. Lindqvist et al., J. App. Phys. 126, 123303 (2019).
［１８］J. Lettry et al., Rev. Sci. Instrum. 85, 02B122 (2014).
［１９］K. Miyamoto et al., AIP Conference Proc. 1655, 050004(2015).
［２０］K. Nagaoka et al., The 7th International Symposium on
Negative Ions, Beams and Sources, 1-11 Sept., videoconference (2020).
［２
１］池田勝則 他：プラズマ・核融合学会誌 95, 327 (2019).

強力なツールとなり得る．PIC-MCC を用 い た 負 イオン
ビームの加速に関するシミュレーション例については，文献［２１］等を参照していただきたい．

参考文献［１］C.K. Birdsall and A. B. angton: PLASMA PHYSICS VIACOMPUTER SIMULATION (Taylor & Francis, 2005).
［２］内藤裕史：プラズマ・核融合学会誌 74, 470 (1998).
［３］石黒静司：プラズマ・核融合学会誌 74, 591 (1998).
［４］K. Miyamoto et al., AIP Conference Proc. 1515, 22 (2013).
［５］高木俊宜：電気学会大学講座 電子・イオンビーム工学，電気学会 p.48 (2005).
［６］G. Fubiani et al., New J. Phys. 19, 015002 (2017).
［７］S. Kuppel et al., J. Appl. Phys. 109, 13305 (2011).
［８］K. Miyamoto et al., Appl. Phys. Lett. 102, 023512 (2013).
［９］T. Takizuka and H. Abe, J. Comput. Phys. 25, 205 (1977).

５９１

J. Plasma Fusion Res. Vol.97, No.10 (2021)５
９２‐５９５

小特集 イオン源の数値シミュレーション技術の進展と応用

６．まとめ
6. Summary
宮 本 賢 治，吉 田 雅 史１）＊，畑 山 明 聖２）
MIYAMOTO Kenji, YOSHIDA Masafumi１）＊ and HATAYAMA Akiyoshi２）
鳴門教育大学，１）山口大学，２）慶應義塾大学（原稿受付：２
０
２
１年６月１
７日）

６．
１ イオン源プラズマシミュレーションモデル
の統合による将来的なイオン源開発の発展
性への貢献

おいては，ITER や原型炉等の将来的な実機負イオン NBI

１９８０年〜１９９０年代における負イオン生成部の解析には，べき課題の１つとして，電子電流を低減する必要がある．

プラズマ中の各粒子密度・温度を空間的に一様であるとし負イオン源ではいわゆる セシウム（Cs）効果 により，負て，それらの時間変化を反応のレート方程式により求める

イオン電流密度の増大と同時に電子電流密度の低下を図

０次元モデル を用いることで，負イオン電流の運転ガス

る．

における重水素負イオン源開発への貢献が挙げられる．安定な負イオン生成や高エネルギー加速を図るために解決す

圧依存性やアークパワー依存性等の実験結果を定性的に説

しかしCs効果にもかかわらず，重水素放電は水素放電に明できた［１，
２］．その後，モンテカルロ法などの導入により比べて，負イオンと同時に引き出される電子電流が多いこ

負イオン生成部の数値計算コードが進展し，電子エネルとが実験的に報告されている［５，６］．負イオンビームの安ギー分布関数や粒子密度の空間分布を求められるように定な加速や長時間運転を実現する上で，引き出される電子なった．その結果，例えば負イオン源プラズマの非一様性

電流の抑制は克服すべき課題の１つである．そのために，などのような，空間分布を要する物理現象への適用が可能

重水素の場合，水素よりも電子電流が増大する物理メカニになった．

ズムを解明することが必須である．現在，重水素実験結果

また負イオンの引出・加速については，従来の数値計算と KEIO-MARC（Kinetic model of Electrons in IOn source

コードは初期のメニスカス形状を仮定して，負イオンビー

with Multi-cusp ARC discharge）コードや PIC-MCC コード

ム軌道を計算していた．しかし PIC-MCC（Particle-In-Cellとの比較・検討による解析が進められている［７］．

with Monte-Carlo Collision）の導入により，負イオンビーム

また重水素の方が水素よりも，Cs 効果の持続時間が短い

軌道を決定するプラズマメニスカス形状を仮定せずに自己ことが実験的に報告されている［８］．これもまた実機の長矛盾なく解けることができるようになった［３］．その結果，

時間運転を実現する上で克服すべき課題の１つである．

従来の数値計算コードでは再現できなかったビームハロを

Cs 効果の持続時間が短くなる理由として，"）引出・加シミュレーションで再現することができるようになった

速部から逆流した正イオンによるプラズマ電極（PG）表面［４］等，負イオンの引出・輸送の物理メカニズムの解明への Cs 吸着層への直接的な剥離［９］や，または#）アーク

貢献することが可能になった．

チャンバー壁表面に吸着した Cs 原子の PG への堆積による

このようにコンピュータの性能向上に伴い，特に運動論

仕事関数の増加［１０］などが挙げられる．シミュレーション

的モデリングを土台としたイオン源プラズマシミュレー

により上述のような Cs 輸送現象等を解析し，Cs 効果の持

ションの手法・モデルの進展により，シミュレーションは続性に関する実験結果についての物理的メカニズムの詳細

負イオン源内で生じる複雑かつ多層的な物理現象を明らかな検討を行う．それにより，Cs の有効利用・最適な導入条にし，実験結果を定性的に十分，説明することが可能に件の導出や引出・加速部への漏洩量低減の手法の確立に貢なってきた．以下では，第１〜５章までで解説したイオン

献ができると考えられる．核融合ダイバータ物理の解析用

源プラズマシミュレーションモデル及びその統合による将に開発した，モンテカルロ計算による不純物輸送コード

来的なイオン源開発の発展性や課題，及び個々の計算技法［１１‐１３］を Cs の挙動解析に適用することが挙げられる．たの他分野への応用性について議論する．

だしCsの質量数が大きく，運動の時間スケールが電子と大

６．
１．
１ 重水素負イオン源開発

きく異なるため，バルクプラズマを含めて運動論的モデリ

核融合炉における負イオンを用いた中性粒子入射装置

ングで自己矛盾なく，同時に解くことは計算コストの観点（負イオン NBI;Neutral Beam Injection）用の負イオン源にから現状では難しい．そのため不純物輸送コードにおいてNaruto University of Education, Naruto, TOKUSHIMA 772-8502, Japan＊現在の所属：宇部工業高等専門学校

corresponding author’s e-mail: kmiyamot@naruto-u.ac.jp

５９２!2021 The Japan Society of Plasma
Science and Nuclear Fusion Research

Special Topic Article

6. Summary

K. Miyamoto et al.

６．
２ 個々のモデリング・計算技法の課題

Cs の運動を計算する際には，バルクプラズマについては流

イオン源プラズマシミュレーションのモデリングにおけ体的モデリング［１４，１５］から計算し，バルクプラズマを背景として与える方法が現実的である．

る課題の１つにとして，空間スケールや時間スケールが，

６．
１．
２ 高周波放電型（RF）負イオン源の負イオンビーム

オーダーで広範囲に渡って変化する物理現象のモデリングの光学特性の改善をどうするか？というマルチスケーリングの問題が挙げら

核融合用，加速器用の RF 負イオン源に共通して，RFれる．空間スケールについて，例えば核融合用負イオン

負イオン源で生成された負イオンビームの発散角は大きい

NBI の負イオン源の実サイズは，実機だと数 10 cm〜数 m

ことが実験的に報告されている［１６］．その原因の１つとしのオーダーに対し，シース長は 10 μm〜数 100 μm のオーて，RF 揺動電場によるメニスカス形状の揺動が指摘されダーである．時間的スケールについては，シース内の正・ている［１７］．この物理現象の解明と負イオンビームの光学

負イオンの運動は数 ns のオーダーだが，表面生成過程は数的特性の改善に向けた指針を示す上で，例えば PIC-MCC

ps〜数 fs のオーダーである．
そのため例えば第５章で述べたように，PIC-MCC では

コードによる検討が有効であると考えられる．
６．
１．
３ Cs フリー負イオン源の開発

スケーリングにより引出領域をモデリングする．しかしこ

負イオン電流値を増大させると共に電子電流値を抑制すの場合，バルクプラズマ中のポテンシャルと引出電圧をスるために，負イオン源に Cs を添加し，PG 表面での表面生ケーリングして求めたポテンシャルとの整合性が必ずしも成を促進させ る 方法 が 一 般 的 で あ る．し か し こ の方法

取れてない．そのためにimplicit PIC法
［２５］などの適用によは，６．
１．
１でも述べたようなプラズマや引出・加速部からる実サイズでのモデリングを進めることが考えられる．実の逆流イオンによる Cs 効果の劣化や引出・加速部への Cs

際 RF イ オ ン 源 に つ い て は，implicit PIC 法 を CERN の流入による耐電圧の劣化等の解決すべき問題がある．そのLinac4負イオン源のモデリングへ既に適用した事例がある

ため現在，負イオン源へCs添加をしないで高電流密度の負［２６］．しかし，空間スケールや時間スケールがオーダーでイオン生成を行うために， Cs フリー負イオン生成 の研

広範囲に渡って変化する物理現象を，運動論的方法のみで究が勢力的に行われている．現存するCsフリー負イオン生モデリングするのは難しい場合もあり得る．その場合は流成には，国内で以下の５つの方法が提案されている．

体モデルとの併用などが考えられる．

（a）シートプラズマを用いた体積生成法［１８］

また他の課題として，バルクプラズマ‐表面間の相互作（b）Ta フィラメントを用いた体積生成法［１９］用のモデリングが挙げられる．例として，1〜2 eV 程度の（c）C12A7 エレクトライド製プラズマグリッドを用いた運動エネルギーでフランクコンドン原子が PG 表面に入射して，H−が表面生成される過程をモデリングする際に，こ表面生成法［２０］（d）Al 製プラズマグリッドを用いた表面生成法［２１］のような低いエネルギー領域では量子力学的な効果が重要（e）ヘリコン波負イオン源［２２］になる．そのため従来の二体衝突近似（BCA:Binary Colli-（a）については，厚さ数ミリメートルのシート状の狭い

sion Approximation）に基づく数値シミュレーションコー領域で高密度のプラズマや励起子を生成して，シート近傍

ド（代表的なものとして，TRIM［２７］や ACAT
［２８］がある）での体積生成により負イオン生成する．
（b）については，の適用妥当性が問題となる．バルクプラズマ‐表面間の相Ta フィラメントによる壁表面への蒸着により，水素振動互作用のモデリングにおいては，分子動力学シミュレー

励起分子の高い振動励起準位（!!!）のポピュレーション

ションの適用が考えられる［２９］．

を増大させて，負イオンの体積生成を増大させる方法であ６．
３ モデリング・計算技法の他分野への応用性

る．（c）に つ い て は，C12A7 electride 表 面 の 仕 事 関 数 を
2.4 eV 近くまで下げることで負イオンを生成する方法で，第１章で述べられているように，イオン源は，核融合プ

Cs を用いた表面生成と原理的には近い．
（d）については，

ラズマ加熱，大型加速器，がん治療など医療機器，探査衛

Al 表面上で負イオンを生成する方法で，仕事関数に強く依星のエンジン，イオン注入など半導体製造など，多分野に存せずに生成が可能である．現在は生成機構の解明を図っ

おける粒子源として用いられており，本特集で紹介したシている．（e）についてはヘリコン波で高密度プラズマを生ミュレーションもまたこれらの分野での研究開発で適用さ成し，体積生成により負イオンを生成する．

れている．他分野への応用性については，例えば質量分析

これらの Cs フリー負イオン生成による負イオン源の開器のイオン源やイオン光学系への適用が考えられる．質量発においても，今後シミュレーション側からの貢献が期待分析の分野では，イオンビーム軌道計算用に Simion というされる．

市販の計算コードがよく用いられる［３０］．Simion はラプ

６．
１．
４ その他

ラス方程式を差分法で解いて，ポテンシャルを求める．イ
オンビームの軌道計算では，クーロン衝突や中性粒子とのその他としては，大面積負イオン源の更なるビーム一様衝突を考慮しない．

性の改善で，長手端部のビーム電流向上［２３］や，アーク放質量分析においてイオンビームの収束性を良くするため

電式でのフィラメントアーキング低減と放電効率の改善に，イオン源に Xe などのガスを外部から導入して，ガス分［２４］への貢献が考えられる．

子とイオンとの衝突によりイオン温度を冷却する方法がある［３１］．ただし，ガス流量が多すぎると，イオンを散乱さ
５９３

Journal of Plasma and Fusion Research Vol.97, No.10 October 2021

せて，かえってイオンビームの発散角が大きくなる．その理メカニズムを抽出することが大事である．０次元モデル

ため適切なガス流量が存在するが，設計計算を行う際にはモデルとしては簡易ではあるが，これまでの負イオン生

Simion ではシミュレーションできない．そのため，例えば成部の解析において有効性を示してきた．他方，３次元のPIC-MCC の適用により，ガス分子とイオンとの衝突過程

シミュレーションコードはモデリングが複雑なため，得らを考慮したシミュレーションを行うことが考えられる．

れた結果が真に物理的に意味のある結果なのか，あるいは
数値計算上のエラーから発生した

また質量分析では，土壌に含まれるダイオキシンや生物

人工的な

結果を識別

の体内における代謝物などの絶縁物を測定試料としてすしにくくなる．計算結果の解釈を誤らないためには，０次る．質量分析器を長時間使用すると，分解能や感度といっ元モデルの併用による計算結果の物理的な検討は必須であた性能が劣化するため，イオン源や電極を洗浄しなければるといっても過言ではない．また，計算コストがかかる２ならない．この性能劣化の原因として，絶縁物がイオン源

次元や３次元モデルのシミュレーションの解析方針の策定

や電極に付着すると帯電し，局所的に不適切な電場を形成

や計算ケースの選択などのためにも，０次元モデルによるすることが考えられるが，今のところ物理メカニズムは解幅広いパラメータ・サーベイは極めて有効な手段になる．
２つ目として，数値計算技法にとらわれた

明されていない．PIC-MCC の適用による解析が考えられ

シミュレー

るが，測定試料となる有機化合物の原子・分子データの拡ションのためのシミュレーション

ではなく，実験との比張や表面での誘電体の帯電現象のモデリングを行う必要が較・検討によるモデリングの検証・改善である．そのためある．

にも，国内・国外との実験研究者と連携して研究を進める

その他の他分野への応用性としては，正・負荷電粒子数

べきである．また，例えば表面物理などこれまでとは異なの比が１に近いペアプラズマにおけるプラズマ振動周波数

る専門分野の研究者との連携により，今後のイオン源プラ

やプラズマの分散関係などの物性解明が挙げられる．銀河

ズマのシミュレーション手法・モデリングのさらなる発展

中心核，パルサー磁気圏に存在する電子・陽電子プラズマが望まれる．

の理解は，宇宙プラズマ物理の重要な課題であり，ペアプ
ラズマの分散関係など諸特性はこれまでほとんど未解明で参考文献ある．一方，負イオン源のプラズマ電極近傍には，表面生［１］O. Fukumasa, J. Phys. D: Appl. Phys., 22, 1668 (1989).
［２］M. Ogasawara et al., Fusion Eng. Des. 26, 507 (1995).
［３］K. Miyamoto et al., Appl. Phys. Lett. 100, 233507 (2012).
［４］K. Miyamoto et al., Appl. Phys. Lett. 102, 023512 (2013).
［５］H. Nakano et al., Jpn. J. App. Phys. 59, SHHC09 (2020).
［６］B. Heinemann et al., New J. Phys. 19, 015001 (2017).
［７］加藤 凌：NIFS 重水素負イオン源における電子輸送の数
値解析，第３
７回プラズマ・核融合学会年会
（２
０
２
０．
１
２．
１
‐
１
２．
４ 中継本部：愛媛大学常北キャンパス）．
［８］P. Franzen et al., Nucl. Fusion 55, 053005 (2015).
［９］A. Mimo et al., AIP Conf. Proc. 2052, 040009 (2018).
［１０］M. Yoshida et al., Rev. Sci. Instrum. 91, 033313 (2020).
［１１］Y. Homma, A. Hatayama, J. Comput. Phys. 231, 3211
(2012).
［１２］Y. Homma, A. Hatayama, J. Comput. Phys. 250, 206 (2013).
［１３］S. Yamoto et al., Comp. Phys. Comm. 248, 106979 (2020).
［１４］T. Mizuno et al., Rev. Sci. Instrum. 75, 1760 (2004).
［１５］S. Lishev et al., Plasma Sources Sci. Technol. 27, 125008(2018).
［１６］W. Kraus et al., AIP Conf. Proc. 2052, 040004 (2018).
［１７］Y. Haba et al., Jpn. J. Appl. Phys. 59, SHHA01 (2020).
［１８］H. Kaminaga et al., Rev. Sci. Instrum. 91, 113302 (2020).
［１９］S. Masaki et al., AIP Conf. Proc. 2052, 020002 (2018).
［２０］M. Sasao et al., Appl. Phys. Express 11, 066201 (2018).
［２１］M. Yoshida et al., Phys. Plasma, 28, 033512 (2021).
［２２］負イオン研究会２０１８「負イオンの基礎と実践の最前線」
https://workshop.nifs.ac.jp/n-ion2018/［２
３］M. Yoshida et al., Fusion Eng. Des. 96-97, 616 (2015).
［２４］M. Ichikawa et al., Plasma Fusion Res. 11, 2405108 (2016).
［２５］J.U. Brackbill and D. W. Forslund, J. Comput. Phys. 46, 271(1982).
［２６］S. Mattei et al., J. Comput. Phys. 350, 891 (2017).
［２７］W. Eckstein et al., J. Appli. Phys. 51, 1194 (1980).
［２８］W. Takeuchi, Y. Yamamura, Rad. Effects 71, 53 (1983).

成により多量の負イオンが存在し，正・負イオンから成る
ペアプラズマ（イオン性プラズマ）が存在する［３２］．そのため，負イオン源の開発への直接的な寄与ではないが，実験室プラズマ研究から宇宙プラズマ研究への新たな貢献という観点から，基礎研究への貢献が期待される．

６．
４ 最後に
シミュレーションのモデリングにおいては，次の２点が重要である．まず１つ目は，適材適所なモデリングを行うことである．冒頭で述べたが，最近の負イオン源のシミュ
レーションにおいては特に運動論的モデリングを土台として，シミュレーションの手法・モデリングが進展してきた．しかしながら 6.2 でも述べたように，イオン源プラズマには空間的・時間的にマルチスケーリングの物理現象が多く存在する．そのため，運動論的モデリングでは上で述べたようなスケーリングなどの工夫なしには，計算コストの観点から実質的にシミュレーションが不可能な場合も多くある．したがって運動論的な効果は犠牲にし，流体モデル
による解析を選択することが望ましい場合もある．例えば，!）先に述べたように，背景となるバルクプラズマには，
流体モデルを適用し，対象粒子（例えば Cs のみ）の運
動は運動論的モデルを用いる"）いわゆるheat flux limit［３３］を磁力線方向の輸送係数に採用するなど，できる限り運動論的な効果を取り入れた流体モデルを開発する，など柔軟な対応が求められる．
また３次元のシミュレーションコードと０次元モデルとの併用により，複雑かつ多層的な物理現象から本質的な物
５９４

Special Topic Article

6. Summary［２９］M. Warrier et al., J. Nucl. Mater. 467, 457 (2015).
［３０］SIMION, https://simion.com/［３１］A.V. Tolmachev et al., Rapid Commun. Mass Spectrom.
14, 1907 (2000).

K. Miyamoto et al.

［３２］K. Tsumori et al., Rev. Sci. Instrum. 83, 02B127 (2012).
［３３］P.C. Stangeby, The Plasma Boundary of Magnetic Fusion
Devices (Taylor & Francis Group, New York), Chap. 26,
pp.656-666 (2000).

５９５

!!
しば

た

たか

小特集執筆者紹介

のり""わ

柴 田 崇 統高エネルギー加速器研究機構（KEK）加速器
研究施設・助教．２
０１
４年慶應義塾大学大学院
理工学研究科・後期博士課程修了，博士（工
学）
．２
０
１
５年 KEK 特別助教を経て，現職．主として J-PARC 加速器用のイオン源開発に従事しています．最近
インターンシップ用にイオン源の紹介・組立て動画を制作したのですが，BGM 選定が難しく，どうやっても YouTube のクッキング動画風になる問題に直面しています．

だ

もとい

和 田

元

同志社大学理工学部・教授．１
９８
３年カリフォ
ルニア大学バークレー校 Ph. D.取得，日立製作
所を経て１
９
８
７年同志社大学着任．専門は弱電離プラズマと固体壁の相互作用．講義が多す
ぎて自分の研究は殆どさせて貰えないので，学生さんの研究のお手伝いをしています．感染症対策で毎日がビデオ教材づくり，ゴールデンウィークもお盆休みもありませんでした．
けん

よし

だ

まさ

吉 田 雅 史

だ

けい

みや

じ

ます．
かず

もと

けん

じ

鳴門教育大学大学院学校教育研究科 教授．
博士（工学）（慶應義塾大学大学院理工学研究科）．旧 日 本 原 子 力 研 究 所（現 QST）で は，
ITER等におけるNBI用負イオン源の開発に従事してきました．大学では負イオン源の数値計算シミュレー
ションの他，学校教育に関する授業や研究も行っています．今
年６歳になる息子がいますが，最近，将棋を覚えたり，プール
で水に顔つけができるようになったりと日々の成長に驚かされます．

信州大学工学部教授．主な研究分野はプラズ
マ分光学．しばらくぶりに野球と空手で動けるようにしたいなあと思いつつ，鍛錬のため
近所の善光寺とか戸隠とかのお散歩をしていのひろ

宮 本 賢 治

澤 田 圭 司ほしたか同志社大学生命医科学部 教授 博士（理学）ここ１
０年ほど，PWI研究から離れて，生命科学分野で研究を進めて お り ま す．昨 今 の エ ネ ル
ギー問題や SDGs に関連して，核融合エネル
ギーに対する社会的な期待が高まっておりますので，核融合
コミュニティの研究者の皆様方の研究成果が身を結ぶことを願っております．

宇部工業高等専門学校 電気工学科 准教授．
２
０１
１年に九州大学大学院総合理工学府・博士（工学）．日本原子力研究開発機構，量子科学
技術研究開発機構，山口大学を経て，今年９月より現職．核融合プラズマ壁相互作用に始まり，現在はプラズ
マ加熱用 NBI の負イオン生成とその引出加速に関わっていま
す．４年ほど前からソフトバレーボールをはじめました．今はサーブとレシーブにはまっています．
さわもつ

剣 持 貴 弘ふみお

星 野 一 生はた

慶應義塾大学理工学部・准教授．２
００
８年に慶應義塾大学大学院理工学研究科・後期博士課
程を修了（博士（工学））．日本原子力研究開発機構，量子科学技術 研 究 開 発 機 構 を 経 て
２
０
１
８年４月より現職．周辺・ダイバータプラズマおよびイオ
ン源プラズマのモデリング・シミュレーション研究を行っています．昨年からの新型コロナで在宅が増え，家中のスマート化がだいぶ進みました（ただし体型を除く）．最初は興味本位でしたが，今や頼りっきりの生活です．

やま

あき

よし

畑 山 明 聖慶應義塾大学 名誉教授，主な研究分野：プラズマ物理，ダイバータプラズマ，負イオン源プ
ラズマ，数値シミュレーション．私事恐縮ですが，退職直前（２０１
９年１月）不覚にも心不全で倒れ，ICU で生死の淵を彷徨いました．これは，そろそろ若手に任せてフェードアウトしろとのサインだったようにも思いますが，家族，そして皆さんのお蔭で何とか再起することができました．現在，生きていることの喜びを感じています．この小特集についても，執筆者の皆様，特に若手研究者の皆様との議論を通していっぱい元気をもらいました．大変，感謝しています．

５９６