ALDセミナー

反応性力場分子動力学法および密度汎関数法による化学気相成長および原子層成長プロセスの解析
ReaxFF Molecular Dynamics and Density Functional Theory Studies
for Analysis of Chemical Vapor Deposition and Atomic Layer Deposition Processes

上根直也1,2、馬渕拓哉2,3 、Jin Yong 4 、財津優4 、安原重雄4 、○徳増崇1
1東北大学
2東北大学
3東北大学

流体科学研究所
大学院工学研究科

ファインメカニクス専攻

学際科学フロンティア研究所

4株式会社ジャパン・アドバンスト・ケミカルズ

2022/12/22 15:30-16:10

1/43

Contents

反応性力場分子動力学法および密度汎関数法による化学気相成長および原子層成長プロセスの解析

1. 研究背景2. 方法と理論
3. CVD法における気体分子表面反応機構の反応性力場分子動力学解析
4. PECVD法におけるSiGe成膜機構の反応性力場分子動力学解析
5. ALD法におけるBN成膜機構の実験および反応性力場分子動力学解析
6. 結論

博士論文予備審査会 [2022.12.27]

第1章 | 半導体材料と成膜工程

携帯端末

2/43

自動運転

IOTデバイス

ビッグデータ

人工知能

半導体集積回路（IC）の製造工程

01 基板

09 平坦化処理

02 酸化/窒化膜形成

10 ゲート形成

03 レジスト塗布

11 パターン形成

04 露光

12 イオン注入

05 現像⋯

06 エッチング

17 配線形成

07 剥離/洗浄08 絶縁膜埋込

6/17 成膜工程

半導体材料の成膜は今後の材料開発の中でも重要なトピック[1] NTT DOCOMO, inc., KY-o1L, 2018-11-22[2] Toyota Motor corp., e-Palette, 2018-01-09[3] Amazon.com, inc.，Amazon echo，2020-09-25

博士論文予備審査会 [2022.12.27]

第1章 | CVD（化学気相成長）法とはプラズマCVD成膜装置
原料ガス

成膜過程（気相反応・表面反応）
SiH4 + H2

解離

投入口

プラズマ

排気口
基板

3/43

電子付着

電子衝突
電離

ラジカル

ポリマー

負イオン

SiH3 H
SiH2 Si2Hm
SiH SinHm
Si

SiHm-

H+

SiHm+

Si3H8

Si2Hm-

H2+

Si2Hm+

SinH2n+2

SinHm-

H3+

SinHm+

Si2H6

表面反応

正イオン

[21]

一次気相反応（電子衝突）

二次気相反応（化学種衝突）表面反応

反発（表面衝突）
反射

脱離
化学吸着
DB

加熱機構

薄膜成長

CVD法では気相及び表面における反応動力学現象を経て薄膜が成長[21] J. Plasma Fusion Res. Vol.91, No.53,14‐316 (2015)

博士論文予備審査会 [2022.12.27]

第1章 | ALD（原子層堆積）法とは

4/43

ALD法はCVD法の一種であるが，原子層単位での膜厚制御を可能とする．これは自己終端性によるものであり，以下の４つのステップを繰返すことで実現される．

Chemisorption

Substrate

Substrate

Feed Purge
Purge Reaction

Reactive gas

Growth per cycle (GPC)

ALD法による成膜プロセス(3) Self-limited growth

ALD Window

Uncontrollable
Substrate

Substrate

Controllable

Uncontrollable

Substrate temperature

ALD法は自己終端性を満たす必要があるため反応動力学現象は更に複雑[22] Chem. Rev. 2010, 110 (1), 111–131

博士論文予備審査会 [2022.12.27]

第1章 | CVD/ALDにおける実験的先行研究と数値シミュレーション手法
実験による先行研究（1855年~）
入力データ
原料組成
流入速度

出力データ
相関関係

印加電圧
堆積温度
•

•

•

10 0

結晶構造

10 -3

結合形態
不純物濃度

etc.

Length scale [m]

成膜速度

etc.

Siを中心とした様々な材料のプロセス
パラメータと成膜速度や膜質との相関関係が調査されたプロセスパラメータは膨大な数が存在するため，多変数パラメータの最適化は多大な労力を伴う
本質的な要因を解明することは難しい

数値シミュレーション研究

I
II
III
IV
V

10-9
10-12

V-FEM, FVM:Electrons:Atoms:Microstructure:Film:Reactor

10 -6

5/43

III-KMCII-Classical MD

IV-Geometrical,
Potts, Level-set

薄膜成長

結晶成長

表面反応

I-DFT, QM
10-15

10-12

10-9

10-6

Time scale [sec]

10-3

100

時空間スケールの異なる数値シミュレーションがCVD/ALDプロセスが活用[46] Surf. Coat. Technol. 2004, 177–178, 245–251

博士論文予備審査会 [2022.12.27]

第1章 | 研究目的
研究目的
第2章
第3章
第4章
第5章

6/43

CVD/ALD法による表面近傍での反応動力学現象を正確に理解し，プロセスパラメータが反応機構や薄膜成長機構に及ぼす影響を解明する．
手法と理論
計算化学の手法を概観した後に，本研究で使用する計算手法について理論を述べる．

CVD法における気体分子表面反応機構の反応性力場分子動力学解析

Si薄膜を成長させるPECVDプロセスにおいて生成される気体分子を単体で表面に衝突させることによって，反応動力学現象を
表面素過程に分解して解明する．

PECVD法におけるSiGe成膜機構の反応性力場分子動力学解析

SiGe薄膜を成長させるPECVDプロセスにおいて生成されるラジカルを連続的に表面に衝突させて成長薄膜の組成および結晶性を解析し，反応動力学的過程を含む反応機構および成膜機構を解明する．

ALD法におけるBN成膜機構の実験，密度汎関数および反応性力場分子動力学解析
BN薄膜を成長させるALDプロセスにおいて，ALDウィンドウの存在を実験的に確認する．

また，密度汎関数法および反応動力学的手法によって，基板温度が成膜機構に及ぼす影響を解明する．

第6章

結論
本研究の結言を述べる

博士論文予備審査会 [2022.12.27]

7/43

Contents

反応性力場分子動力学法および密度汎関数法による化学気相成長および原子層成長プロセスの解析

1. 研究背景2. 方法と理論
3. CVD法における気体分子表面反応機構の反応性力場分子動力学解析
4. PECVD法におけるSiGe成膜機構の反応性力場分子動力学解析
5. ALD法におけるBN成膜機構の実験および反応性力場分子動力学解析
6. 結論

博士論文予備審査会 [2022.12.27]

第2章 | 計算科学の全体像

8/43
10 0

マクロスケール
モデル

ミクロスケール
モデル

10 -3
10 -6
10-9

ナノスケールモデル

10-12化学反応と動力学の両方を考慮可能な
数値シミュレーション手法が必要

Length scale [m]
I
II
III
IV
V

:Electrons:Atoms:Microstructure:Film:Reactor

V-FEM, FVM
IV-Geometrical,
Potts, Level-set

III-KMC

薄膜成長

結晶成長

II-Classical MD

ReaxFF MD
表面反応

I-DFT, QM
10-15

⇒ 動力学

⇒ 化学反応＋動力学⇒ 化学反応

10-12

10-9

10-6

Time scale [sec]

10-3

100

比較的高精度なDFT法と反応動力学現象を考慮可能なReaxFF MD法を選択
博士論文予備審査会 [2022.12.27]

第2章 | 密度汎関数法（Density Functional Theory: DFT）

9/43

電子密度が求まれば系の基底状態における電子エネルギーが確定するというHohenberg-Kohnの定理に基づいたKohn-Sham方程式を解くことで，電子密度と全電子エネルギーを計算する方法．

• Schrödinger Equation

Hohenberg–Kohn theorems• Kohn-Sham Equation

ZA
1 2[− ∇ − ∑+ V eff ]ϕk =ε kϕ k
2
A rA

AccuracyHψ = Eψ

for ground state method

DFT
HF
|

MP2

|

CCSD CCSD(T)
SAC
EOM
MP4
CISD

|

|

Full CI|

N7 … N!
N6
N3~N4 N5
Cost (N= number of electrons)

計算コストと精度の観点からDFT法を採用
博士論文予備審査会 [2022.12.27]

第2章 | Molecular Dynamics (MD) simulations

10/43

原子・分子を質量を持つ粒子として扱い、古典的な運動方程式を
逐次数値的に解くことで，多体系の時間発展を模擬する手法
N
d 2 rimi 2  
ij  Fi
Fi  （ニュートンの運動方程式）dtj 1,i  jri t  t   ri t   vi t  t Fi t 2mi

t 2  O t 3 

1 dFi t  2vi t  t   vi t  t t +O t 3 2m
2 dt
Fi t Fi t  t   Fi t  ri t  t   ri t   vi t  t vi t  t   vi t  Fi t dt

t +O t 2 Fi t t 2  O t 3 2m
Fi t  t   Fi t tt  O t 3 

ϕki

ϕij
j

ϕjk

k

ri : i原子の位置
mi : i原子の質量
φij : i原子に作用するエネルギーNd 2 rimi 2    ij  Fidtj 1,i  j

自身が対象とする系に合わせて適切な力場（Force-Field：FF）を選定する必要
博士論文予備審査会 [2022.12.27]

第2章 | ReaxFF MD法の理論

11/43

Esystem = Ebond + Eangle + Etors + Econj + Eover + Eunder+ Ecoa + E pen + Elp + EC 2 + Etriple + EH −bond + EvdWaals + ECoulomb
結合次数依存項

結合次数（Bond order: BO）

結合性/反結合性軌道の電子数の差•

化学結合の強さを定量的に表す値•

ReaxFFではσ，π，ππ結合に分割した連続値として扱う
エネルギ

↑

低

結合性→

Na

1s

↑

↑

Nb

↑

↑

↑

1s

↑

↑

反結合性→



BOij  BOij  BOij  BOij

•

高



結合次数非依存項

BO=(Nb-Na)/2



pbo 2 
rij
 exp  pbo1      σ結合
 ro  
pbo 4 
rij
 exp  pbo 3      π結合
 ro  
pbo 6 
rij
 exp  pbo 5      ππ結合
 ro  

ReaxFF MDは結合次数依存の関数形式を採用し化学反応を再現可能
博士論文予備審査会 [2022.12.27]

12/43

Contents

反応性力場分子動力学法および密度汎関数法による化学気相成長および原子層成長プロセスの解析

1. 研究背景2. 方法と理論
3. CVD法における気体分子表面反応機構の反応性力場分子動力学解析
4. PECVD法におけるSiGe成膜機構の反応性力場分子動力学解析
5. ALD法におけるBN成膜機構の実験および反応性力場分子動力学解析
6. 結論

博士論文予備審査会 [2022.12.27]

第3章 | PECVDプロセスにおける気体分子の表面反応
原料ガス

成膜過程（気相反応・表面反応）
SiH4 + H2

解離
ラジカル
SiH3 H
SiH2 Si2Hm
SiH SinHm
Si

投入口

プラズマ

排気口

ポリマー

電子付着

電子衝突
電離

負イオン

正イオン

SiHm-

H+

SiHm+

Si3H8

Si2Hm-

H2+

Si2Hm+

SinH2n+2

SinHm-

H3+

SinHm+

Si2H6

表面反応

反発

1022
1021
1020
1019

Density [ m-3 ]プラズマCVD成膜装置

13/43

1017
10

1015

脱離

Si2H6

13

13.56 MHz~ 0.1 W/cm2~ 0.1 Torr

Generatedstable species
SiH 3

Neutralradicals

H
SiH 2
SiH
Si

Si2Hy
Si3Hz

High-order
radicals
SiH x+

Charged speciese
SiH x-

1012

化学吸着
DB

加熱機構

16

1014

反射

H2

1018

10

基板

SiH 4 Parent molecule

・親分子であるSiH4・主要成膜前駆体であるSiH3・膜質劣化の原因として知られるSiH2
これらの主要分子の表面反応機構解明はPECVDプロセスのモデル化に向けて重要

H原子が被覆したSi表面上で基板温度がSiHx (x = 2-4)の表面反応に及ぼす影響を解明
博士論文予備審査会 [2022.12.27]

第3章 | 反応性力場分子動力学法（ReaxFF MD）

14/43

Esystem = Ebond + Eangle + Etors + Econj + Eover + Eunder+ Ecoa + E pen + Elp + EC 2 + Etriple + EH −bond + EvdWaals + ECoulomb
結合次数依存項( (− Deσ + BOijσ ⋅ exp  pbe1 1 − BOijσ
Ebond =

− Deπ ⋅ BOijπ − Deππ ⋅ BOijππ
•
•

)

pbe 2

)

120




Energy [kcal/mol]

結合エネルギー項

結合解離エネルギー：De
結合エネルギーパラメータ：pbe

角度エネルギー項E=
f 7 ( BOij ) ⋅ f 7 ( BO jk ) ⋅ f8 ( ∆ j )
angle

{

}

2
⋅ pval1 − pval1 exp  − pval 2 (θ o ( BO ) − θijk ) 



•

角度エネルギーパラメータ：pval100

結合次数非依存項
20

ReaxFF (This study)
ReaxFF (Ref.)
DFT (GGA/PBE)

15

80

10

60
40

5

20
0
1

ReaxFF (This study)
ReaxFF (Ref.)
DFT (GGA/PBE)

2
3
4
Si–H bond length in SiH 4 [Å]

5

0

80

100
120
H–Si–H angle in SiH 4 [°]140

DFT法から得られる結果を再現するように既存のポテンシャルパラメータを修正
博士論文予備審査会 [2022.12.27]

第3章 | 計算系および計算条件

15/43
 系全体

×2000

×2000

×2000

10Å

サイズ：30.72×26.88×25.00 [Å3]
時間：15 ps（0.25 fs × 60000）
試行回数：各条件2000回

Θ = 0.00

 気体分子

化学種：SiH4 /SiH3 /SiH2
入射速度：300K （最確速度）
位置：X, Y = 無作為，Z = 固定

計算系の正面図
H（気体分子）Si（気体分子）Si（基板原子）H（被覆原子） 水素被覆基板
原子種/結晶面：Si(100)-(2×1)
原子数：224個（Si）
H被覆率：0.00-1.00（0.25刻み）
基板温度：300K-1200K（100 K刻み）

Θ = 1.00

基板温度とH被覆率を独立して変化させた時に熱運動が表面反応に及ぼす影響を解析
博士論文予備審査会 [2022.12.27]

第3章 | 衝突後の気体分子の挙動

16/43
表面拡散なし1. 反射（Reflection）

ー 衝突後に表面拡散（Z方向の速度が2回以上変化）せずに直接的に計算系外へ出た

2.
3.
4.

Z 方向

気体分子が表面に衝突した後の挙動に応じて，4種類のイベントに分類した．

表面拡散あり反射 ＠基板温度Ts = 300K, 被覆率Θ = 1.00

脱離（Desorption）

ー 衝突後に表面拡散を経て，計算系外へ出た

化学吸着（Chemisorption）

ー 計算終了時にSisub-Sigas結合が形成された物理吸着（Physisorption）

ー 計算終了時に基板上に存在し Sisub-Sigas結合が形成されない（Si-H結合および物理吸着を含む）

特定の瞬間の幾何学的構成と運動状態に基づいて表面反応を分類
博士論文予備審査会 [2022.12.27]

第3章 | 衝突後の気体分子の挙動

17/43
表面拡散なし1.

反射（Reflection）

ー 衝突後に表面拡散（Z方向の速度が2回以上変化）せずに直接的に計算系外へ出た

Z 方向

気体分子が表面に衝突した後の挙動に応じて，4種類のイベントに分類した．

表面拡散あり脱離 ＠基板温度Ts = 1200K, 被覆率Θ = 0.50

2. 脱離（Desorption）

ー 衝突後に表面拡散を経て，計算系外へ出た

3.
4.

化学吸着（Chemisorption）

ー 計算終了時にSisub-Sigas結合が形成された物理吸着（Physisorption）

ー 計算終了時に基板上に存在し Sisub-Sigas結合が形成されない（Si-H結合および物理吸着を含む）

特定の瞬間の幾何学的構成と運動状態に基づいて表面反応を分類
博士論文予備審査会 [2022.12.27]

第3章 | 衝突後の気体分子の挙動18/43
表面拡散なし1.

2.

反射（Reflection）

ー 衝突後に表面拡散（Z方向の速度が2回以上変化）せずに直接的に計算系外へ出た

Z 方向

気体分子が表面に衝突した後の挙動に応じて，4種類のイベントに分類した．

表面拡散あり化学吸着 ＠基板温度Ts = 1200K, 被覆率Θ = 0.00

脱離（Desorption）

ー 衝突後に表面拡散を経て，計算系外へ出た

3. 化学吸着（Chemisorption）

ー 計算終了時にSisub-Sigas結合が形成された

4.

物理吸着（Physisorption）

ー 計算終了時に基板上に存在し Sisub-Sigas結合が形成されない（Si-H結合および物理吸着を含む）

特定の瞬間の幾何学的構成と運動状態に基づいて表面反応を分類
博士論文予備審査会 [2022.12.27]

第3章 | Si基板温度依存性（SiH4，θ=0.00）
100

Probability [%]

低温

80

高温

19/43
@300 K

Reflection
Desorption
Chemisorption
Physisorption@1200 K

60
40
20
0
300

900
600
Substrate temperature [K]

1200

低基板温度

高基板温度

Sigas－Hgas結合

気体分子

解離優位

脱離優位

化学吸着確率

化学吸着確率40%程度増加30%程度減少

SiH4は気体分子内の結合解離効果と脱離効果のバランスで化学吸着確率が決定する．
⇒ 従来モデルでは熱運動による脱離効果は未考慮のため精度の大幅な改善が期待
博士論文予備審査会 [2022.12.27]

第3章 | H被覆Si基板温度依存性（SiH2，θ=1.00）

Probability [%]

100
(c)
80

20/43
@300 K

Reflection
Desorption
Chemisorption
Physisorption@1200 K

60
40

低基板温度

20
0
300

600
900
Substrate temperature [K]

1200

基板

価電子なし

気体

価電子２個

高基板温度
基板

価電子なし

気体

価電子２個

H原子を取り込んでH原子を取り込んでSiH3として化学吸着

SiH3として脱離被覆H原子に衝突したSiH2は低基板温度では化学吸着し，高基板温度では脱離
博士論文予備審査会 [2022.12.27]

21/43

Contents

反応性力場分子動力学法および密度汎関数法による化学気相成長および原子層成長プロセスの解析

1. 研究背景2. 方法と理論
3. CVD法における気体分子表面反応機構の反応性力場分子動力学解析
4. PECVD法におけるSiGe成膜機構の反応性力場分子動力学解析
5. ALD法におけるBN成膜機構の実験および反応性力場分子動力学解析
6. 結論

博士論文予備審査会 [2022.12.27]

第4章 | SiGe薄膜の産業的応用例と課題

22/43

SiGe薄膜はプロセスパラメータを変化させることで，バンドギャップなどの物性を段階的に制御可能であることから，産業的に数多くの応用例が見られる．

PECVD装置
Sourcegas

Region 2(SiH2+GeH3)

Inlet

SiH2

PlasmaOutletSubstrate

Heater

PECVDプロセスで生成される主要な気体分子種密度の圧力依存性 [38]Region 3(SiH2+GeH2)

GeH3
2 1
3 4

Region 1(SiH3+GeH3)
SiH3
Region 4(SiH3+GeH2)

GeH2

気体分子組成が薄膜組成および構造に及ぼす影響と成膜機構の解明を目的[1] J. Alloys Compd. Vol. 803, No. 30, 260–264 (2019)[2] J. Electron. Mater., Vol. 48, No. 11, 7360–7365 (2019)[3] Phys. Status Solidi, Vol. 215, No.12, 1700736 (2018)[38] MRS Online Proceedings Library volume 1426, 403–408 (2012)

博士論文予備審査会 [2022.12.27]

第4章 | 計算系詳細

23/43
 系全体

系大きさ：30.72×26.88×45.00 [Å3]緩和時間：7.5 ps / 60 ps（前/後）
時間刻み：0.25 fs

 気体分子分子組成：SiH3/SiH2 + GeH3/GeH2を5%刻み
位置：X, Y = 無作為，Z = 1-2 Å/ns
入射間隔：5.5 ps
入射総数：成長原子 (Si+Ge)が400個（±5%）: 追加H原子: 基板Si原子: 追加Si 原子: 追加Ge原子

 基板原子原子種/結晶面：Si(100)-(2×1)
基板温度：1500K

気体分子を表面に対して連続的に照射することで薄膜を成長
博士論文予備審査会 [2022.12.27]

第4章 | 解析方法

24/43

薄膜内の組成および構造は公開可視化ツール（OVITO）を用いて解析を行った．

 組成（原子含有量：atomic content）
Si content =# of Si atoms in film
×100
# of total atoms in film

 構造（結晶性：Crystallinity）
Cd =





Cubic + 1st neighbor + 2nd neighbor
N depoCubic diamond
1st neighbor2nd neighbor

Cubic：第1近接原子および第2近接原子の全てが立方体ダイヤモンド格子位置に配置されている
1st neighbor： Cubicと判定された原子の最初の隣接原子．4つの隣接原子は格子位置に配置されているが，2番目に近い隣接原子の少なくとも1つは配置されていない．
2nd neighbor ：Cubicと判定された原子の二番目の隣接原子．原子自体は格子位置に配置されているが，隣接の少なくとも1つが欠落しているか格子に配置されていない．

気体分子を連続的に入射させて得られた成長薄膜の組成および構造を定量的に算出[20] Comp. Phys. Commu., Vol. 205, 13–21 (2016)

博士論文予備審査会 [2022.12.27]

第4章 | 成長薄膜のスナップショット
Region 1

25/43

No. 01

No. 02

No. 03

No. 04

No. 05

No. 06

No. 07

No. 08No. 09No. 10No. 11

No. 12

No. 13

No. 14

No. 15

No. 16

No. 17

No. 18

No. 19No. 20

0%

5%

10%

15%

20%

（SiH3+GeH3）

Region 2（SiH2+GeH3）

Region 3（SiH2+GeH2）

Region 4（SiH3+GeH2）

Low

High

博士論文予備審査会 [2022.12.27]

第4章 | 組成の気体分子組成依存性

26/43
膜内Ge含有量の気体分子組成依存性• 気相Ge化学種の増加に対して，膜内Ge含有量は線形に増加する傾向が確認された• 気相中の化学種分布が特定されれば，膜内組成が推定可能であることを示唆している• SiH3 から構成される薄膜に対するGeH2 からのGe原子取り込み量はGeH3と比べて3倍程度大きい• SiH2 およびGeH2 は100％化学吸着するとしてKMCモデ
ルでは再現できない，Ge原子の優先的な取り込みを確認

気体分子組成が成長薄膜の組成に及ぼす影響を解明[21] Comp. Phys. Commu., Vol. 205, 13–21 (2016)

博士論文予備審査会 [2022.12.27]

第4章 | 組成の気体分子組成依存性

27/43

SiH3化学吸着@300 K

192.50 ps
傾きの変化が確認された時の表面状態

SiH3脱離吸着@1200 K

SiH2による被覆H原子の取り込み反応が熱運動の脱離効果によって阻害されたことで，SiH2は必ず化学吸着することを仮定した従来のモデルとの差異が確認[21] Comp. Phys. Commu., Vol. 205, 13–21 (2016)

博士論文予備審査会 [2022.12.27]

第4章 | 結晶性の気体分子組成依存性

28/43

結晶性の気体分子組成依存性• SiH3から成長した薄膜の平均結晶性は73.1%，SiH2は58.8％• Ge気体分子種比率の増加に伴い結晶性は減少する傾向が確認• SiH3から成長した薄膜はGeH3の供給では高い結晶性を維持するが，GeH2の供給は結晶性を急激に低下

従来のKMC法では困難な気体分子種が薄膜の結晶性に及ぼす影響を解明
博士論文予備審査会 [2022.12.27]

第4章 | XH3およびXH2ラジカルの反応機構

29/43

XH3ラジカルの反応機構••

1つのSiH3ラジカルが表面で拡散し，表面のH原子を1つ引き抜き，
1つのSiH4分子と1つのダングリングボンドを生成する．
ダングリングボンドに別のSiH3ラジカルが化学吸着する．

SiH 3  g   H  s   SiH 4  g   DB  s 

SiH 3  g   DB  s   SiH 3  s 

2段階反応

Eley-Rideal型機構

XH2ラジカルの反応機構
SiH 2  g   H  s   SiH 3  g   DB  s or SiH 2  g   2H  s   SiH 4  g   2DB  s 

SiH 2  g   DB  s   SiH 2  s 
•

2段階反応

Eley-Rideal型機構

SiH2ラジカルは表面のH原子を1つ引き抜いた後，SiH3として衝突した
H原子の下にあるSi原子と新たにSi-Si結合を形成する．

SiH 2  g   H  s   SiH 3  s 

1段階反応

H捕獲型機構

XH3ラジカルはER型機構，XH2ラジカルはER型機構とH捕獲機構の両方を確認
博士論文予備審査会 [2022.12.27]

第4章 | XH3およびXH2ラジカルの成長機構

30/43
衝突シミュレーションの計算結果
Surface events

SiH3吸着表面上のSiH3

SiH2吸着表面上のSiH3

XH3ラジカル
Eley-Rideal型機構
水素含有量が少なく，高い結晶性の薄膜を成長

各条件で1000回の衝突

SiH2 desorption
SiH3 desorption
SiH4 desorption
Chemisorption

Gaseous species /Adsorbed species on Si(100) surface
SiH3 / SiH2
SiH2 / SiH3
SiH2 / SiH2
SiH3 / SiH3
0.0
3.5
20.0
0.0
13.2
6.5
27.7
3.5
83.5
36.7
43.3
0.0
3.3
53.3
9.0
96.5

XH2ラジカル

Eley-Rideal型機構＋
H捕獲型機構
水素含有量が多く，低い結晶性の薄膜を成長

XH3ラジカルとXH2ラジカルの異なる成長機構が膜質におよぼす影響を解明
博士論文予備審査会 [2022.12.27]

31/43

Contents

反応性力場分子動力学法および密度汎関数法による化学気相成長および原子層成長プロセスの解析

1. 研究背景2. 方法と理論
3. CVD法における気体分子表面反応機構の反応性力場分子動力学解析
4. PECVD法におけるSiGe成膜機構の反応性力場分子動力学解析
5. ALD法におけるBN成膜機構の実験および反応性力場分子動力学解析
6. 結論

博士論文予備審査会 [2022.12.27]

第5章 | BNの広域な産業応用範囲

32/43

ホウ素と窒素から成るホウ素窒化物（Boron Nitride: BN）は0次元から3次元まで様々な構造を持つ化合物であり，幅広い産業的な応用先が存在している魅力的な材料である．

BNナノ材料の産業的応用先

0D ナノ粒子=>がんの予防および治療薬

1D ナノチューブ=>機能性複合材料の補強材

2D ナノシート

Fig. Structural models of representative BN nanomaterials:=>デバイス用材料誘電体ゲート層，コーティング (a) 0D B N fullerene, (b) 1D BN single-walled nanotube,
24 24(c) 2D mono-layered BN nanosheet. [1]
3D ナノポーラス=>水素貯蔵，汚染物質処理，触媒担体

BNは様々な産業的応用が可能であり魅力的な材料として期待されている[1] Journal of Materials Science & Technology 31 (2015) 589-598

博士論文予備審査会 [2022.12.27]

第5章 | 開発されてきたBN用のALDプロセス

33/43

ALD法はコンフォーマルで極薄の誘電体を成長させることができるため，BN薄膜を成長させる最も効果的な方法の一つとして期待されている．

Reactant

Boron precursors

BCl3

BBr3

TEB*

B2H6

TDMAB*NH3

ALDウィンドウの温度領域および基板温度が成膜機構に与える影響は未解明[2] Matter 2, 587–630, March 4, 2020
*TEB: Triethylboron, TDMAB: Tris(dimethylamino)boron

博士論文予備審査会 [2022.12.27]

第5章 | BCl3+NH3使用したBN-ALDプロセスのDFTシミュレーション
Cl
H
Cl

B

2x

Cl

Cl
H

H

B

O

O

O

Si

Si

2x

B

Si

Si

Si

Start

Reaction (a)Reaction (b)

H

2HCl

2HCl

H

H

Cl

N

N

H

Si

B

N

4x

Cl

H

HCl
Cl

Si

N

H

Cl

34/43

H

Reaction (c)

H

4HCl

Cl

Cl

B Cl B

H

N

O

Si

H

Si

B

N

Cl

H

O

Si
Reaction (d)4つの反応についてDFTシミュレーションを実施し，エネルギーダイアグラムを取得
博士論文予備審査会 [2022.12.27]

第5章 | 反応(a)および(b)に対するDFTシミュレーション結果
Reaction (a) : BCl3 + Si(OH)4 → BCl2OSi(OH)3 + HCl

35/43
SiliconChorine

OxygenBoron

Hydrogen

Nitrogen

Eact =15.4 kcal/mol
Eact =17.0 kcal/mol(1) Initial state: IS (2) Reactant: R (3) Transition (4) Product: P
state: TS
Reaction (b) : NH3 + Cl2BNH2 → BCl(NH2) + HCl(5) Final state: FS(1) Initial state: IS (2) Reactant: R (3) Transition (4) Product: P
state: TS(5) Final state: FS反応(a), (b)の反応障壁は15.4 kcal/mol，17.0 kcal/molと予測
博士論文予備審査会 [2022.12.27]

第5章 | 反応(c)および(d)に対するDFTシミュレーション結果
Reaction (c) : BCl3 + H2BNH2 → H2BNHBCl2 + HCl

36/43
SiliconChorine

OxygenBoron

Hydrogen

Nitrogen

Eact =30.0 kcal/mol
Eact =16.9 kcal/mol(1) Initial state: IS (2) Reactant: R (3) Transition (4) Product: P
state: TS
Reaction (d) : NH3 + BCl2OH → NH2BClOH + HCl(5) Final state: FS(1) Initial state: IS (2) Reactant: R (3) Transition (4) Product: P
state: TS(5) Final state: FS反応(c)の反応障壁は30.0 kcal/molと予測され，ALDウィンドウの下限を律速博士論文予備審査会 [2022.12.27]

第5章 | BN-ALDを実現するための新たなReaxFFの構築

37/43

BCl3とNH3を用いたOH終端Si(100)表面へのBN成長をシミュレーションする力場は存在しないため，
ALDプロセスにおける反応動力学現象を適切に模擬する新たな力場を構築する必要がある．
新規力場構築のベースとした既存の力場

300

40

DFT
ReaxFF

250

DFT
ReaxFF

35

200

Energy (kcal/mol)

Energy (kcal/mol)

30

150

100

50

1.5

2.0

2.5

3.0

3.5

4.0

4.5

B-Cl distance (Å)

B/N/H/Si/O/Cl を含む初期力場を構築

15
10

0
1.0

Si/C/O/Cl 力場

20

5

0

B/N/H 力場

25

BCl3のB‐Cl距離に対するエネルギー変化

5.0

5.5

-5
90

100

110

120

130

140

150

Cl-B-Cl angle (°)

BCl3のCl-B‐Cl距離に対するエネルギー変化

BCl3の分子内ポテンシャルが一致するように力場を訓練構築した初期力場をベースにBCl3の分子内ポテンシャルが一致するように力場を訓練[33] The Journal of Physical Chemistry B 2018 122 (14), 4083-4092[34] Nanoscale, 2019, 11, 7447-7456

博士論文予備審査会 [2022.12.27]

第5章 | BN-ALDを実現するための新たなReaxFFの構築•

BCl3の幾何学的構造に加えて反応経路も訓練に追加•

BN-ALDプロセスにおけるOH-Si(100)およびBN表面でのBCl3およびNH3の反応をモデル化

4BCl 3 +Si ( OH )4 → 4HCl + Si ( OBCl2 )4ー 反応 (a)BCl 3  NH 3  HCl  Cl2 BNH 2ー 反応 (b), (c), (d)
15
10

38/43

ReaxFF w/o train.
ReaxFF w train.
DFT

Energy (kcal/mol)

5
0
-5
-10
-15
-20
-25

Reaction coordination

OH-Si(100)およびBN表面でのBCl3およびNH3の反応経路を模擬できるようにモデル化
博士論文予備審査会 [2022.12.27]

第5章 | BCl3とNH3を用いた薄膜成長シミュレーション

39/43

開発した力場がBCl3とNH3を用いたBN薄膜の成長を適切にモデル化できることが検証された．しかし，成膜機構を理解するためには，表面での反応動力学現象を含む成膜シミュレーションが必要である．
(1) BCl3 Exposure

BCl3(2) 1st Purge(3) NH3 Exposure(4) 2nd Purge
ChorineNH3

Boron
Nitrogen

OH-Si(100)

Hydrogen
Oxygen
Silicon×5 ALD Cycles

1500 Kにおいて5ALDサイクルを含む薄膜成長シミュレーションを実行
博士論文予備審査会 [2022.12.27]

第5章 | 薄膜原子組成の時間発展
(a)

40/43
(b)

After 1st ALD cycle

After 2nd ALD cycle

After 4th ALD cycle

After 5th ALD cycle

After 3rd ALD cycle

3次元的なクラスター成長の後に2次元的な薄膜成長に変化する様子が確認
博士論文予備審査会 [2022.12.27]

第5章 | 薄膜原子組成の時間発展
1st ALD cycle 2nd ALD cycle 3rd ALD cycle 4th ALD cycle 5th ALD cycle

41/43• BCl3の供給ステップにおいて，H原子の減少，Cl原子およびB原子の増加が確認⇒ BCl3によるHCl分子の脱離およびBClxの吸着• NH3の供給ステップにおいて，Cl原子の減少，H原子
およびN原子の増加が確認⇒ NH3によるHCl分子の脱離およびNHxの吸着• ALDサイクルを繰り返すことで，原子数が段階的に増減する様子が確認

連続的な気体分子の入射により段階的にBN薄膜が成長する様子が確認
博士論文予備審査会 [2022.12.27]

42/43

Contents

反応性力場分子動力学法および密度汎関数法による化学気相成長および原子層成長プロセスの解析

1. 研究背景2. 方法と理論
3. CVD法における気体分子表面反応機構の反応性力場分子動力学解析
4. PECVD法におけるSiGe成膜機構の反応性力場分子動力学解析
5. ALD法におけるBN成膜機構の実験および反応性力場分子動力学解析
6. 結論

博士論文予備審査会 [2022.12.27]

第6章 | 結論

43/43

CVD/ALD法による表面近傍での反応動力学現象を正確に理解し，プロセスパラメータが反応機構や薄膜成長機
構に及ぼす影響を解明することを目的として，ReaxFF MD法およびDFT法を用いた数値シミュレーションを実施．
第3章：CVD法における気体分子表面
反応機構の反応性力場分子動力学解析

Probability [%]100
80

第4章：PECVD法におけるSiGe成膜機構の反応性力場分子動力学解析

第5章：ALD法におけるBN成膜機構の実験および反応性力場分子動力学解析
低温度領域

Reflection
Desorption
Chemisorption
Physisorption

中温度領域

高温度領域

60
40
20
0
300

600
900
Substrate temperature [K]1200

ReaxFF MD法による反応動力学現象を考慮した数値シミュレーションはCVD/ALD法における反応機構や成膜機構を解明する上で有効 ⇒ 多くのプロセスに適応可能
博士論文予備審査会 [2022.12.27]

Thank you for your kind attentions• 本 研 究 の 一 部 は JSPS 科 研 費 （ JP20J20915 お よ び
JP21H01243）の助成を受けたものである．
• 本研究の数値計算は， 東北大学流体科学研究所未来流体情
報創造センターの次世代融合研究システム(AFI-NITY)を用いて行った．ここに記して， 厚く謝意を示す