平成２１年度
九州大学大学院 総合理工学府
先端エネルギー理工学専攻

修 士 論 文論文名
宇宙線中性子起因
ソフトエラーシミュレーションに関する研究

氏

名

指導教員名

安部 晋一郎

渡辺 幸信 准教授

目次
第1章
1.1
1.2

1
研究背景 ・・・・・・・・・・・・・・・・・・・・・・・・・・1
研究目的 ・・・・・・・・・・・・・・・・・・・・・・・・・・3
序論

第 2 章 ソフトエラー発生原理および電荷付与に関わる現象について
4
2.1 高エネルギー中性子起因ソフトエラー発生原理 ・・・・・・・・・42.2

電荷付与過程に関わる物理現象 ・・・・・・・・・・・・・・・・5

第 3 章 ソフトエラーシミュレーション手法
7
3.1 BGR とソフトエラー率の定義 ・・・・・・・・・・・・・・・・8
3.2 PHITS ・・・・・・・・・・・・・・・・・・・・・・・・・・・8
3.3 核反応モデル ・・・・・・・・・・・・・・・・・・・・・・・・93.3.1 INC モデル ・・・・・・・・・・・・・・・・・・・・・・・93.3.2 QMD モデル ・・・・・・・・・・・・・・・・・・・・・・10
3.3.3 GEM モデル ・・・・・・・・・・・・・・・・・・・・・・13
3.3.4 イベントジェネレーターモード ・・・・・・・・・・・・・・14
3.4 阻止能計算コード・・・・・・・・・・・・・・・・・・・・・・15
3.5 ストラグリング分布関数・・・・・・・・・・・・・・・・・・・17
第 4 章 実デバイス体系を模擬した計算結果および考察
20
4.1 計算体系および計算条件・・・・・・・・・・・・・・・・・・・20
4.2 計算体系サイズに関する計算結果の収束性の確認・・・・・・・・24
4.3 中性子の入射エネルギーに関する確認・・・・・・・・・・・・・26
4.4 二次イオン種の寄与・・・・・・・・・・・・・・・・・・・・・28
4.5
4.6
4.7
4.8

核反応モデルの違いに関する計算結果の比較・・・・・・・・・・30
阻止能計算コードの違いに関する計算結果の比較・・・・・・・・31
デバイスの材質の違いに関する計算結果の比較・・・・・・・・・32
有感領域サイズを変更したときの計算結果の比較・・・・・・・・35

第 5 章 電荷付与過程に関わる物理現象の影響調査結果および考察
36
5.1 ストラグリングによる影響の調査結果・・・・・・・・・・・・・36
5.2 イオントラック構造による影響の検討・・・・・・・・・・・・・39i41
6.1 核反応モデルに関する実験値との比較および考察 ・・・・・・・・41
6.1.1 Si(n,p)X に関する実験値との比較 ・・・・・・・・・・・・・42
6.1.2 Si(n,4He)X に関する実験値との比較 ・・・・・・・・・・・・47
6.2 阻止能計算コードと実験値との比較および考察 ・・・・・・・・・52

第6章

各計算モデルの妥当性の確認

第7章

まとめおよび結論

55

参考文献

57

謝辞

59

ii

1. 序論1.1 研究背景
近年、電子機器の信頼性問題として二次宇宙線中性子によるソフトエラーが注目されている。ソフトエラーとは、放射線が半導体デバイスに入射したことによって生じる一時的な誤動作・故障のことである。特に大型計算機やペースメーカーなど、高い信頼性が要求される機器への影響が懸念されている。また半導体デバイスの微細化に伴いソフトエラー発生率は上昇しており、将来的には一般家電や自動車などに関しても影響が現れると考えられている。[1,2]図 1.1-a に示すように、二次宇宙線中性子は宇宙から降り注ぐ一次宇宙線と大気との核反応によって生成される。地上での中性子スペクトルは時期や位置によって異なるが、日本ではおおよそ図 1.1-b に示すようなスペクトル[3]である。また 1MeV 以上の中性子は 1 時間に 1cm2 あたり 25 個程度降り注いでいる。
ソフトエラーに関する研究の歴史は古く[1,2,4]、1975 年に高エネルギー宇宙線（陽子、He イオン、Fe イオンなど）によって人工衛星の電子機器に異状が発生していると報告されたのが始まりである[5]。この当時は、高エネルギー宇宙線が直接地上まで届く
ことはほとんどないため地上で使用する電子機器にとってソフトエラーは深刻な問題ではないと考えられていた。しかし 1978 年、Intel 社の May らによって地上においてもソフトエラーが生じているとの指摘がなされた[6]。このときの原因はパッケージ
材料中に含まれていた微量の放射性不純物から放出されるアルファ粒子であった。
4

Neutrons [/MeV/cm2/hour]

10

103
102

二次宇宙線中性子スペクトル@ 東京 2009/01/01(cal. by EXPACS ver.2.17)

101
100
10-1
10-2
10-3
10-4
0.001 0.01

0.1

1

10

100

Neutron Energy [MeV]

図 1.1

1000

(a) 二次宇宙線生成過程の模式図 および (b) 東京での二次宇宙線中性子のエ

ネルギースペクトル（EXPACS[3] により計算）
1表 1.1

ソフトエラー研究の歴史[1,2,4]

年代

内容

1975 年

宇宙空間において高エネルギー宇宙線によるソフトエラーが確認される。

1978 年

地上においてアルファ粒子によるソフトエラーが確認される。

1979 年
1980 年代～
1990 年代前半
1990 年代後半
2002 年
2005 年

核反応を介したソフトエラーの発生が報告され、中性子入射であってもソ
フトエラーが発生することが指摘される。
アルファ粒子起因のソフトエラーに関して研究・対策が進む。
複数の研究グループから中性子起因ソフトエラーの危険性が指摘される。
半導体デバイスの微細化に伴いソフトエラー率が上昇すると報告され、この問題の重要性が広く認知される。
JEDEC や JEITA がソフトエラーの評価手法の標準文書をまとめる。

1979 年には IBM 社の Ziegler らによって二次宇宙線の陽子や He イオン、また電荷を持たない中性子であっても、核反応過程を介することでソフトエラーを引き起こす
ことが可能であると示された[7]。そして同年に Guenzer らによって行われた半導体チップへの陽子・中性子照射実験によって、核反応が原因でソフトエラーが生じることが確認された[8]。さらに 1983 年に Dicello によって実施された実験では π 中間子によって LSI にソフトエラーが発生することが確認され[9]、全ての放射線はソフトエラーを引き起こすことが可能であることが判明した。
1980 年代に入るとアルファ粒子起因のソフトエラーに関する研究が多く為され、パッケージ材料中の放射性不純物の除去やチップの耐性向上など対策が進められた。
1990 年代後半になると二次宇宙線中性子によるソフトエラーの重要性が指摘され始め、DRAM のソフトエラーの主要因は中性子であることが確認された[10]。また同時期に論理回路での中性子起因ソフトエラーが確認され、その影響は無視できないと
指 摘 さ れ た [11] 。 2000 年に開催された半導体の信頼性に関する国際会議『IEEEInternational Reliability Physics Symposium(IRPS)』では、Texas Instruments 社（TI 社）の Baumann らによって BPSG 膜中の 10B が熱中性子を捕獲することによって放出されるアルファ粒子がソフトエラーを引き起こすとの報告があった[12]。
2002 年には TI 社の Baumann や日立製作所の伊部らによって半導体デバイスの微細化によりソフトエラーの発生確率が上昇すると報告され[13,14]、ソフトエラーは切実な問題として広く認識されるようになった。2005 年には米国の半導体技術協会（JEDEC）
や日本の電子情報技術産業協会（JEITA）によってソフトエラーの評価手法の標準文書が作成され、現在もソフトエラーに関する研究は活発に行われている。
2

ソフトエラーの評価手法としては、フィールド試験、加速器試験、シミュレーショ
ンの 3 つが挙げられる。フィールド試験とは、多数のデバイスを数ヶ月から数年間連続稼働させてソフトエラーの発生率を直接測定する方法である。この方法は現実に則したデータが得られる点では非常に有効であるが、実験期間が非常に長いため安易に行うことができないのが欠点である。加速器試験はデバイスに放射線を照射してソフ
トエラー発生率の測定を行うため、フィールド試験と違い短期間で測定が行えるのが利点である。しかし中性子に関しては連続スペクトル中性子源施設が限られており、またコストが非常にかかる面も問題となる。
これら実測的な方法はデバイスの品質を示すためにも重要であるが、実際にデバイスを製造しなければ実験が行えないために、測定結果をデバイスの研究開発へ反映するためには一定の期間を要する。一方シミュレーションは設計段階でソフトエラー率を評価でき、その結果を早期に設計にフィードバックすることでデバイスの研究開発
を効率よく進めることに役立つ。シミュレーション結果の信頼性を確保するためには、
ソフトエラー発生の物理過程を精度よくモデル化することが重要となる。
当研究グループでは、宇宙空間における陽子起因ソフトエラーおよび地上における宇宙線中性子起因ソフトエラーについて、半経験的モデル（McSEE-Q code：Monte
Carlo simulation in Single-Event Effects-Kyudai version）を作成し、ソフトエラー発生確率の評価などを行ってきた。[15,16]
現在は設計ルール 45nm 以降の最先端デバイスを対象とした宇宙線中性子起因ソフ
トエラーシミュレーションの高精度化を目標として研究を行っており、半導体理工学
研究センター（STARC）および筑波大学と共同研究を進めている。

1.2 研究目的
本研究ではソフトエラーの発生過程のうち、核反応過程および電荷付与過程において重要となる要因を解明することを目的とする。
まず、実デバイスを模擬した計算体系への宇宙線中性子入射について粒子・重イオ
ン輸送計算コード PHITS[17]を用いてシミュレーションを行い、種々の要因がシミュレ
ーション結果へ与える影響について定量的に評価を行う。
また生成二次イオンによる電荷付与量計算の精度向上を図るため、微細化に伴い影響が顕在化すると予測される現象（ストラグリングおよびデルタ線によるイオントラ
ック構造形成）による、ソフトエラーシミュレーションへの影響について調査を行う。
さらに、PHITS に搭載されている核反応モデルおよび阻止能計算コードの精度を確認するため、計算結果と実験値との比較を行う。
3

2. ソフトエラー発生原理及び電荷付与に関わる現象について
2.1 高エネルギー中性子起因ソフトエラー発生原理
中性子起因ソフトエラーの発生に関する模式図と時間・空間発展の様子を図 2.1 に示す。まず、中性子が半導体デバイス内に入射したとき、ある確率でデバイス材料と核反応を起こす。この核反応によってさまざまな二次イオンが生成される。生成された二次イオンはデバイス内の電子や原子核との相互作用を介してエネルギーを付与し、その結果二次イオンの飛程に沿ってデバイス内に電子・正孔対が発生する。空乏層内やその近傍に発生した電荷はドリフト、ファネリング *、拡散などによってメモリセルに収集される。この収集された電荷がある臨界量を超えたとき、メモリセル内の情報の反転（0 → 1 または 1 → 0）が起こり、ソフトエラーが引き起こされる。[1,2]このようにソフトエラーの発生には広範囲の空間・時間スケールに亘る物理が関与しているため、シミュレーションの高精度化のためには各段階の物理過程を精密よくモデル化する必要がある。

図 2.1 ソフトエラー発生過程に関する模式図（左上）および時間・空間発展* 二次イオンが空乏層を通過したとき、二次イオンによって発生した電荷が空乏層の電場を緩和し、それにより空乏層外の領域からも電荷が引き寄せられる現象
4

2.2 電荷付与過程に関わる物理現象
半導体デバイスの微細化に伴い電荷付与計算に対して影響が顕在化すると考えられる現象として、エネルギーストラグリング（Energy Straggling）とデルタ線によるイ
オントラック構造の形成が挙げられる。
2.1 節で述べたように、デバイスへの電荷付与は生成二次イオンとデバイス中の電子や原子核との相互作用によって生じる。この相互作用は確率論的に起こるため、実際の電荷付与量はある平均値の周りで揺らいでいる。微小領域では関与する電子や原子核が尐なくなるため、相互作用の回数が減尐し揺らぎの幅が大きくなる。この統計的な揺らぎをエネルギーストラグリングと呼ぶ。これまでのソフトエラーシミュレー
ションの多くはストラグリングを無視しており、単位長さ当たりのエネルギー付与量である阻止能（dE/dx）の平均値を用いてソフトエラー率を算出していた。しかし実際には、平均の電荷付与量が臨界電荷量を下回っている場合でも、ストラグリングの影響によって臨界電荷量を上回るイベントがいくつか発生している可能性がある。今後デバイスの微細化がさらに進んだとき、ストラグリングがシミュレーション結果へ与える影響が顕在化することが予想される。
二次イオンがデバイス中を通過するとき、その飛跡に沿って大量の電子・正孔対を生成する。その多くは飛跡近傍に生成されるが、デルタ線と呼ばれる数 keV 程度の比
較的高いエネルギーを持つ電子の放出が確率的に起こり、イオンの飛跡から離れた場
所まで電子輸送が起こる。またデルタ線によって更に電子・正孔対が生成される。その結果、電荷の空間的な広がり（イオントラック構造）が形成される。次世代デバイ
スの微細化は有感領域の縮小をもたらすため、イオントラック構造を厳密に考慮したシミュレーションが必要となると考えられる。

図 2.2

エネルギー付与量の揺らぎ（左図、文献[18]より転載）
デルタ線によるイオントラック構造の形成（右図、文献[19]より転載）
5図 2.3 にストラグリングとイオントラック構造が電荷付与過程に与える影響について予測された結果を示す[20]。平均の阻止能を用いた計算で十分な領域（領域 1）、ス
トラグリングを考慮すべき領域（領域 2）、ストラグリングに加えてイオントラック構造を考慮すべき領域（領域 3）で構成されている。またそれぞれの領域の境界は、個々の現象による効果がエネルギー付与計算に対して最大 10%の変化を引き起こす場合に対応している。デバイスの微細化が進むに従って、後述する有感領域のサイズも小Sensitive volume
Mean Chord Length (nm)

さくなるため、両者の影響を考慮する必要性が示唆されている。

dE/dx only
(region 1)

103

102

dE/dx &
straggling(region 2)

dE/dx & straggling& track structure(region 3)

1

10

1

10

100

Incident ion energy (MeV/amu)

図 2.3

1000

ストラグリングとイオントラック構造に関する影響予測（文献[20]より）

6

3. ソフトエラーシミュレーション手法
本研究で用いたソフトエラーシミュレーション手法について述べる。前述したようにソフトエラーの発生過程は、(1) 核反応過程、(2) 電荷付与過程、(3) 電荷収集過程の 3 つに分けられる。核反応過程および電荷生成過程に関しては PHITS[17]と呼ばれる
汎用粒子輸送計算コードを採用する。また電荷収集過程に関しては、図 3.2 に示すよ
うにデバイス内部に有感領域を定義し、領域内に付与された電荷がすべてメモリセルに収集されるという簡易的なモデルを用いる。

図 3.1

ソフトエラーシミュレーション過程

図 3.2

有感領域モデルの概念図

7

3.1

BGR とソフトエラー率の定義

シミュレーションでは、まず PHITS を用いて有感領域に付与されたエネルギー分布

N ( Ed ) を計算する。ここで、シリコン中で 1 対の電子・正孔対を生成するために必要なエネルギーが 3.6eV であることより、エネルギー付与量 E d と電荷付与量 q には以下の関係が成り立つ。

q[fC] e
 Ed [ MeV]  44.5  Ed [ MeV]
3.6  10 6(3.1)

ここで、 e は電気素量（1.602×10-19[C]）である。これよりエネルギー付与分布を電荷付与分布に換算し、入射中性子フラックス n と有感領域の体積 VSV から次式で与えられる BGR（Burst Generation Rate）を算出する。


BGR(Qc )  N (q)dq

Qc

(3.2)

n  VSVBGR(Qc ) と 1 ビット当たりのソフトエラー率 SER(Qc ) は以下の関係にある。
SER(Qc )  C  BGR(Qc )  n  VSV

3.2

(3.3)

PHITS

PHITS とは Particle and Heavy Ion Transport code System の略称で、ほぼすべての粒子と重イオンの核反応および粒子輸送を計算する汎用三次元モンテカルロシミュレー
ションコードのことである[17]。ここで、モンテカルロ法とはフォンノイマンにより開発された方法で、物理現象を確率的にとらえ、確率過程を記述する確率変数を導き、その過程を乱数を発生させて追跡する手法である。
PHITS は旧 日 本 原 子 力 研究 所 （ 現 在 の原子 力 研 究 開 発 機 構 ）で 開 発 さ れ た
NMTC/JAM を基に、東北大学および高度情報科学技術研究機構（RIST）において開発された。その利用範囲は加速器や原子炉の設計および安全評価、宇宙や高高度での被ばく線量評価などの航空宇宙分野、がんの放射線治療などの医療分野と、多岐にわたる分野で利用されている。
PHITS には核反応の動的過程に関する計算モデルとして Bertini による INC（Intra
Nuclear Cascade）モデル[21]、JAM（Jet AA Microscopic Transport Mode）モデル、JQMD（Jaeri Quantum Molecular Dynamics）モデル[22-24]を、蒸発過程に関する計算モデルとして SDM（Statistical Decay Model）、GEM[25]（Generalized Evaporation Model）が搭載されている。また、核データを用いた計算も行うことが可能である。さらに 20MeV
以下の中性子核反応については、核データを利用した計算をイベントジェネレーター
8

として模擬できる計算オプションであるイベントジェネレーターモード（e-mode）[26,27]も搭載されている。

また荷電粒子の阻止能計算に関しては SPAR（Stopping Power And Range）コード[28]および ATIMA（Atomic Interaction with Matter）コード[29]が搭載されている。さらに
PHITS にはストラグリングを再現するための分布関数として Gauss 分布、Vavilov 分
布、Landau 分布を適宜用いるオプションが搭載されている。

3.3 核反応モデル
3.3.1

INC モデルINC とは Intra Nuclear Cascade の略で、原子核内の核子をサンプリングして抽出し、抽出された核子と入射粒子との衝突を計算するものである[21]。INC モデルでの核子の取り扱いは古典的なもので、原子核内の核子間の相互作用は考慮されていない。
原子核内の核子の密度分布は、次式の Hofstadter の分布関数を用いている。

 (r) 1e

r c  Z1

(3.4)

1

ここで

r  核中心からの核内での距離c  密度が 1 2 になる r の値  r0 A1 3
r0  1.07  10

13

Z1  0.545  10

； A  質量数

cm

13

cm

である。計算上では原子核を同心球殻領域に分割した階段関数で近似している。各領域の外半径は次式の関係から決定する。
 ( r )   i   (0) ,i  1,2,3

1  0.9 , 2  0.2 ,(3.5)

 3  0.01

原子核内の陽子と中性子の分布は厳密には同一でないが、計算上では同一としている。
核子の運動量分布については、中性子や陽子が各領域で絶対零度で縮退したフェルミ分布を持つと仮定し、次式を用いている。

f ( p )  ci p 2(3.6)

ci は規格化因子で、上式の規格化をpF

f ( p )dp  領域 i 内の中性子（または陽子）の数 N i(3.7)とすると、 ci  3N i p F3 となる。 p F はフェルミエネルギーに対応する運動量である。
フェルミエネルギーは核子密度、核子の種類に依存する。一般に i 領域に対して
9

p Fi  ci i1 3 である。ここで、 ci は定数である。
ポテンシャルエネルギー分布は、各領域内のフェルミエネルギーと結合エネルギーの和で与えられる。計算上では、結合エネルギーは全ての核種および全ての領域に対して一律に 7MeV としている。
次に、飛程・衝突点・衝突反応の種類および衝突された粒子の運動量のサンプリング法について述べる。入射粒子が原子核と衝突し原子核内に侵入すると、核内領域におけるポテンシャルエネルギー分が入射粒子の運動エネルギーに加算される。核内での粒子の運動は相対論的二体衝突理論によって記述され、自由核子-核子散乱断面積つまり  nn ,  n p ,  p p が用いられる。
核内における衝突点と衝突反応の種類は同時に決定される。反応が散乱のとき、散乱核が微分散乱断面積を用いて標本抽出される。反応が π 中間子発生反応であれば、反応生成物は重心系で等方に放出される。核内の座標系は、実験室系、標的粒子静止系と重心系を適宜用いている。
衝突に関与したエネルギーと運動量ベクトルが定まると、これらは一時記憶され、その衝突点から衝突に関与した全ての粒子について粒子を出発させて、同じ計算を繰り返す。このようにしてカスケードが展開される。
衝突後の粒子のエネルギーが衝突のあった領域に対するフェルミエネルギーより小さい場合は、パウリの排他律によりその衝突は禁止されているので衝突はなかったこととする。衝突後の粒子のエネルギーが原子核外部場を基準にして測ってカットオ
フエネルギー以下になったら、その粒子の計算は打ち切る。カットオフエネルギーとは、以後の衝突に対しては二体衝突近似が適用できない場合に、そこで計算を打ち切るために設定されたものである。その値としては、すべての粒子について原子核表面におけるクーロンポテンシャルの半分の値が採用されている。これは、クーロンエネ
ルギー以下の陽子が中性子と衝突したとき、陽子エネルギーの大部分は中性子に与えられるため、中性子はさらにカスケードを続けることが可能となるためである。全クーロンポテンシャルをカットオフエネルギーとすると、このような過程を切り捨ててしまうことになるので、便宜的に半分の値を用いている。

3.3.2

QMD モデルQMD モデルとは、核子の状態を量子的に取り扱うことで核子間の相互作用を考慮した半古典的なモデルである[22-24]。INC モデルと違い、系全体の時間発展を追っている。核子は幅 L のガウス型波束で次のように表される。

i (r ) 

 (r  R ) 2 i

1
exp

 r  Pi 34
( 2L)
4L



10

(3.8)

ここで、Ri と Pi はそれぞれ i 番目の核子の位置と運動量の中心値である。全系の波動
関数は、これらの波束の直積でかけると仮定する。
Ri と Pi の時間発展は、ニュートン方程式と二体の衝突頄で記述される。ニュートン
方程式は時間に依存する変分原理から求められ、次のようになる。

  HRi
Pi

HP i 
R i

,

(3.9)

ハミルトニアン H は質量頄を含んだ粒子エネルギーと二体の相互作用エネルギーの和で与えられる。二体の有効相互作用としては、クーロン力、対称頄を含んだ Skyrme
型の相互作用を用いる。ガウス型の波束を用いるので、それらは次式のようになる。

H   Ei
i



1 A
1 B
 i     i 
2 0 i1    0 i



1
e2
c
c
 i j R  R erf R i  R j
2 i , j ( i )
i
j



Cs
 1  2 ci  c j  ij
2  0 i , j ( i )





4L



(3.10)



ここで2

Ei  mi  Pi

2

(3.11)

また erf はエラー関数、ci の値は陽子に対しては 1、中性子に対しては 0 である。<ρi>は核子同士の密度のオーバーラップ関数で、
  i     ij    dr i (r )  j (r )
j i

j i



  ( 4L) 3 2 exp  R i  R j  4 L
j i

2



(3.12)

ここで i (r )  dp
f i (r , p )( 2L) 3 2

 ( 2L) 3 2 exp  ri  R i  2 L2

(3.13)

計算では、 A  219.4 MeV、 B  165.3 MeV、   4 3 が使われている。これらの値は、非圧縮率 K  237.3 MeV、飽和密度    0  0.168 fm-3、束縛エネルギー16MeV となるように設定されている。また、対称パラメータ Cs  25 MeV、ガウス型波束の幅L  2.0 fm2 で計算を行っている。

核反応の時間発展はニュートン方程式に加えて、二体の衝突頄でも記述される。ま11

ず、二体の衝突頄は終状態のパウリ排他律因子（ 1  f (r, p, t ) ）を含む。2 粒子の衝突の判定は、ローレンツ共変な衝突変数と、そのエネルギーでの断面積 σ から定義される bmax    の大小によってなされる。二体衝突は瞬時の相互作用と考え、衝突イ
ベントは時空間の 2 点で定義される。
今回用いた QMD モデルには、高エネルギーの核反応を記述するために N（核子）、Δ（1232）、N*（1440）、π 中間子が含まれる。Δ と N*は、対称頄を除いて核子と同じ相互作用をし、π 中間子はクーロン相互作用だけをすると仮定する。そしてこれらの粒子の生成、消滅は二体衝突頄で扱う。以下に QMD モデルで扱われている衝突チャンネルを示す。ここで、B はバリオン、N は核子を表す。
1. Bi  B j  Bi  B j
2.

N N N 

3.

N N N

4.

NN NN

5.

N  N*  N  N

6.

N   

7.

N   N*(3.14)

   N*
衝突頄の他に、Δ, N*の共鳴状態からの崩壊による π 中間子の発生も考慮する。
8.

相対論的修正に関して、今回用いた QMD モデルでは運動学とエネルギーの表式に相対論を採用しているが、運動方程式自体は非相対論的である。この問題は高エネル
ギーで重要となってくるため、エネルギーの表式における変数にローレンツスカラー
量を適用し、またローレンツ収縮を初期状態のブーストされた原子核に適用すること
により、相対論的効果の重要な部分を近似的に取り込んでいる。
基底状態を作る方法としては、次のようなランダムパッキングの方法を採用している。まず半径 R0  1.124 A1 3  0.5 、表面パラメータ a  0.2 fm の Woods-Saxon 型の分布を用意し、その中にガウス型波束の中心 Ri をランダムに分布させる。次に、この空間分布より、密度分布、ポテンシャルエネルギーを計算する。局所的に Thomas-Fermi
近似を用いてランダムに抽出する。最後に、全エネルギーが以下の質量公式の値を満たしているかどうかチェックする。

Ebin  15.56 A  17.23 A

232

N  Z
3 Z 2e2 46.57


2A

5 1.24 A1 3

(3.15)

ここで、A は質量数、N は中性子数、Z は陽子数である。もし核子当たりの束縛エネ
ルギーが Ebin  0.5 に入っていれば、このサンプリングを基底状態として採用する。

12

3.3.3

GEM モデルGEM は蒸発過程と核分裂を取り扱うモデルで、n から 28Mg までの合計 66 個の核種を放出粒子として扱う[25]。親核を i、放出粒子を j、娘核を d とする。また親核の励起エネルギーを E、反応の Q 値を Q、質量数をそれぞれ Ai, Aj, Ad、電荷数を Zi, Zj, Zd とする。ε から ε+dε の間のエネルギーをもつ粒子 j の放出確率 Pj は、Pj ( )d 

2S  1 m  ( )  ( E  Q   ) d
j

j

 
2

2

d

inv

i (E )

(3.16)

ここでS j ：粒子 j のスピン
m j ：粒子 j の質量

 inv ( ) ：逆反応断面積、すなわち、残留核にエネルギーε の粒子 j が衝突して複合核が形成される反応に対する断面積

 i (E ) ,  d (E ) ：親核、娘核の準位密度
逆反応断面積は次式で与えられる。
 inv ( )   g cn 1  b  :中性子 g c j 1  V  



:荷電粒子

 g 1    (3.17)

ここで  g  R 2 は半径 R の核に対する幾何学的な断面積、 V  Z j Z d e 2 Rc はクーロン
障壁、cn および cj は定数である。この計算モデルでの逆反応断面積に関するパラメー
タは、n, p, d. t. 3He, α については Dostrovsky らによるパラメータを、その他の粒子については Matsuse らによるパラメータを用いている。
準位密度についてはフェルミガスモデルを用いることで、次式で表される。

 (E) e 2 a  E  
12 a 1 4 E   5 4

E  Ex(3.18)

ここでa  Ad 8 ：準位密度パラメータ

 ：娘核の対エネルギーである。Ex は Gilbert と Cameron によって定義された次式が用いられている。
Ex  U x   ,U x  2.5  150 Ad(3.19)
一方 E が Ex 以下のとき、次式の経験式が用いられる。

 (E)  1
12 T

e  E  E0  T

E  Ex

ここで13(3.20)

1 T  a U x  1.5 U x

E0  E x  T log T  0.25 log a  1.25 log U x  2 aU x

である。
粒子 j の全放出確率は(2.14)式をクーロン障壁 V から（ E  Q ）まで積分することで
得られる。

j 

g j g  (E) 
i

E Q

V




 1   ( E  Q   )d d(3.21)

上 式 は E  Q  V  Ex の と き 解 析 的 に 計 算 す る こ と が 可 能 で あ る 。 ま た 、E  Q  V  E x のときは近似式が用いられている。
放出粒子の種類は次式の確率分布で得られる。

p j  j  j j
放出粒子と娘核の運動エネルギーは(3.22)
式によって選ばれる。その運動の角度分布は、重心系で等方性の乱数で決定する。また、娘核の励起エネルギーは次式で計算される。

Ed E  Q  

3.3.4(3.23)

イベントジェネレーターモード（e-mode）

核反応を模擬する方法としては、上記のようなモデル計算のほかに、核データを用
いて Boltzmann 方程式を解く方法がある[26,27]。核データとは原子核の反応、壊変およ
びその構造に関わるデータの総称である。
しかし、核データを用いた計算で求められる量は平均値であり、統計的な揺らぎは求めることが出来ない。これは Boltzmann 方程式が位相空間での一体の分布関数の時
間発展を記述する方程式であり、また核データにも生成二次粒子の相関に関する情報が入っていないため、2 つ以上の粒子を放出するチャンネルが開いた際、エネルギー
および運動量保存則が成立しないことが原因である。また、核データを用いて中性子によるエネルギー付与量を計算する際に使用するカーマ近似は、エネルギー付与量の平均値を導出するための近似であり、それを用いて統計的な揺らぎを計算することはできない。
そのような揺らぎを計算するためには、特定の二次粒子の断面積を独立に計算するだけでなく、ひとつの核反応で生成される全ての二次粒子と反跳核を、エネルギーおよび運動量保存則の制約のもとでモンテカルロ的に生成する、すなわち、核反応から
14

の全ての放出粒子と反跳核のエネルギーと運動量を記述できるモデルを用いる必要がある。このようなモデルをイベントジェネレーターと呼ぶ。上記の INC モデルなど
もイベントジェネレーターである。
PHITS には、20MeV 以下の中性子反応に関する核データを用いた計算をイベント
ジェネレーターとして模擬できる計算オプションが搭載されている。この計算オプシ
ョンのことを e-mode と呼ぶ。

3.4 阻止能計算コード
PHITS には荷電粒子の阻止能の計算を行うために SPAR コード[28]、ATIMA コード[29]が搭載されている。両者とも Bethe-Bloch の式[18]に基づいて平均の阻止能の計算を行うコードである。Bethe-Bloch の式は次式で与えられる。
dE
Z z 2  1  2me  2 2 v 2Wmax    2  L 4N a re2 me c 2 
ln 
2 
2
dx
A  2 I(3.24)

ここで、

N a ：アボガドロ数 = 6.022×1023 [mol-1]

re ：古典電子半径 = 2.817×10-13 [cm]me ：電子質量 ：吸収物質の密度Z ：吸収物質の原子番号

A ：吸収物質の原子質量z ：入射粒子の素電荷基準の電荷

I ：平均励起ポテンシャル v c 11  2

Wmax ：一回の衝突における最大エネルギー伝達、

Wmax 

s  me M （ M ：入射粒子の質量）  

2me c 2 2
1  2s 1   2  s 2

また L は補正頄を表す。一般的には密度効果からの補正  と殻補正 C を考慮して計算を行う。このときの L は次式となる。

L  


2



C
Z

(3.25)SPAR コードでは、入射粒子のエネルギーおよび核電荷によって異なる 3 つの計算方法を用いている。高エネルギー領域（   0.04 z 2 3 ）に関しては(3.25)式の補正を用15

いた Bethe-Bloch の式が採用されている。中エネルギー領域（ 0.0046 z1 3    0.04 z 2 3 ）に関しても同様に Bethe-Bloch の式を用いるが、電子捕獲による入射粒子の正味の電荷の低下が無視できなくなるため、入射粒子の核電荷 z の替りに次式で計算される有効核電荷 z eff を使用している[30]。

zeff  z 1  exp 125z 2 3 (3.26)

図 3.3 に示すように Bethe-Bloch の式は粒子のエネルギーが低い領域では精度が悪い
ため[31]、低エネルギー領域（ 0    0.0046 z1 3 ）では J.Lindhard, M.Sharff, H.E.Schiottらによる理論（LSS 理論）[32]が採用されている。
ATIMA コードもまた入射粒子のエネルギーによって異なる 3 つの計算方法が適用されている。30MeV/amu 以上の高エネルギー領域では(3.25)式の補正に加えて Barkas
補正、LS 補正を用いた Bethe-Bloch の式で計算が行われる[33]。10MeV/amu 以下の低エネルギー領域では、実験値を基にした半経験式で計算を行う SRIM[34]を採用している。またこれらの中間のエネルギー領域では両者の補間式を使用している。

Stopping Power [keV/(g/cm2)]

1.5
Exp. data
Bethe formula
LSS theory
1

0.5

図 3.3

0

1

10

100

1000

Proton Energy [keV]104

シリコンへの陽子入射に関する阻止能の比較[31]163.5 ストラグリング分布関数
PHITS にはストラグリングを考慮するために Gauss 分布、Vavilov 分布、Landau 分
布を適宜用いるオプションが搭載されている。厚さがそれぞれ 1.0μm、10μm 、1.1mmの Si セルへの 20MeV の陽子入射によるエネルギー付与分布について、ストラグリン
グを考慮した PHITS 計算の結果を図 3.4～3.6 に示す。
適用される分布関数は次式で計算される  によって決定される。

  0.01 : Landau分布 0.01    10 : Vavilov分布
 10  : Gauss分布   Wmax

(3.26)

ここで  2N a re2 me c 2 

Z z2
x
A 2

Wmax 

2me c 2 2
1  2s 1   2  s 2

 はターゲット材質：  , Z , A 、ターゲット厚さ： x 、入射粒子の種類： z , s および入射エネルギー： E に依存することが分かる。例えばターゲットを厚さ 500nm のシ
リコンと仮定したとき、適用される分布関数は図 3.7 に示すように入射粒子の種類とエネルギーによって決まる。
0

10

Landauf(E) [1/keV]

10-1

proton@20MeV
into 1.0m Si( = 0.009)

-2

10

10-3
-4

10

10-5

E
0

10

20

30

図 3.4

Landau 分布

40

Deposit Energy [keV]

17

50

0

10

Vavilov-1f(E) [1/keV]

10

proton@20MeV
into 10m Si( = 0.09)

-2

10

10-3
-4

10

10-5

E
0

50

100

150

Deposit Energy [keV]
図 3.5

Vavilov 分布

0

10

Gauss

f(E) [1/keV]

10-1

proton@20MeV
into 1.1mm Si( = 10.1)

-2

10

10-3
-4

10

10-5

E
0

4

4

2000 4000 6000 8000 1 10 1.2 10

Deposit Energy [keV]
図 3.6

Gauss 分布18

1000

Projectile Energy [MeV]Landau100

Vavilov

10

Gauss

1

Si target thickness = 500 [nm]
0.1

1
図 3.7

3

5

7

9

11

Atomic Number

適用されるストラグリング分布関数

19

13

15

4. 実デバイス体系を模擬した計算結果および考察
4.1 計算体系および計算条件
計算体系の模式図を図 4.1 に示す。本計算では、設計ルール 45nm のデバイスを模擬するために有感領域の面積を 0.1μm×0.1μm とし、有感領域の深さを 0.5μm と定義した[35]。また統計精度を上げるために、反応領域中に有感領域を 1μm 間隔で二次元格子状に配置した。Si セルの上部には SiO2 で構成された厚さ 0.35μm の絶縁層およびCu と SiO2 を 1:1 の比で混合した金属層を置いた。金属層の厚さについては文献[36]に記載されている配線幅およびアスペクト比より推定した結果、十分な厚さを取るために 3.0μm とした。金属層の上部にはパッケージ層を配置した。パッケージ材料はエ
ポキシ樹脂と SiO2 や Al2O5 などの充填剤で構成されているが、充填剤が全体の 50～
90 重量%を占めているため、ここでは SiO2 のみ考慮した。反応領域のサイズおよびパッケージ厚さの決定については次節で述べる。
入射中性子スペクトルに関しては EXPACS[3] を用いて計算した中性子スペクトル（図 1.1-b 参照）を使用した。また中性子はパッケージ上部から垂直に入射した。
以下の計算では、核反応モデルは QMD+GEM、20MeV 以下の中性子核反応については核データとして JENDL-3.3[37]を用いた e-mode 計算オプションを用いた。また阻止能計算は ATIMA コードを採用した。
臨界電荷量 Qc はデバイスの種類や設計ルールなどによって異なるが、最近のデバ
イスでは 1fC を下回るという報告がある[38]。そこで本研究では Qc = 1fC 前後について議論を行った。

図 4.1

計算体系の模式図20

表 4.1

配線幅、アスペクト比および配線厚さ（2005 ITRS JEITA 和訳 配線 10 頁 表 80a[36]）配線名

配線幅 [nm]

アスペクト比配線厚さ [nm]Metal 1

90

1.8

162

Intermediate (via)

90

1.6

144

Intermediate

90

1.8

162

Semi-global (via)

180

2.2

396

Semi-global

180

2.4

432

Global (via)

200

2.2

440

Global

200

2.4

480

Total

図 4.2

2216

デバイスの断面図および配線名（2005 ITRS JEITA 和訳 配線 7 頁 図 71[36]）2110Critical Charge [fC]8
180 nm6
4

130 nm
90 nm

2

65 nm
45 nm

0
200
図 4.3

160

120

80

design rule [nm]40設計ルールと臨界電荷量の相関（SRAM）[38]

シミュレーションの都合上、有感領域 1 つ 1 つを区別してエネルギー付与分布を出力することが困難であった。そのため、1 つの荷電粒子が複数の有感領域を通過した
場合、最初に通過した有感領域へのエネルギー付与のみ考慮して計算を行った。
上記の条件の妥当性について確認した結果、1 つの荷電粒子が複数の有感領域を通過するようなイベントは全体の数%程度しか見られなかった。またこのようなイベン
トによる影響について確認するため、(a) 全ての有感領域を同一として取扱った計算と、(b) 有感領域に 1～9 までのセル番号を割り振ってそれぞれの領域を区別して取扱った計算について比較を行った。その結果、図 4.5 に示すように両者の計算結果はほぼ一致したため、前述の仮定は妥当であると結論付けた。

22

図 4.4 有感領域の取扱い方

BGR [cm2/m3]

10-10

全ての有感領域を同一として取扱った場合
9つの有感領域を区別して取扱った場合

-11

10

-12

10

0

0.5

1

1.5

Critical Charge [fC]
図 4.5

計算条件の妥当性の確認

232

4.2 計算体系サイズに関する計算結果の収束性の確認
臨界電荷量の低下に伴い、阻止能の低い軽イオンによってもソフトエラーを引き起
こすことが可能となる。軽イオンは飛程が長いため、有感領域から離れた場所で生成された場合でも有感領域へ到達することができる。そのため、軽イオンによる影響を正確に評価するためには計算領域を大きく取る必要がある。本節では計算体系のサイズによる計算結果の収束性について確認を行った。ここで反応領域の深さを x, 一辺の長さを y, 上部パッケージ層の厚さを z としている。
図 4.6～4.8 にそれぞれ反応領域の厚さ、面積サイズおよびパッケージ層の厚さに関する収束性の確認結果を示す。反応領域の厚さに関しては 500μm で計算結果が収束していることがわかる。また反応領域の面積に関しては多尐の揺らぎがあるものの、一辺 1000μm で概ね収束していると判断した。パッケージ厚さに関しては厚さ 500μmで収束性を確認した。
以上の結果より、本章では反応領域 1.0×1.0×0.5mm、パッケージ層厚さ 0.5mm と定義して計算を行った。

BGR [cm2/m3]

10-10

RV depth: x(y = 1000m, z = 0m)
5m
10m
100m
500m
1000m

10-11

10-12

0

0.5

1

1.5

Critical Charge [fC]
図 4.6

反応領域の厚さに関する収束性の確認
242

BGR [cm2/m3]

10-10

RV side: y(x = 500m, z = 0m)
1m
10m
100m
1000m2000m
4000m

10-11

10-12

0

0.5

1

1.5

2

Critical Charge [fC]
反応領域の面積サイズに関する収束性の確認

10-10package thickness: z(x = 500m, y = 1000m)
0m
10m
100m
500m
1000mBGR [cm2/m3]

図 4.7

-11

10

-12

10

0

0.5

1

1.5

Critical Charge [fC]
図 4.8

上部パッケージ層の厚さに関する収束性の確認
252

4.3 中性子の入射エネルギーに関する確認
二次宇宙線中性子のエネルギースペクトルは広範囲に渡っているが、シミュレーシ
ョンにおいては低エネルギー中性子の存在は無視される場合が尐なくない。そこで本節ではシミュレーションにおいて考慮すべき入射中性子のエネルギー範囲を明確にするため、表 4.2 のように入射エネルギーを 6 群に分割して計算を行い、各エネルギ
ー群の寄与について確認した。
計算結果を図 4.9、4.10 に示す。ここで Total とは各エネルギー群を用いて計算した
電荷付与分布 N i (q) と、各エネルギー群の入射中性子フラックス i と全入射中性子フラックス  の比を用いて、次式で計算した値である。
6

total   N i ( q) 
i 1

i


(4.1)

1MeV 以下の中性子は、入射する割合は全体に対して半数を占めるものの、BGR へ
ほとんど寄与しないことがわかった。また 3.5GeV 以上の入射中性子による影響もほ
とんどないことが確認された。
以上の結果より、本計算では 1MeV から 3.5GeV までの入射中性子スペクトルを用いてシミュレーションを行うこととした。

表 4.2
入射中性子群Group 1

入射中性子エネルギー群
エネルギー領域
0.01eV ~ 1eV

入射割合
0.202

Group 2

1eV ~ 1MeV

0.315

Group 3

1MeV ~ 20MeV

0.231

Group 420MeV ~ 200MeV

0.179

Group 5200MeV ~ 3.5GeV

0.072

Group 6

3.5GeV ~ 90GeV

0.001

26

BGR [cm2/m3]

10-10

Total
Group 1
Group 2
Group 3
Group 4
Group 5
Group 6

10-11

10-12

10-13

10-14

0

0.5

1

1.5

2

Critical Charge [fC]
図 4.9

各入射エネルギー群に関する BGR への寄与

Group 1
Group 4

Group 2
Group 5

Group 3
Group 6

1

Fraction

0.8
0.6
0.4
0.2
0

0

図 4.10

0.5

1.0

1.5

Critical Charge [fC]2.0

各入射エネルギー群に関する BGR への寄与度

27

4.4 二次イオン種の寄与
入射中性子とデバイス材料との核反応によって様々な二次イオンが生成される。本節ではいずれの二次イオン種がソフトエラーの主要因となるか調べた。
計算結果を図 4.11、4.12 に示す。Qc = 1.2fC 前後で He イオンの寄与が残留核（A > 4のイオン）と同等であることが分かる。また臨界電荷量が低下するにつれて He イオ
ンが主要因となることが分かる。さらに、Qc = 0.2fC 前後まで低下したとき、ソフト
エラーの主要因が He イオンから陽子へ替わることが予測される。今後は臨界電荷量が低下すると考えられているため、核反応からの軽イオン生成について精度よくモデル計算を行うことが重要であることがわかった。

BGR [cm2/m3]

10-10
total
H
He
A>4

-11

10

-12

10

10-13

0

0.5

1

1.5

Critical Charge [fC]
図 4.11

各二次イオン種による BGR への寄与

28

2

H
F - Si1

He
P - CuLi - O
othersFraction

0.8
0.6
0.4
0.2
0

0
図 4.12

0.5

1.0

1.5

Critical Charge [fC]

各二次イオン種による BGR への寄与度

29

2.0

4.5 核反応モデルの違いに関する計算結果の比較
PHITS には核反応の計算を行うモデルが複数搭載されている。本節では核反応モデ
ルの違いが計算結果へ与える影響について調査を行った結果を示す。なお、核反応モデルの妥当性の検証は後述する 6.1 節で行っている。
比較を行う核反応の動的過程モデルとして、Bertini による INC モデル、QMD モデ
ルを用いた。蒸発過程に関しては GEM モデルを採用した。また 20MeV 以下の中性子
核反応について、核データとして JENDL-3.3 を用いた e-mode オプションを適用した
場合についても比較を行った。
計算結果を図 4.13 に示す。QMD+GEM と e-mode を併用した計算結果と比べて、INC＋GEM による結果は 30～40%程度低い値となった。また INC+GEM と e-mode を併用
した結果は 20～30%程度低い値となった。
以上の結果より、核反応モデルの選択もまた計算精度へ大きく影響を及ぼすことが確認された。

BGR [cm2/m3]

10-10

INC
INC & e-mode
QMD
QMD & e-mode

-11

10

-12

10

0

0.5

1

1.5

Critical Charge [fC]
図 4.13

核反応モデルの違いに関する BGR の比較
30

2

4.6 阻止能計算コードの違いに関する計算結果の比較
PHITS には阻止能計算を行うコードとして SPAR コードおよび ATIMA コードが搭載されている。本節では阻止能計算コードの違いが計算結果へどのような影響を与えるかについて調査を行った。なお、阻止能計算コードの妥当性の確認は後述する 6.2
節で行っている。
計算結果を図 4.14 に示す。Qc = 1fC 以上の領域では、SPAR を用いた計算結果がATIMA を用いた計算結果に対して 5～15%程度低い値を示した。イオン種の寄与について比較すると、これは He イオンによる高電荷付与イベントが減尐したためである
ことがわかる。一方、Qc = 1fC 以下の領域では両者の結果にはほとんど違いは見られなかった。
現在のデバイスに関する臨界電荷量は 1fC 以下であり、次世代のデバイスでは臨界
電荷量はさらに下がると予測されているため、阻止能計算の違いによる影響はそれほど大きくないと考えられる。

-10

BGR [cm2/m3]

10

ATIMA
SPAR
H
He
A>4

-11

10

10-12

10-13

0

0.5

1

1.5

Critical Charge [fC]
図 4.14

阻止能計算コードの違いに関する BGR の比較

31

2

4.7 デバイスの材質の違いに関する計算結果の比較
実際のデバイスには様々な材料が用いられている。その一方で、ソフトエラーシミュレーションを行う際にデバイスの材質を Si のみとして計算を行う場合もある。そこで本節ではデバイスの材質の違いがどのような影響を与えるか調べるため、以下の4 つの計算体系を用いて計算を行い、BGR について比較を行った。図 4.15 に計算体
系の模式図を示す。
(a)1×1×1mm の Si 単一セル中に有感領域を配置した計算体系(b)

有感領域上部に存在する層を全て SiO2 とした計算体系(c)4.1 節で定義した計算体系(d)

計算体系(c)に対し、有感領域の上部にタングステン配線を配置した計算体系

ここでは 4.1 節で定義した計算体系(c)を基準として比較を行った。それぞれの計算体系に関する BGR の比較結果を図 4.16 に示す。
Si のみを考慮した計算体系(a)と、計算体系(c)による結果を比較すると Qc = 0.3～
1.5fC の領域で BGR が 10～20%程度低い値となった。図 4.17 を見ると、両者の違い
は He イオンによる寄与が低下したため生じたと考えられる。
基準となる計算体系から Cu を取り除いた計算体系(b)と計算体系(c)による計算結果を比較すると、(b)の BGR の方が Qc = 0.5fC 以上の領域で BGR が 10～20%程度高い値となった。図 4.18 を見ると、He イオンより重いイオンによる寄与が増加していることが確認できる。これは Si 近傍のイオンや O 近傍のイオンは、より重い Cu 近傍のイオンと比べて飛程が長いため、有感領域へ到達するものが増加したためであると予想される。
タングステン配線を考慮した計算体系(d)と計算体系(c)を比較すると、BGR の値に変化が見られなかった。デバイス内に存在するタングステン配線の量が極微量である
ため、シミュレーションへの影響がほとんど表れなかったと考えられる。
以上の結果より、デバイスを構成する材料の違いによって計算結果が 10～20%程度変動することが確認された。これより、ある程度現実に則した計算体系を定義する必要があると考える。一方で、タングステン配線の存在による影響は見られなかったた
め、デバイスの微細な構造を考慮する必要はないことが示唆された。

32

図 4.15

それぞれの計算体系の模式図

BGR [cm2/m3]

10-10(a) only Si(b) without Cu(c) with Cu(d) with W

-11

10

-12

10

0

0.5

1

1.5

Critical Charge [fC]
図 4.16

それぞれの計算体系に関する BGR の比較
33

2

-10

BGR [cm2/m3]

10(a) only Si(c) with Cu
H
He
A>4

-11

10

10-12

10-13

0

0.5

1

1.5

2

Critical Charge [fC]
図 4.17

計算体系(a), (c)に関する BGR への二次イオンの寄与の比較-10

BGR [cm2/m3]

10(b) without Cu(c) with Cu
H
He
A>4

-11

10

10-12

10-13

0

0.5

1

1.5

2

Critical Charge [fC]
図 4.18計算体系(b), (c)に関する BGR への二次イオンの寄与の比較
34

4.8 有感領域サイズを変更したときの計算結果の比較
デバイスの微細化・高集積化によるソフトエラー率の変化を予測するために、文献[38]を参考にし、有感領域サイズを 0.137×0.137×0.5μm としたもの（設計ルール 90nmに相当）、また設計ルールが 45nm から更に縮小したものを模擬するために有感領域サ
イズを 0.07×0.07×0.5μm としたものについて計算を行った。
計算結果を図 4.19 に示す。ここでは電荷収集効率 C  1 と仮定し、ソフトエラー発
生確率 SER を入射中性子フラックス n で割ったもので比較した。また設計ルール
90nm、45nm に関しては図 4.3 で示した臨界電荷量について黒点をプロットした。
ソフトエラー率の変化は臨界電荷量に強く依存するため一概には言えないが、設計ルールが 90nm から 45nm へ縮小したときのソフトエラー率の増加は 2 倍以内であると考えられる。また更に微細化が進んだとしても、1 ビット当たりのソフトエラー率の変化は 2 倍以内に収まると予測される。
一方で 1 デバイス当たりのソフトエラー率は今後増加すると考えられる。例えば設計ルールが 90nm から 45nm へと 50%縮小したとき、1 ビット当たりの面積は 25%の縮小となる。よって単位面積当たりのビット数は単純に考えて 4 倍となり、ソフトエ
ラー率も更に 4 倍増加すると考えられる。

-120.137×0.137×0.5m (90nm)
0.1×0.1×0.5m (45nm)
0.07×0.07×0.5m

SER/n [cm2/bit]

10

-13

10

-14

10

0

0.5

1

1.5

2

Critical Charge [fC]
図 4.19

有感領域サイズの違いに関するソフトエラー発生確率の比較
35

5 電荷付与過程に関わる物理現象の影響調査および考察
5.1 ストラグリングによる影響の調査結果
4 章で用いた計算体系について、PHITS 内蔵のストラグリング分布関数を使用して計算を行った。計算結果を図 5.1 に示す。
ストラグリングを考慮した結果、陽子が主要因の領域となっている領域で BGR が増加した。一方で、He イオンが主要因となる領域および残留核が主要因となる領域
では計算結果にほとんど変化が見られなかった。
この原因を探るため、有感領域へのエネルギー付与に関与した二次イオンの初期エ
ネルギー分布を出力し、適用されるストラグリングの分布関数について調査を行った。
ここで有感領域へのエネルギー付与に関与したイオンとは、外部から有感領域へ入射するイオンと、有感領域内で生成されたイオンを加えたものである。（図 5.2）
図 5.3 に陽子、He イオン、O イオン、Si イオンに関する初期エネルギー分布を示す。
また図 5.4 に対して、それぞれのイオンの初期エネルギーの平均値、上限値、下限値をプロットしたものを図 5.4 に示す。

-10

BGR [cm2/m3]

10

without straggling
with straggling
H
He
A>4

-11

10

10-12

10-13

0

0.5

1

1.5

Critical Charge [fC]
図 5.1

ストラグリングを考慮したときの計算結果の変化
362

有感領域へのエネルギー付与に関与したイオン

10-5

10-6

Number of Particles [1/source]

Number of Projectiles [1/source]

図 5.2

H
-6

10

E = 44 [MeV]

10-7
10-8
-9

10

10-10

0

200

400

600

800

He
10-7

10-8

10-9

10-10

1000

E = 6.2 [MeV]

0

10-7

O
E = 0.5 [MeV]
-8

10

-9

10

10-10

0

0.5

1

1.5

2

2.5

3

Energy [MeV]

図 5.3

5

10

15

20

25

30

Particle Energy [MeV]

Number of Projectiles [1/source]

Number of Projectiles [1/source]

Projectile Energy [MeV]
10-7

Si
E = 0.15 [MeV]
-8

10

-9

10

10-10

0

0.5

1

1.5

2

2.5

Energy [MeV]

有感領域へのエネルギー付与に関与した各イオンの初期エネルギー分布
37

3

1000

Projectile Energy [MeV]Landau100

Si thickness= 500 [nm]Vavilov

10

Gauss

1

0.1

図 5.4

1

3

5

7

9

11

Atomic Number

13

15

各イオンへ適用されるストラグリング分布関数

図 5.4 を見ると、陽子に関しては主に Landau 分布が適用されており、これによって高エネルギー付与イベントが増加した。一方 O イオンや Si イオンに関しては揺らぎの小さい、平均付与エネルギー周りで左右対称な分布をしている Gauss 分布が適用されている。また He イオンに関しては Gauss 分布に近い形の Vavilov 分布が適用されていることが分かる。このためストラグリングを考慮したとしても高エネルギー付与イ
ベントの増加が見られなかったと考えられる。
以上の結果より、臨界電荷量が 1fC 前後の領域では He イオンが主要因であるため、ストラグリングによる影響は現在のところ顕在化していないと予想される。しかし臨
界電荷量が低下して主要因が陽子となったとき、ストラグリングによる影響を考慮する必要がある。

38

5.2 イオントラック構造による影響の検討
本節では、二次イオントラック内部の初期電荷の動径方向分布について検討を行った。二次イオンによって生成される二次電子（デルタ線、オージェ電子、プラズモン崩壊による生成電子など）の輸送をモンテカルロ法で計算した結果から、初期エネル
ギー付与の動径方向分布が次式で与えられている[39]。
C1 (T ) f ( r )r 2.308 (0.16 nm  r  10 nm)
2
D( r )  Z eff
 C 2 (T )r 2.055(10 nm  r  100 nm)
 C (T )r (100 nm  r  1000 nm)
 3

(5.1)

ここで、

T ：イオンの運動エネルギー [MeV/amu]
Z eff ：有効電荷、

r ：動径方向の位置 [um]Z eff  Z  1  exp 125Z 2 3  

C1 (T )  1.421  105 T 0.913

C2 (T )  9.137  10 4 T 0.9404

C3 (T )  1.402  10 6 T 0.887

f ( r )  1.124  0.223r  0.03656 r 2  0.00156 r 3

  2.7T 0.005

T T  1876.61 2
T  938.3

上式を用いて、Si への 44MeV 陽子入射によって形成されるイオントラック構造に関して計算を行った結果を図 5.5 に示す。また半径 r の微小厚円柱内に存在する電荷の、全体に占める割合を次式で計算した結果を図 5.6 に示す。
r

f (r )  2rD(r )dr

0.16


(5.2)

 2rD(r )dr0.16

図 5.6 を見ると、半径 100nm 以内に存在する電荷量は全体の 85～95%程度であることがわかる。最近のデバイスのセルサイズもまた 100nm オーダーである。そのため、
イオントラック構造によるシミュレーションへの影響があると考える。今後イオント
ラック構造を考慮したシミュレーションを行い、計算結果への影響について定量的に評価する必要がある。

39

Si への 44MeV 陽子入射によって形成されるイオントラック構造

1
0.8
0.6

f (r)

図 5.5

0.4
proton 44MeV
He 6.2MeV
O 0.5MeV
Si 0.15MeV

0.2
0

0

図 5.6

50

100

150

radial distance [nm]

動径方向に関する線電荷密度の割合
40200

6 各計算モデルの妥当性の確認
6.1 核反応モデルに関する実験値との比較および考察
4.5 節では、核反応モデルの違いによる計算結果への影響は小さくないことが確認された。また 4.4 節ではソフトエラーの主要因が軽イオンであることが示された。そこで本章では、軽イオン生成に関するそれぞれの核反応モデルの精度について検討するため、中性子と半導体デバイスの主要構成材料である Si との核反応からの陽子、アルファ粒子生成に関する実験データとの比較を行った。
比較する核反応モデルとして、Bertini による INC モデル、QMD モデルを用いた。
また 20MeV 以下の中性子核反応に関して、核データとして JENDL-3.3 を用いた
e-mode 計算オプションについても計算、比較を行っている。
表 6.1 に比較を行った実験データの一覧を示す。比較に用いた実験値は生成断面積（CS）、エネルギー微分断面積（DX）および二重微分断面積（DDX）である。

表 6.1

比較を行った実験データ一覧[40-45]

生成イオン種入射エネルギー [MeV]

測定者

測定年

比較内容

proton14.4

F.L.Hassler+

1962

CS

14.1

Y.Terada+

2002

DX, DDX

5.0~50.0

F.B.Bateman+

1999

CS, DX, DDX

62.7

S.Benck+

2002

CS, DX, DDX

95.6

U.Tippawan+

2004

CS, DX, DDX

14.8

D.W.Kneff+

1986

CS

5.0~58.0

F.B.Bateman+

1999

CS, DX, DDX

62.7

S.Benck+

2002

CS, DX, DDX

95.6

U.Tippawan+

2004

CS, DX, DDX

4

He

41

6.1.1Si(n,p)X に関する実験値との比較

n+Si 反応による陽子生成に関して、核反応モデルによる計算結果と実験値とを比較した結果を図 6.1～6.9 に示す。
INC+GEM による結果は CS に関して非常に良い一致を示した。また DX に関しても良い一致を示したが、高エネルギー陽子放出に関して実験値を過大評価する傾向が見られた。また DDX に関して比較すると、陽子の放出角度が後方になるに従って実験値を過小評価する傾向が見られ、実験値の傾向と異なる結果を示した。
QMD+GEM による結果は CS、DX、DDX ともに良好な一致が得られた。また低エ
ネルギー中性子核反応に関してもある程度の一致を示すことを確認した。
e-mode による結果は CS に関して 10MeV 以下の中性子入射の場合、実験値を過大評価する傾向を示した。DX、DDX に関しては QMD+GEM による結果と同程度の一致を示した。
以 上 の結 果よ り、 陽子 生成 に 関し ては QMD+GEM あ る い は QMD+GEM にJENDL-3.3 を用いた e-mode 計算を併用するのが妥当であると結論づけた。

1

10

proton0

Cross Section (b)

10

-1

10

F.B.Bateman+ 1999
F.L.Hassler+ 1962
S.Benck+ 2002
U.Tippawan+ 2004
QMD
INC
e-mode

-2

10

10-3
10-4

10

Incident Energy (MeV)
図 6.1

n+Si 反応からの陽子生成に関する CS 比較

42

100

4

10

Ein= 14MeV

3

d/dE (mb/MeV)

10

2

10

101
100
10-1

図 6.2

Y.Terada+ 2002
QMD
INC
e-mode0

5

10Proton Energy (MeV)

15

n+Si 反応からの陽子生成に関する DX 比較（Ein = 14MeV）

104

Ein= 38MeV

d/dE (mb/MeV)

103

F.B.Bateman+ 1999
QMD
INC

2

10

1

10

0

10

-1

10

-2

10

-3

10

0

10

20

30

40

50

60

Proton Energy (MeV)
図 6.3

n+Si 反応からの陽子生成に関する DX 比較（Ein = 38MeV）
43

4

10

Ein= 63MeV

d/dE (mb/MeV)

103

S.Benck+ 2002
QMD
INC

2

10

101
100
10-1
10-2
10-3

図 6.4

0

10

20

30

40

50

60

Proton Energy (MeV)70

80

n+Si 反応からの陽子生成に関する DX 比較（Ein = 63MeV）

104

Ein= 96MeV

d/dE (mb/MeV)

103

U.Tippawan+ 2004
QMD
INC

2

10

1

10

0

10

-1

10

-2

10

-3

10

図 6.5

0

20

40

60

80

Proton Energy (MeV)

100

n+Si 反応からの陽子生成に関する DX 比較（Ein = 96MeV）
44

図 6.6n+Si 反応からの陽子生成に関する DDX 比較（Ein = 14MeV）

図 6.7

n+Si 反応からの陽子生成に関する DDX 比較（Ein = 38MeV）
45

図 6.8n+Si 反応からの陽子生成に関する DDX 比較（Ein = 63MeV）

図 6.9

n+Si 反応からの陽子生成に関する DDX 比較（Ein = 96MeV）
46

6.1.2Si(n,4He)X に関する実験値との比較

n+Si 反応によるアルファ粒子生成に関して、核反応モデルによる計算結果と実験値について比較した結果を図 6.10～6.18 に示す。
INC+GEM による結果は CS、DX、DDX ともに実験値を大きく過小評価する結果となった。
QMD+GEM による結果は CS に関して 25MeV 以上の領域で良好な一致を示した。
DX、DDX に関して、高エネルギー領域では INC+GEM と同様に実験値を大きく下回
ったが、低エネルギー領域に関しては実験値と同等の結果を示した。低エネルギー4Heは核反応の蒸発過程から生成されたもので、QMD は蒸発過程に至るまでの過程を精度よく記述できるため、実験と良い一致を示したと考える。
e-mode による結果は CS に関して 9MeV 以上の中性子入射に関してよい一致を示した。DX に関しても実験値との良好な一致を示したが、DDX に関しては実験値を上回る結果となった。
以上の結果より、4He 生成に関しても QMD+GEM と JENDL-3.3 を用いた e-mode 計
算の併用が妥当であると考えられる。しかし高エネルギー4He 放出が行われていない
点や中性子の入射エネルギーが 8MeV 未満での 4He 生成について記述できていない点など、改善すべき点がいくつかあることが確認された。
1

10

4

He

0

Cross Section (b)

10

-1

10

F.B.Bateman+ 1999
D.W.Kneff+ 1986
S.Benck+ 2002
U.Tippawan+ 2004
QMD
INC
e-mode

10-2
10-3
10-4

10

Incident Energy (MeV)
図 6.10100

n+Si 反応からの 4He 生成に関する CS 比較
47

4

10

Ein= 15MeV

d/dE (mb/MeV)

103

F.B.Bateman+ 1999
QMD
INC
e-mode

2

10

101
100
10-1
10-2
10-3

図 6.11

0

5

4

10

15

20

He Energy (MeV)

25

30

n+Si 反応からの 4He 生成に関する DX 比較（Ein = 15MeV）

104

Ein= 38MeV

d/dE (mb/MeV)

103

F.B.Bateman+ 1999
QMD
INC

2

10

1

10

0

10

-1

10

-2

10

-3

10

図 6.12

0

10

4

20

30

40

50

He Energy (MeV)

n+Si 反応からの 4He 生成に関する DX 比較（Ein = 38MeV）
48

4

10

Ein= 63MeV

d/dE (mb/MeV)

103

S.Benck+ 2002
QMD
INC

2

10

101
100
10-1
10-2
10-3

図 6.13

0

10

20
4

30

40

50

60

70

80

He Energy (MeV)

n+Si 反応からの 4He 生成に関する DX 比較（Ein = 63MeV）

104

Ein= 96MeV

d/dE (mb/MeV)

103

U.Tippawan+ 2004
QMD
INC

2

10

1

10

0

10

-1

10

-2

10

-3

10

図 6.14

0

20

4

40

60

He Energy (MeV)

80

100

n+Si 反応からの 4He 生成に関する DX 比較（Ein = 96MeV）
49

図 6.15n+Si 反応からの 4He 生成に関する DDX 比較（Ein = 15MeV）

図 6.16n+Si 反応からの 4He 生成に関する DDX 比較（Ein = 38MeV）
50

図 6.17n+Si 反応からの 4He 生成に関する DDX 比較（Ein = 63MeV）

図 6.18n+Si 反応からの 4He 生成に関する DDX 比較（Ein = 96MeV）
51

6.2 阻止能計算コードと実験値との比較および考察
PHITS に内蔵されている SPAR コード、および ATIMA コードを用いて、Si セルへの各イオン入射に関する阻止能について計算を行い、実験値と比較を行った。また2008 年に更新された SRIM（Stopping and Range of Ion in Matter）コード[46]による計算結果も合わせてプロットした。比較結果を図 6.19～6.22 に示す。
SPAR による計算結果は重イオンに関しては実験値と良好な一致を示したが、陽子や He イオンなどの軽イオンに関しては 100keV 前後の低エネルギー領域で実験値を下回った。
ATIMA による計算結果はいずれのイオンに関しても良好な一致が得られた。また
ATIMA は 10MeV/amu 以下の領域で旧バージョンの SRIM を利用して計算を行っているため、SRIM-2008 による計算結果と近い値となっていることも確認された。
現在の段階でソフトエラーの主要因は He イオンであり、臨界電荷量の低下に従って陽子による寄与も増加すると考えられるため、軽イオンの阻止能に関して実験値と良好な一致を示した ATIMA を採用することが妥当であると考える。
また 4.6 節の図 4.14 で示したように、SPAR による計算結果は ATIMA による結果と比べて軽イオンによる高電荷付与イベントが尐なかった。これは図 6.19 および 6.20に示したように、Si への軽イオン入射の阻止能のピーク値が、ATIMA と比べて SPARによる結果が低い値となっているためであることもわかった。

52

H into Si

2

Stopping Power [keV/(g/cm )]1.5

図 6.19SH.Z.Izmailov (1980)
E.I.Sirotinin (1984)
P.Mertens (1988)
A.Ikeda (1996)

1

SRIM-2008
PHITS (ATIMA)PHITS (SPAR)

0.5

0

1

10

100

1000

Incident Ion Energy [keV]104

Si への陽子入射に関する阻止能計算モデルと実験値[47-50]との比較

He into Si

2

Stopping Power [keV/(g/cm )]2.5

図 6.20

W.K.Lin (1973)
D.C.Santry (1980)
J.D.Pearce (1981)
D.C.Santry (1984)

2

SRIM-2008
PHITS (ATIMA)PHITS (SPAR)

1.5
1
0.5
0

1

10

100

1000

Incident Energy [keV]104

Si への He イオン入射に関する阻止能計算モデルと実験値[51-54]との比較
53

O into Si

2

Stopping Power [keV/(g/cm )]15

図 6.21

10

J.G.Kelley (1973)
D.C.Santry (1992)
W.Jiang (1999)
J. of ICRU (2005)

SRIM-2008
PHITS (ATIMA)PHITS (SPAR)

10

104

5

0

1

100

1000

105

Incident Energy [keV]106

Si への O イオン入射に関する阻止能計算モデルと実験値[55-58]との比較

Si into Si

2

Stopping Power [keV/(g/cm )]25

図 6.2220

J.Keinonen (1992)
H.J.Whitlow (2002)
J. of ICRU (2005)

SRIM-2008
PHITS (ATIMA)PHITS (SPAR)

10

104

15
10
5
0

1

100

1000

105

Incident Energy [keV]106

Si への Si イオン入射に関する阻止能計算モデルと実験値[58-60]との比較
54

7. まとめおよび結論
本研究では、ソフトエラー発生過程のうち、核反応過程および電荷付与過程において重要となる要因を解明するために実際の半導体デバイスを模擬した計算体系を用いてシミュレーションを行うと共に各計算モデルの精度について確認した。また二次イオンによる電荷付与量計算の精度向上を図るために、半導体デバイスの微細化に伴い影響が表面化すると予測される現象による、ソフトエラーシミュレーションへの影
響について調査を行った。
まず、実デバイス体系を模擬した計算体系への中性子入射に関するシミュレーショ
ン結果について、以下に整理する。


二次宇宙線中性子は広範囲のエネルギースペクトルを有しているが、ソフトエラーに関与するものは 1MeV から 3.5GeV の間の入射中性子がほとんどで、その寄与は全入射中性子に対して約 98%であった。



臨界電荷量が 1fC 程度の領域でソフトエラーの主要因は He イオンであり、今後
臨界電荷量が減尐するに従って陽子による寄与も現れてくると予測される。また
陽子や He イオンなどの軽イオンは飛程が長いため、その影響を正しく評価するためには尐なくとも一辺 1mm 程度の大きさを持つ計算体系が必要となる。



核反応モデルの違いによるシミュレーション結果への影響は大きいことが確認された。PHITS に搭載されている核反応モデルの精度について調査した結果、軽イオン生成に関しては QMD+GEM と、核データとして JENDL-3.3 を用いた e-mode
計算オプションを併用することが望ましいとわかった。しかし動的過程からの高エネルギーHe 放出や低エネルギー中性子入射反応からの He 生成について、本モデル計算では十分に記述できていないため、これらについて改善する必要がある。



デバイスの材質を変えてシミュレーションを行った結果、計算体系はある程度現実に則して定義する必要があるが、タングステン配線の存在など詳細な構造の違いによる影響は小さいことが判明した。



有感領域サイズを変更して計算した結果、設計ルールが 90nm から 45nm へ縮小
したとき、1 ビット当たりのソフトエラー発生確率の変化は 2 倍程度増加することがわかった。一方で 1 デバイス当たりのソフトエラー発生確率はデバイスの集積度に比例するため、4 倍程度上昇すると考えられる。
55

次に、半導体デバイスの微細化に伴い影響が現れると予測される現象として、ストラグリングとデルタ線によるイオントラック構造の形成による、ソフトエラーシミュ
レーションへの影響について調査した結果を以下に整理する。


ストラグリングによる影響を定量的に評価した結果、陽子がソフトエラーの主因となる領域で BGR の増加が見られた。将来さらに微細化が進み臨界電荷量が減尐し、陽子による影響が現れたとき、ストラグリングを考慮する必要があると考
えられる。



デルタ線によるイオントラック構造の形成について、動径方向の電荷密度分布の計算を行った結果、半径 100nm 程度までイオントラック構造を考慮する必要性が確認された。現在のデバイスセルサイズも 100nm 程度であり今後さらに縮小が進むと考えられるため、イオントラック構造による影響を考慮する必要があると予測される。
最後に今後の展望について述べる。今回は電荷収集過程について有感領域モデルを使用したが、回路シミュレーションなどを行うためには過渡電流応筓が必要となる。
過渡電流について計算するためには、核反応からの生成二次イオンによってデバイス内に付与される初期電荷分布を出力し、これを入力データとしてデバイスシミュレー
ションを行う必要がある。そのためには初期電荷分布を再現する三次元メッシュの生成、特に斜めイオン入射の取り扱いが重要となる。
共同研究の枠組みで、筑波大学の研究グループによって高エネルギー粒子の任意の侵入方向と二次イオントラックの微細構造に対応できるように、メッシュ生成プログ
ラムが新規開発された[61]。今後はこれを用いて PHITS と汎用三次元デバイスシミュ
レータ ENEXSS[62]（Environment for NEXt Simulation System）をリンクさせたソフトエ
ラーシミュレータを構築し、デバイスへの中性子入射による過渡電流応筓について計算を行う予定である。

56

参考文献
[1] T.Horikiri, 日経エレクトロニクス, No.903, 63 (2005).
[2] Y.Tosaka, 日経エレクトロニクス, No.903, 145 (2005).
[3] T.Sato, EXPACS ver.2.17 (2009); http://phits.jaea.go.jp/expacs/jpn.html[4] J.F.Ziegler et al., IBM J. Res. Develop., 40, No.1, 3 (1996).
[5] D.Binder et al., IEEE Trans. on Nucl. Sci., 22, 2675 (1975).
[6] T.C.May and M.H.Woods, IEEE Trans. Electron Devices, 26, 2 (1979).
[7] J.F.Ziegler and W.A.Lanford, SCIENCE, 206, 776 (1979).
[8] C.S.Guenzer et al., IEEE Trans. on Nucl. Sci., 26, 5048 (1979).
[9] J.F.Dicello, IEEE Trans. Nucl. Sci., 30, 4613 (1983).
[10] W.McKee et al., IEEE International Reliability Physics Symposium (IRPS) Proceedings,
1 (1996).
[11] Y.Tosaka et al., IEEE Electron Device Letters, 18, 99 (1997).
[12] R.Baumann et al., IEEE International Reliability Physics Symposium (IRPS)
Proceedings, 152 (2000).
[13] R.Baumann et al., Technical Digest of International Electron on Devices Meeting(IEDM), 329 (2002).
[14] E.Ibe et al., International Conference on Information Technology and Applications(ICITA), No.273, 21 (2002).
[15] 児玉章裕, 2005 年度九州大学先端エネルギー理工学専攻修士論文 (2006).
[16] 西嶋康太, 2007 年度九州大学先端エネルギー理工学専攻修士論文 (2008).
[17] H.Iwase, K.Niita, and T.Nakamura, J. Nucl. Sci. and Technol., 39, No.11, 1142 (2002);
http://phits.jaea.go.jp/index.html[18] W.R.Leo, "Techniques for Nuclear and Particle Physics Experiments", Springer-Verlag(2005), pp.49.
[19] R.A.Weller et al., IEEE Trans. on Nucl. Sci., 50, No.6, 2265 (2003).
[20] M.A.Xapsos, Radiation Res. 132, 282 (1992).
[21] H.W.Bertini, ORNL-3383, Oak Ridge National Laboratory (1963).
[22] K.Niita et al., Phys. rev. C, 52, No.5, 2620 (1995).
[23] K.Niita et al., JQMD code, JAERI-Data/Code 99-042 (1998).
[24] 池内丈人, 2001 年度九州大学先端エネルギー理工学専攻修士論文 (2002).
[25] S.Furihata, Nucl. Instr. and Meth. in Phys. Res. B 171, 251 (2000).
[26] T.Satoh, RIST ニュース, No.40, 3 (2005).
[27] T.Nakamura, et al., 日本原子力学会誌, 48, No.12, 949 (2006).
[28] H.W.Armstrong, K.C.Chandler, ORNL-4869, Oak ridge National Laboratory (1973).
57

[29] ATIMA code; http://www-linux.gsi.de/~weick/atima/[30] W.H. Barkas, "Nucl. Res. Emulsions", Vol.1, Academic Press, New York (1963) p.372.
[31] S.P.Ahlen and K.Kinoshita, Phys. Rev. D 26, 2347 (1982).
[32] J. Lindhard et al., Mat. Fys. Medd. Dan. Vid. Selsk, 33, No.14 (1963).
[33] J. Lindhard,

A.H. Soerensen, Phys. Rev. A 53 (1996) 2443.

[34] J.F. Ziegler, J.B. Biersack, U. Littmark, "The Stopping and Range of Ions in Solids",
Vol.1, Pergamon Press (1985).
[35] D.Lambert et al., IEEE Trans. on Nucl. Sci., 51, No.6, 3435 (2004).
[36] ITRS 2005, "配線" (2005); http://strj-jeita.elisasp.net/strj/[37] K.Shibata et al., J. of Nucl.Sci. and Technol., 39, 1125 (2002).
[38] H.Kobayashi et al., IEEE International Reliability Physics Symposium (IRPS) (2009).
[39] A.Akkerman and J.Barak, Nucl. Instr. and Meth. in Phys. Res. B 227, 319 (2005).
[40] F.L.Hassler and R.A.Peck Jr, Phys. Rev., 125, No.3, 1011 (1962).
[41] Y.Terada et al., J. Nucl. Sci. and Technol., Sup. 2, 437 (2002).
[42] F.B.Bateman et al., Phys. Rev. C 60, No.6, 064609 (1999).
[43] S.Benck et al., Nucl. Sci. Eng., 141, No.1, 55 (2002).
[44] UTippawan et al., Phys. Rev. C 73, No.3, 034611 (2006).
[45] D.W.Kneff et al., Nucl. Sci. Eng., 92, No.4, 491 (1986).
[46] J.Ziegler, SRIM-2008; http://www.srim.org/[47] SH.Z.Izmailov et al., Nucl. Instr. and Meth., 168, 81 (1980).
[48] E.I.Sirotinin et al., Nucl. Instr. and Meth. in Phys. Res. B 4, 337 (1984).
[49] P.Mertens and P.Bauer, Nucl. Instr. and Meth. in Phys. Res. B 33, 133 (1988).
[50] A.Ikeda et al., Nucl. Instr. and Meth. in Phys. Res. B 115, 34 (1996).
[51] W.K.Lin et al., J. Appl. Phys., 44, 3631 (1973).
[52] D.C.Santry and R.D.Werner, Nucl. Instr. and Meth., 178, 523 (1980).
[53] J.D.Pearce and R.R.Hart, J. Appl. Phys., 52, 5056 (1981).
[54] D.C.Santry and R.D.Werner, Nucl. Instr. and Meth. in Phys. Res. B 1, 13 (1984).
[55] J.G.Kelley et al., Phys. Rev. B, 8, No.1, 103 (1973).
[56] D.C.Santry and R.D.Werner, Nucl. Instr. and Meth. in Phys. Res. B 69, 167 (1992).
[57] W.Jiang et al., Phys. Rev. B 59, No.1, 226 (1999).
[58] J. of the ICRU, 5, No.1 (2005).
[59] J.Keinonen et al., Appl. Phys. Rett., 60, 228 (1992).
[60] H.J.Whitlow et al., Nucl. Instr. and Meth. in Phys. Res. B 190, 84 (2002).
[61] 柴野望己, 2009 年度筑波大学電子・物理工学専攻修士論文 (2010).
[62] 3D TCAD simulator ENEXSS (http://www.tcad-international.com/product/index.html)
and SELETE (http://www.selete.co.jp/)58

謝辞
本研究を行うにあたって、多くの方々にお世話になりました。
中島秀紀教授には研究の場を提供していただき、深く感謝いたします。渡辺幸信准
教授には学部 4 年から 3 年間、研究を進めるにあたって熱心に御指導頂き、また多くの貴重な機会を与えて頂きましたこと、心から感謝致します。本研究分野の知識が皆無だった私がここまでやってこられたのも、先生の御指導の賜物だとおもっておりま
す。まだまだ至らない点が多々ありますが、今後とも御指導御鞭撻のほど、よろしくお願い申し上げます。山本直嗣助教には日頃から研究および日常生活の面で多くの助言を頂きまして、大変感謝しております。いつも明るく研究室の雰囲気を和やかにして頂きました大神めぐみ秘書に感謝致します。
半導体理工学研究センターの今村健氏、筒井将史氏、戸坂義春氏、古田博伺氏には研究を進める上で専門的な知見から様々な助言を賜りました。この場を借りて御礼申し上げます。また筑波大学の佐野伸行教授、柴野望己氏にはデバイス物理について門外漢である私に多くのアドバイスを頂きました。ここに感謝致します。
日本原子力研究開発機構の坂本幸夫氏および高度情報科学技術研究機構の仁井田浩
二氏には、夏期実習で輸送計算コード PHITS の高度な使用法や詳細な内容などについて御教示頂きました。ここに謝意を表します。
同じ研究グループの皆様には大変お世話になりました。渡邊健人氏には研究に取り組む上での姿勢、心構えを学ばせていただきました。内藤有紀氏には本論文を作成する上で様々な助言を頂き、また日々の研究生活においても色々と楽しませて頂きました。平山嵩祐氏には実験グループの視点から様々な指摘を頂きました。また持ち前の笑顔で研究室の雰囲気を明るくしていただきました。厳しい状況ですが、就職活動の方も頑張って下さい。卒業された先輩方、林真照氏、西嶋康太氏、本岡親英氏、叶涛氏には大変お世話になりました。また KadrevDimtre Nikolov 氏、蘭長林氏には短い間でしたが英語に触れる機会をいただきました。深く感謝しております。
研究分野は異なっておりましたが、大宅将史氏、小谷優介氏、新谷將氏、前野旭弘
氏、江崎徹氏、栗田知明氏、西山貴章氏、花屋倫生氏、下川床潤氏、杉田健策氏、安永将人氏、矢野貴洋氏には研究生活で色々とお世話になりました。大変感謝しております。
最後に、これまでの人生を支えて頂いた両親、姉に心から深く感謝いたします。

59