修士論文

分子動力学シミュレーションによる

OLC の摩擦メカニズム評価

指導教員： 屋代

西村

如月

英晃

2013 年 2 月
神戸大学大学院 工学研究科

博士課程前期課程

機械工学専攻Master thesis

Molecular Dynamics Studies
on Friction Mechanism of OLC

Hideaki NISHIMURA
Feburary 2013
Department of Mechanical Engineering,
Graduate School of Engineering,
Kobe University, Kobe, Japan要約
OLC は球殻構造を有しているために，固体潤滑剤としての応用が期待される．本研究では，OLC 薄膜の摩擦特性について原子レベルからの知見を得るため，Stone-Wales
欠陥を導入した炭素原子数が 60n2 (n=1〜6) で表わされるフラーレン，並びに，それらを入れ子構造にした OLC を作成し，それぞれ単体の圧縮・摩擦特性を分子動力学
シミュレーションにより検討した後，平面上に並べた薄膜構造のフラーレン・OLC に対して摩擦シミュレーションを行った．頂点保持による圧縮シミュレーションにより，多層構造を有することで OLC の方がフラーレンよりも高い抵抗力を示すこと，OLC
は内部のフラーレンほど強い力を受けることなどを明らかにした．平滑なダイヤモン
ド壁を用いて圧縮を行った場合では，フラーレンより OLC の方が，径の大きなものよ
り小さなものの方が早い段階で高い圧縮応力を示すこと，OLC は中心部の C60 が崩壊したときに応力が急上昇することなどがわかった．単一のフラーレン・OLC に対して平滑なダイヤモンド壁を用いて摩擦シミュレーションを行った場合では，摩擦係数は
摩擦時のつぶれ具合に関係なく 10−2 程度のオーダーを示すこと，半分つぶした ＠C960と ＠C1500 以外は転がることなく滑るように移動することなどを示した．またフラーレ
ンの摩擦係数は径の大きさによる整理はできなかったが，OLC は径の大きなものほど低い摩擦係数を示した．最後の薄膜構造のフラーレン・OLC に対する摩擦シミュレー
ションでは，圧子・基板の表面に凹凸を与えると，実際に測定されるような 10−1 オー
ダーの摩擦係数となったが，このときフラーレン・OLC 共に圧子・基板の凹凸に保持され回転するような変形を受けていた．フラーレンは起伏に入り込みやすいために OLC
より高い摩擦係数を示していた．

Summary
OLC is expected as new solid lubricant because of its spherical structure. For a new
insight on the friction properties of OLC thin ﬁlm, various molecular dynamics (MD)simulations are performed on the compression and friction of isolated fullerene/OLC
and thin ﬁlms composed of the array of fullerenes/OLCs. First, we have made spherical fullerenes based on the polyhedral rule of 60n2 (n=1〜6) atoms with Stone-Wales
defects. OLCs were also made by nesting these fullerenes. Then we performed various compression simulations on isolated fullerene/OLC. In the point compression byholding the ﬁve-membered ring at the top and bottom of fullerene/OLC, the OLCs
showed higher strength than the fullerenes because of its multi-layer structure. Detailobservation revealed that the internal fullerene in the OLC receives higher force than
the isolated fullerene. According to the compression by a ﬂat diamond wall, we also
found the following facts; (1) the OLCs show higher stress than the fullerenes at theearly stage of compression, (2) the smaller fullerene/OLC shows higher stress than thelarger ones, and (3) the OLC shows drastic stress increase when the center C60 was
collapsed. Then we performed various scratch simulations on isolated fullerene/OLC
by a ﬂat diamond wall changing the indentation depth. The friction coeﬃcient showed
the order of 10−2 regardless of the indentation depth or collapse morphologies of thetarget (fullerene/OLC). Here, the target carbons glide without rotation under scratch,
except the OLCs of half crushed ＠C960 and ＠C1500 . We cannot ﬁnd out the relationship between the friction coeﬃcient and the size of fullerenes; however, the large OLC
showed lower friction coeﬃcient than the small OLC. Finally, we performed scratch
simulations on fullerene/OLC thin ﬁlm. When we introduced the surface roughnessof a indenter/substrate, the friction coeﬃcient increased to the order of 10−1 . Thismagnitude agrees with experimental result. All the fullerene/OLC rotate under thescratch, because they were held by the serrated surface. We found that the fullereneshows much higher friction coeﬃcient than the OLC since they can deform to ﬁt theirshape to the surface roughness.

目次
第1章 緒 論1

第 2 章 解析手法の基礎

42.1

分子動力学法 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4

2.2

原子間ポテンシャル. . . . . . . . . . . . . . . . . . . . . . . . . . . .

5

2.2.1

Brenner ポテンシャル . . . . . . . . . . . . . . . . . . . . . . .

5

2.2.2

Lennard-Jones ポテンシャル . . . . . . . . . . . . . . . . . . . .

9

2.2.3

力の表式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

10

2.3

高速化手法

2.4

速度スケーリング法

2.5

フラーレンモデルの幾何形状. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . . . . . . . . . . .
. . . . . . . . . . . . . . . . . . . . . . .

第 3 章 単一フラーレンの解析

3.1

3.2

3.3

15
16
18

初期構造緩和シミュレーション . . . . . . . . . . . . . . . . . . . . . .

18

3.1.1

解析条件 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

18

3.1.2

解析結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

18

圧縮シミュレーション . . . . . . . . . . . . . . . . . . . . . . . . . . .

20

3.2.1

頂点保持による圧縮 . . . . . . . . . . . . . . . . . . . . . . . .

20

3.2.2

ダイヤモンド壁による圧縮 . . . . . . . . . . . . . . . . . . . . .

27

ダイヤモンド壁による摩擦シミュレーション . . . . . . . . . . . . . . .

35

3.3.1

解析条件 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

35

3.3.2

解析結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

36

第 4 章 単一 OLC の解析

4.1

13

43

初期構造緩和シミュレーション . . . . . . . . . . . . . . . . . . . . . .

43

解析条件 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

43

4.1.1

i

目次

ii

解析結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

43

圧縮シミュレーション . . . . . . . . . . . . . . . . . . . . . . . . . . .

45

4.2.1

頂点保持による圧縮 . . . . . . . . . . . . . . . . . . . . . . . .

45

4.2.2

ダイヤモンド壁による圧縮 . . . . . . . . . . . . . . . . . . . . .

52

ダイヤモンド壁による摩擦シミュレーション . . . . . . . . . . . . . . .

56

4.3.1

解析条件 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

56

4.3.2

解析結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

56

4.1.2
4.2

4.3

第 5 章 薄膜構造での摩擦シミュレーション

5.1

5.2

63

圧子の表面起伏が摩擦係数に与える影響 . . . . . . . . . . . . . . . . .

63

5.1.1

解析条件 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

63

5.1.2

解析結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

66

基板の表面起伏が摩擦係数に与える影響 . . . . . . . . . . . . . . . . .

75

5.2.1

解析条件 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

75

5.2.2

解析結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

75

第6章 結 論

86

参考文献

89

学術論文・学術講演

91

謝辞

99

第1章
緒論
炭素の同素体としては黒鉛やダイヤモンドが古くから知られていたが，近年の科学技術の発達に伴い，ナノメートルレベルで C60 [1] やカーボンナノチューブ (CNT)[2] といった新たな同素体が次々と発見され，現在ではそれらの形態や物性の研究が盛んに行われている．その一つに 1992 年の Nature 誌で Ugarte より初めて報告されたオニオンライクカーボン (Onion-Like Carbon:OLC)[3] がある (図 1.1)．OLC は球殻構造を有し，フラーレンが分子量の小さいものから順に入れ子状に積層した物質であり，透過型電子顕微鏡 (Transmission Electron Microscope:TEM) で観察するとタマネギのよ
うに見えることから命名された．その直径は数 nm〜数十 nm にも及び，C60 の直径約0.7[nm] と比較すると非常に大きな値である．こうした形状から，グラファイトの性質である熱伝導性や電気伝導性を有しており [4] ，特に固体潤滑剤 [5] としての応用が期待されている．
実験によるアプローチでは，Kuznetsov は OLC が発見された当初から，その物理的特性や電気的特性といった基礎的な特性の評価を行っている [4],[6] ．Hirata，Igarashi ら
はシリコンウェハ上に散布させた OLC に対してボールオンディスク試験機を用いて摩擦試験を行い，大気中，真空雰囲気中においても摩擦係数が 0.1 以下であったことを報告している [7] ．Joly-Pottuz，Ohmae らのグループは OLC を基油に添加すること
で，その潤滑性能を向上させることが可能であることを示している [8]−[10] ．この他にも，特にその生成方法については様々な研究が行われてきた [11]−[13] ．
計算材料科学の分野では，C60 や CNT などは分子動力学を初めとする電子・原子シ
ミュレーションの格好のターゲットとして盛んに研究がなされてきた．山口，丸山ら

1

第1章

緒論

2

は，C60 生成における温度の影響やその生成メカニズムを調べており [14] ，CNT の多様な物性の一つである熱物性に着目し，伝熱特性の評価も行っている [15] ．尾方，渋谷らは，無欠陥の CNT ではなく，五員環や七員環のような六員環以外の欠陥を有す CNT
についての熱伝導率の評価を分子動力学法により行っている [16] ．Hirai, Nishimaki らは, ピンホール欠陥を持つ CNT の機械的特性を調べている [17] ．Deguchi，Yamaguchiらは単層 CNT に高温化で引張応力を与えることで Stone-Wales 欠陥を発生させ，引張時におけるそうした欠陥の挙動について評価している [18] ．西村らは多層 CNT の解析として，二層構造を持つ CNT の圧縮特性を分子動力学法により検討している [19] ．また，簡易的な原子間ポテンシャルでの検討だけでなく，Tight Binding 法や第一原理計
算など，電子状態をも考慮した解析も行われている [20],[21] ．しかしながら，高次のフラーレンや OLC になると原子数が極めて多くなるため，これらの高精度なシミュレー
ションによる検討は限られている．

TEM で観察される OLC の中心のフラーレンの直径は C60 の直径に近いこと，そうして重なり合うフラーレンの層間距離はグラファイトの層間距離 (0.334[nm]) とほぼ等しく，この値は n 層目の炭素原子数を 60n2 とした場合の層間距離に近いこと，などが報告されている [22] ．また，一般的にフラーレンは五員環と六員環から構成されると考えられており [23] ，この場合単体での安定構造は正二十面体形状とされている．

OLC についてもそのように考えられているが，実際 TEM で観察される OLC の多く
は球殻形状を有しているため，五員環と六員環だけでなく七員環等の欠陥を有しているものと考えられる．このように原子レベルでの構造が未だ解明されていないために，

OLC についてのシミュレーションはその多くが構造安定性に関するものである．Saito
らは，高次のフラーレンにおける欠陥の遷移過程について，いくつかのパターンを提案している [24] ．Terrones らは，60n2 という原子数にとらわれず，正二十面体形状やStone-Wales 欠陥を導入して球殻形状を実現した様々なフラーレンのカイラリティや構造について研究を行っている [25] ．Wang，Chang らは，Stone-Wales 欠陥を元に新たな欠陥構造を提案し，これを導入することで安定な球殻構造の作成に成功している [26] ．
他にも構造安定性に関する研究はいくつか行われているが [27],[28] ，そうして作成した

OLC を対象としたシミュレーションは依然として少なく，行われていても C60 ，C240 ，C540 の三層からなる OLC を対象としたものがほとんどである [29],[30] ．

第1章

緒論

3

本研究では OLC を平面上に並べた仮想的な薄膜構造に対して，押し込みやスクラッ
チを分子動力学法シミュレーションを用いて行い，OLC 薄膜の摩擦特性について評価することを主たる目的とする．単一フラーレン・OLC の違い，OLC の直径，圧子形状，圧子表面粗さ等に着目した様々なシミュレーションを行う．
第 2 章では解析手法の基礎として，分子動力学法を簡単に説明し，分子動力学計算で最も重要となるポテンシャルエネルギーについて述べる．また，大規模計算を行うための高速化手法を示し, 最後に OLC を構成するフラーレンの幾何形状について説明する．
第 3 章では径の大きさの異なる単一フラーレンに対して，初期構造緩和，圧縮，摩擦のシミュレーションを行う．初期構造緩和シミュレーションにより緩和後のフラーレンの構造を示し，その半径について実験値との比較を行う．圧縮シミュレーション
では頂点保持による圧縮とダイヤモンド壁による圧縮の 2 通りを行い，フラーレンの圧縮特性や圧縮方法による力学応答の違いを議論する．摩擦シミュレーションではフラットな面を持つダイヤモンド壁を用いて，フラーレンの径の大きさや押し込み量による摩擦係数の違いについて議論する．
第 4 章では径の大きさの異なる単一 OLC に対して，3 章と同様のシミュレーションを行い，それぞれ径の大きさが与える影響やフラーレンとの違いについて議論する．
第 5 章ではフラーレン・OLC を平面上に並べた薄膜構造に対し，圧子や基板の形状を変えた摩擦シミュレーションを行い，対象 (フラーレンか OLC か)，径の大きさ，圧子形状・表面凹凸などが摩擦係数に与える影響について考察する．
最後に, 第 6 章で本研究の総括を述べる.

10 nm
Fig.1.1 TEM image of OLC.[3]第2章
解析手法の基礎2.1

分子動力学法

分子動力学法 (molecular dynamics method，略して MD 法) は，系を構成する各粒子についてニュートンの運動方程式

mid2 r i
= Fi
dt2

(2.1)

を作成し，これを数値積分することにより粒子の軌跡を求める方法である．ここで，mi ，r i はそれぞれ粒子 i の質量および位置ベクトルである．粒子 i に作用する力 F i は，系のポテンシャルエネルギー Φtot の各位置における空間勾配として次式により求められる．

Fi = −

∂Φtot
∂r i

(2.2)

式 (2.1) の数値積分には，Verlet の方法，予測子–修正子法等がよく用いられる．本
研究では，以下に示す Verlet の方法を用いた．
時刻 t + ∆t と t − ∆t での粒子 i の位置ベクトル r i (t ± ∆t) を Taylor 展開すると，()
dr i (t) (∆t)2 d2 r i (t)
3
+
+
O
(∆t)
r i (t + ∆t) = r i (t) + ∆t
dt
2
dt2
2 2
(
)
dr i (t) (∆t) d r i (t)3r i (t − ∆t) = r i (t) − ∆t
+
+
O
(∆t)
dt
2
dt2

(2.3)
(2.4)

となる．ここで，v i を時刻 t における粒子 i の速度とすると，

dr i= v i (t)
dt
4

(2.5)

第2章

解析手法の基礎

5であり，式 (2.1) と式 (2.5) を式 (2.3) と式 (2.4) に代入すると，()(∆t)2 F i (t)+ O (∆t)3
2
mi
2
(
)
(∆t) F i (t)+ O (∆t)3r i (t − ∆t) = r i (t) − ∆tv i (t) +2mir i (t + ∆t) = r i (t) + ∆tv i (t) +

(2.6)
(2.7)

となる．両式の和と差をとると，

r i (t + ∆t) + r i (t − ∆t) = 2r i (t) + (∆t)2()
F i (t)+ O (∆t)4mi(

r i (t + ∆t) − r i (t − ∆t) = 2∆tv i (t) + O (∆t)3

)

(2.8)
(2.9)

が得られる．これより，時刻 t + ∆t での位置ベクトルと t での速度は

r i (t + ∆t) = 2r i (t) − r i (t − ∆t) + (∆t)2v i (t) =()
F i (t)+ O (∆t)4
mi()
1{r i (t + ∆t) − r i (t − ∆t)} + O (∆t)22∆t(2.10)(2.11)と求められる．t + ∆t での座標を求めるには 2 つの時刻 t と t − ∆t での座標が必要である．初期の計算 (t = 0) では，t = ∆t での座標 r i (∆t) は式 (2.6) と初速度から得ることができる．

2.2

原子間ポテンシャル

粒子に作用する力は系のポテンシャルエネルギ−により決定される．そのため，分子動力学法においてはポテンシャルの選定が重要になる．解析の対象となる物質や解析条件に適当となるようなポテンシャル関数を決定しなければならない．

2.2.1

Brenner ポテンシャル

本解析では，シリコンに対する Tersoﬀ 型ポテンシャルが Brenner によりグラファイ
トにフィッティングされたものを用いる [31] ．
シリコンは常温常圧においてダイヤモンド構造を持つが，炭素の場合ダイヤモンドとグラファイトという 2 つの安定構造がある．したがって，グラファイトの sp2 結合と，ダイヤモンド構造の sp3 結合の違いを表現することが重要になる．

第2章

解析手法の基礎

6

系のエネルギーはΦtot =∑ ∑[]fc (rij ) VR (rij ) − B¯ij VA (rij )(2.12)i j(>i)と表される. rij は粒子 i, j 間の距離を示している．ここで VR (rij ) は斥力をあらわす項であり，−B¯ij VA (rij ) は引力をあらわす項である．また fc (rij ) はポテンシャルエネ
ルギーの打ち切り（カットオフ）を滑らかにするための項である．各項はそれぞれ[]√De
exp −β 2S (rij − Re )
S−1


√
S · De2
VA (rij ) =exp −β(rij − Re )
S−1
S

VR (rij ) =



1,

 [


fc (rij ) =






(

1 + cos

π (rij − R1 )R2 − R1(2.13)(2.14)
rij < R1)]

/2, R1 < rij < R20,

(2.15)

rij > R2である. B¯ij は，粒子 i ，粒子 j 以外に粒子 k も含めた 3 個の原子によって定められる項であり,

1B¯ij = (Bij + Bji )2(2.16)Bij , Bji はそれぞれ原子 i, j を中心とした加算となり一般に Bij ̸= Bji である.
Bij = [1 +∑G(θi )fc (rik )]−δ(2.17)

G(θj )fc (rjk )]−δ(2.18)

k̸=i,jBji = [1 +∑k̸=i,j

ここで rik は粒子 i, k 間の距離, rjk は粒子 j, k 間の距離を示している. θi は原子 i を中心とする j − i − k の内角, θj は原子 j を中心とする i − j − k の内角である. 微分のため,∑G(θi )fc (rik ) を ζij ,∑G(θj )fc (rjk ) を ζji とする. G(θ) は
[

c20
c20
G(θ) = a0 1 + 2 − 2
d0 d0 + (h + cos θ)2](2.19)

である.
式 (2.12) ∼ 式 (2.19) で用いられるパラメータは C2 分子，独立したグラファイト
シートおよびダイヤモンドそれぞれの結合エネルギと平衡状態における結合距離，仮第2章

解析手法の基礎7

想的な単純格子および f cc 構造の結合エネルギ，グラファイトの斜方六面体構造から
ダイヤモンド構造への相変態における障壁エネルギの計算等によりフィッティングされている．それらのパラメータを表 2.1 に示す．

Table 2.1 Potential parameters for Brenner potential.
De [eV]
6.0
Re [nm]0.139
β [nm−1 ]
0.21
S
1.22
h
1.0
a0
0.00020813
c0
330.0
d0
3.5
R1 [nm]0.17
R2 [nm]0.20

シリコンに対するパラメータの場合と最も異なることは，h の値が 1 になっている点である．h が 0 の場合はより小さな結合角で安定するのに対し，h を 1 にすることで結合角はより大きく，すなわち 180 [deg.] に近づこうとする．このことが，グラ
ファイトシート内の平面構造を維持する駆動力になる．また Brenner ポテンシャルには, 炭素原子間距離の値に重点を置きクラスタの形成に最適化されたパラメータ 1 と,炭素間に作用する力の値に重点を置き物性の測定に最適化されたパラメータ 2 が存在
する. パラメータは力の再現を重視したパラメータ 2 を用いて計算を行った.
結合距離の変化に対する各結合角でのポテンシャルおよびポテンシャルが最小の値をとる時の結合角の変化に対する結合距離の変化を図 2.1，2.2 に示す．図 2.1 より, ポ
テンシャルの最小値が結合角の減少に伴い結合距離が大きくなる方向へ推移しているのが確認できる. また図 2.2 より結合角 60∼90[deg.] において, 結合距離が増減し変化が滑らかではないのは, 原子 j, k での結合長に対するカットオフ半径による.

第2章

解析手法の基礎20

E [eV]

θ = 30θ = 60θ = 120
θ = 180
10

0

0.1

0.12

0.14

0.16

0.18

0.2

r [nm]

Fig.2.1 Relationship between potential energy and bond length.

0.2

r [nm]

0.18

0.16

0.14
60

120

180

θ [deg]Fig.2.2 Relationship between stable bond length and bending angle.

8

第2章

2.2.2

解析手法の基礎

9

Lennard-Jones ポテンシャル

グラファイトの層間の Van der Waals ポテンシャルは，次式のような Lennard-Jones型ポテンシャルで表される [32] ．
∑ ∑

(

σ
Φtot =4ϵ riji j(̸=i)

)12

(

σ
−
rij

)6 


(2.20)

第 1 項は斥力，第 2 項は引力を表し，式 (2.20) で用いられるパラメータを表 2.2, ポテ
ンシャル曲線を図 2.3 に示す．

Table 2.2 Potential parameters for Van der Waals.

ϵ
σ0.004783 [eV]0.3345 [nm]

–3

×10

E [eV]

5.0

0.0

–5.0
0.3

0.4

0.5

0.6

r [nm]

Fig.2.3 Relationship between Van der Waals potential and atomic distance.

第2章

2.2.3

解析手法の基礎

10力の表式

本解析で用いる各ポテンシャルについて，式 (2.2) に基づいて各ポテンシャル関数
より原子間力を導く．まず，Brenner ポテンシャルについて原子 ij 間のエネルギー寄与は

Φij = VR (rij ) −Bij + Bji
VA (rij )2(2.21)

である. これの r i , r j , r k による微分

Fi = −

∂Φij
∂Φij
∂Φij
,Fj = −
,Fk = −
∂r i
∂r j
∂r k

を, j > i 全ての組について加算すれば原子に働く力が求まる. 以下で個々の原子にかかる力を示す.
原子 i の力は]r[[∂Bij ∂Bji1
F i = B¯ij VA′ (rij ) − VR′ (rij )+ VA (rij )
+
rij
2
∂r i
∂ri

]

ij

(2.22)

である. ここで ′ は r による微分を表す. VR (rij ), VA (rij ) の r ij での微分は
]r
[
√
√
De
ij
VR′ (rij ) = −β 2S
exp −β 2S (rij − Re )S−1
rij

(2.23)

2 S · De2r ij
exp −β(rij − Re )
S S−1
S
rij

(2.24)√VA′ (rij ) = −β





√

である. fc (rij ) の r ij での微分は
(



 −1

∂fc (rij )
=

∂r ij


π
π (rij − R1 )sin
2 R2 − R1
R2 − R1)0: R2 > rij > R1

(2.25): rij < R1 , R2 < rijである. 以降では Bij , Bji の r i での微分を考える.

∂Bij
∂ζij
= −δ(1 + ζij )−δ−1
∂r i
∂r i

(2.26)(

中心∂ζij i∑r ik∂ cos θi=fc′ (rik )G(θi )+ fc (rik )G′ (θi )
∂r i
rik
∂r i
k̸=i,j)(2.27)G′ (θ) は G(θ) を cos θ で微分したもので以下である.
[∂G(θ)
2c20 (1 + cos θ)G (θ) =
= a02
∂ cos θ[d20 + (1 + cos θ)2 ]′](2.28)

第2章

解析手法の基礎11

cos θ の位置ベクトルでの微分を考える. 原子 j, i, k の内角であるから
r ij · r ik
rij rik()()∂ cos θi1
cos θi r ij
1cos θi r ik
=
−
+
−
∂r i
rik
rij
rij
rij
rik
rik
cos θi =(2.29)(2.30)

次に Bji の微分を考える.

Bji = (1 + ζji )−δ より

∂Bji
∂ζji
= −δ (1 + ζji )−δ−1
∂r i
∂r i

(2.31)

中心∂ cos θj∂ζji j∑=fc (rjk )G′ (θj )
∂r i
∂r i
k̸=i,j(2.32)

r ji · r jk
rji rjk∂ cos θj
1 r jk cos θj r ji
=−
+
∂r i
rji rjkrji rjicos θj =(2.33)(2.34)

以上を整理すると[

F i = B¯ij VA′ (rij ) − VR′ (rij )

]r

ij

rij

(

)

θi
中心
(
) i∑
fc (rik )G′ (θi ) r1ik − cos
1
−δ−1
rij )
(
+ VA (rij ) −δ (1 + ζij )cos θi1′′2
k̸=i,j fc (rik )G(θi ) + fc (rik )G (θi ) rij − rik
r ji
θj
中心
(
) j∑
fc (rjk )G′ (θj ) cos
1
−δ−1
rji
rji
+ VA (rij ) −δ (1 + ζji )
r jk
1
′
2
−f
(r
)G
(θ
)
c jk
j rji
k̸=i,jrjk
原子 j の力はFj =[B¯ij VA′ (rij ) − VR′ (rij )

]

(

r ij
−
rij

)

[

1
∂Bij ∂Bji+ VA (rij )
+
2
∂r i
∂r i

r ij

]

(2.35)

である. Bij の微分の ζij について(中心
∂ cos θi∂ζij i∑=fc (rik )G′ (θi )
∂r i
∂r j
k̸=i,j)(2.36)∂ cos θi
1 r ik cos θi r ij
=−
+
∂r j
rij rikrij rij(中心
r jk∂ cos θj∂ζji j∑fc′ (rjk )G(θj )=+ fc (rjk )G′ (θj )
∂r j
rjk
∂r j
k̸=i,jrij

r ik(2.37))(2.38)rik第2章∂ cos θj
=
∂r j

(

1
cos θj−rjkrji

以上を整理すると

Fj =[B¯ij VA′ (rij ) − VR′ (rij )

]

(

r ij
−
rij

)

(

r ji
1
cos θj
+
−
rji
rji
rjk

解析手法の基礎
)

r jk
rjk

12

(2.39))

r ij

中心 f (r )G′ (θ ) cos θi
(
) i∑
1
c ik
i rij
−δ−1
+ VA (rij ) −δ (1 + ζij )1′2
k̸=i,j −fc (rik )G (θi ) rij

rij

r ik
rik

(

)

θj
r ji
fc (rjk )G′ (θj ) r1jk − cos
1
−δ−1
r
rji
ji
(
)
+ VA (rij ) −δ (1 + ζji )r jk
cos θj1′′2
k̸=i,j fc (rjk )G(θj ) + fc (rjk )G (θj ) rji − rjkrjk
原子 k の力は][1∂Bij ∂Bji+(2.40)
F k = VA (rij )
2
∂r k
∂r k

(

中心
) j∑

である. Bij , Bji の微分の ζij , ζji については
(

(

)

中心
∂ζij i∑r ik∂ cos θi=fc′ (rik )G(θi ) −+ fc (rik )G′ (θi )
∂r k
rik
∂r k
k̸=i,j

(

(

中心
r jk
∂ζji j∑fc′ (rjk )G(θj ) −
=
∂r k
rjk
k̸=i,j))∂ cos θj+ fc (rjk )G (θj )∂r k(2.41))′

である. 以上を整理すると
中心
(
) i∑
−fc (rik )G′((θi ) r1ik)
1
F k = VA (rij ) −δ (1 + ζij )−δ−1
cos θi′′2k̸=i,j −fc (rik )G(θi ) + fc (rik )G (θi )
rik
中心
(
) j∑
−fc (rjk )G′ (θj ) r1jk
1
−δ−1
+ VA (rij ) −δ (1 + ζji )cos θj′′2k̸=i,j −fc (rjk )G(θj ) + fc (rjk )G (θj ) rjk
同様に，Van der Waals ポテンシャルについては，

Fi = −∂Φ(rij )
∂r ij
∑



12
=−
4ϵ −
rij
j(̸=i)(

F j = −F i

σ
rij

)12

6
+
rij

(

σ
rij

)6 


r ij
rij

(2.42)

r ij
rij

r ik
rik
r ji
rji
r jk
rjk

(2.43)(2.44)

なお，添字 i，j ，k は原子をあらわし，rij , rik はそれぞれ粒子 i, j 間，i, k 間の距離，

θi は結合 rij , rik 間の角度を表している．

第2章

2.3

解析手法の基礎13

高速化手法

原子数 N の系において粒子間の全相互作用を評価すると，1step に N × (N − 1) 回の計算が必要となり，N が大きくなると極めて膨大な計算量となる．実際には，一定距離以上離れた粒子は影響を及ぼさないので，作用を及ぼす範囲 (カットオフ半径 rc ) 内の粒子からの寄与を効率よく計算することにより高速化できる．従来よく用いられてきた高速化手法に粒子登録法がある．これは，図 2.4 に示したように，rc よりひとまわ
り大きい半径 rfc 内の粒子をメモリーに記憶し，その中で rc 内の相互作用を評価する方法であり，N × (rc 内粒子数 ≪ N − 1) に計算負荷が減少される．しかし，粒子登録法では rfc 半径より外の粒子が rc 内に達すると力の評価が適切でなくなるので，一定のステップ毎に登録粒子の更新 (N × (N − 1) 回の探査) を行わなければならない．このため，系がある程度の規模以上になると，粒子登録による高速化は登録更新の計算負荷により打ち消される．

rfc

rc

Fig.2.4 Schematic of bookkeeping method.

第2章

解析手法の基礎

14

別の高速化手法としてブロック分割法がある．図 2.5 に示すように，シミュレートする系をカットオフ距離程度の格子状に分割し，各ブロックに属する粒子をメモリーに記憶する．着目している粒子に作用する力を評価する際には，その粒子が属するブロックおよび隣接するブロックから相互作用する粒子を探索して行う．粒子が属する
ブロックは，粒子の位置座標をブロックの辺長 bx，by で除した際の整数により判断できるので，ブロック登録時の計算負荷は粒子数 N のオーダーとなる．したがって，粒子登録法では登録更新の負荷が大きくなるような大規模な系でも高速化が可能である．

y

by
0

x

bx

Fig.2.5 Schematic of domain decomposition method.

第2章

2.4

解析手法の基礎15

速度スケーリング法

分子動力学解析における温度制御には一般的には速度スケーリング法が用いられる．
この方法は，統計熱力学より導かれる式 (2.45) を用いて，以下のように制御する．

1 α α α 3m vi vi = k B T
22(2.45)

mα : 粒子αの質量
viα : 温度 T での粒子αの速度kB : Boltzmann 定数 = 1.38 × 10−23 [J/K]
目標の温度 T0 における原子 α の速度を viα0 とおくと viα0 は式 (2.46) のように表される．

(

viα0 =

3kB T0mα)0.5(2.46)

同様に，温度 T の時の原子 α の速度は式 (2.47) のように表される．
(

viα =

3kB T
mα

)0.5

(2.47)

よって，式 (2.46) と式 (2.47) より以下の式が得られる．
(

viα0
T0
=
α
vi
T

)0.5

(2.48)

つまり，系の温度を T から T0 にするには，式 (2.48) の右辺を現在の速度に掛けてやればよい．ただ，これだけでは原子配置に反映されないので，Verlet 法における ∆r αi (t+∆t)√（式 2.49）を

T0 /T ∆r αi (t + ∆t) と置き換える必要がある．

∆r αi (t + ∆t) = r αi (t + ∆t) − r αi (t) = r αi (t) − r αi (t − ∆t) + (∆t)2

F αi (t)(2.49)mα

平衡状態では，能勢の方法 [33] など外部との熱のやりとりをする変数を考慮した拡張系の分子動力学法によって得られるカノニカルアンサンブルに一致することが示されて
いる．

第2章

2.5

解析手法の基礎

16

フラーレンモデルの幾何形状

一般的にフラーレンは五員環と六員環から構成されていると考えられている [23] ．この場合，Euler の定理に従うと五員環の数は 12 個であり，その五員環は孤立五員環則(Isolated Pentagon Rule) と経験的な C60 の五員環の位置から，正二十面体の頂点の位置にあるものと考えられる [23] ．しかし，これまでの検討でこのルールに従って原子数が 60n2 で表されるフラーレンを作成した場合，緩和を行うと五員環の部分が外側に突
き出し，正二十面体形状になることがわかった [34] ．このようにフラーレンの形状は五員環の位置によって決まるため，球殻形状を実現するためには五員環の数を増やす必要があると考える．そこで今回は五員環を含む形状欠陥として，図 2.6 に示すような

6-6-6-6 形状からの遷移で知られる 5-7-5-7 形状欠陥 (Stone-Wales 欠陥) を導入することとした．導入過程においては孤立五員環則を順守し，対称性も考慮した．また，今回は五員環，六員環，七員環のみでフラーレンを構成するものとした．上述のルール
により C240 は七員環を導入することができず，今回 C60 及び C240 には導入を行っていない．導入例として，図 2.7 に C540 ，図 2.8 に C960 の欠陥を導入する前後の原子配置について，半球分を平面上に展開したものを示す．図の赤い部分は各フラーレンの最小構成単位を示したものである．

(a) Before inserting(b) After inserting

Fig.2.6 Snapshots of Stone-Wales defect.

第2章(a) Before inserting

解析手法の基礎(b) After inserting

Fig.2.7 Snapshots of expanded C540 .

(a) Before inserting

Fig.2.8 Snapshots of expanded C960 .

(b) After inserting

17

第3章
単一フラーレンの解析
3.1

初期構造緩和シミュレーション

3.1.1

解析条件

原子数が 60n2 の形で表わされるフラーレンの，n=1〜6 の場合を解析対象とし，全方向自由境界条件 (真空) のもとで 20000[fs] の緩和計算を行った．炭素原子の初期原子間距離は, 今回使用した原子間ポテンシャルの平衡粒子間距離である 0.145[nm] としている. 各原子の初速度は Maxwell Boltzmann 分布に従って乱数で与えており，温度は10[K]，積分の時間ステップは全て 1.0[fs] である．

3.1.2

解析結果

緩和後の各種のフラーレンについて，n の値，半径 r ，そして，OLC にした場合の層間距離に当たる，n の値が 1 少ないフラーレンとの半径の差を表 3.1 に示す．ここで，半径 r は，原点から各原子までの距離の平均値で求めている．表 3.1 を見ると，C60 の半径は実験値の約 0.35nm とほぼ等しく，半径の差 (層間距離) も多少のばらつきはあるもののグラファイトの層間距離 (0.336nm) に近い．図 3.1 に C1500 の緩和計算後の構造を例として示す．図 3.1 を見ると，欠陥の導入により五員環の数が増えたことで，多少の凹凸はあるものの球殻に近い形状となっているのがわかる．他のフラーレンでも同様の傾向が見られた．

18

第3章

単一フラーレンの解析

Table 3.1 Fullerene parameter n, radius and radius diﬀerence after relaxation.
fullerenenradius [nm]radius diﬀerence [nm]

C60
C240
C540
C960
C1500
C2160

1
2
3
4
5
6

0.364
0.721
1.078
1.455
1.805
2.170

−
0.357
0.357
0.377
0.350
0.365

z

x

1

y

2

(a) top view(b) side view ( 1 )(c) side view ( 2 )Fig.3.1 Snapshots of C1500 after relaxation.

19

第3章

3.2

圧縮シミュレーション

3.2.1

頂点保持による圧縮

単一フラーレンの解析

20

解析条件
前節の n=3〜6 のフラーレンを対象に，上下頂点の五員環をつかみ部として固定しながら圧縮するシミュレーションを行った．圧縮はひずみ制御で行い，原子間の距離を均等に縮めることでひずみを与えた．圧縮ひずみは毎ステップ増加させており，ひずみ速度に換算すると 5.0 × 10−6 [/fs] となる．各原子の初速度は Maxwell Boltzmann 分
布に従って乱数で与えており，温度は 10[K]，積分の時間ステップは全て 1.0[fs] である．
解析結果
図 3.2 に C540 の圧縮シミュレーション中の応力−ひずみ関係と変形の様子をまとめて示す．図 3.2(ii) のスナップショットは図 3.2(i) 中に矢印で示した各点にそれぞれ対応している．図 3.2(i) の ε < 0.05 の圧縮初期を見ると，(a)ε = 0.005 でわずかに上昇した応力は (b)ε = 0.03 までわずかに下降し，その後は (b)ε = 0.03 から (c)ε = 0.05 まで
あまり上昇せずに進む．このときの様子を図 3.2(ii) で確認すると，(b) では下から見た
図中に緑の丸印をつけて記したつかみ部の五員環が内側に向かってへこんでいること，また (b)，(c) 間ではあまり変化がないことがわかる．図 3.2(i) の 0.05 < ε < 0.22 の圧縮中期を見ると，応力は (c)ε = 0.05 から (d)ε = 0.17 に向かって著しく上昇し，その後ピークを示して急激に下降する．この応力ピーク前後の構造を図 3.2(ii) で確認すると，(e) を下から見た図中に青丸をつけて示した五員環群が球の内側にへこむことで座屈していた．
図 3.3 に C960 の圧縮シミュレーション中の応力−ひずみ関係と変形の様子をまとめて示す．図 3.3(i) の応力−ひずみ関係を見ると，圧縮直後は応力上昇するものの，(a) →(b) で低下し，(b) → (c) 間ではあまり応力上昇せずに進む．このときの様子を図 3.3(ii)で確認すると，先の C540 と同様に，(b) を上部から見た図に示したつかみ部の五員環が内側に向かってへこんでいた．図 3.3(i) の 0.095 < ε < 0.3 の圧縮中期では，応力は(c) → (d) に向かって上昇するが，C540 ほど顕著ではない．(d) の応力ピーク後の構造

第3章

単一フラーレンの解析

21を図 3.3(ii) の (e) を上部から見た図で確認すると，青丸で示したつかみ部の周囲の五員環群ではなく，赤丸で示したさらに一つ外側の五員環群が球の内側にへこむことで座屈していた．
図 3.4 に C1500 の圧縮シミュレーションの結果を示す．図 3.4(i) を見ると，圧縮直後にわずかに上昇→低下→再上昇して応力ピーク，という全体的な傾向は同じであるが，のこ歯状のゆらぎが大きく表れている．(a) → (c) の初期応答を図 3.4(ii) で確認すると，やはりつかみ部の五員環が内側に向かってへこんでいる．また (d) → (f) の応力ピーク前後の変化についても，C960 同様頂点から第二周目に相当する赤丸で示した五員環群で変形を生じていた．
図 3.5 に C2160 の圧縮シミュレーションの結果を示す．図 3.5(i) の初期応答 (a) → (c)は長く複雑であり，これまでのフラーレンに比べ，あまり応力上昇せずに進む．また
圧縮直後でもほとんど応力上昇せず，これまで述べてきたようにつかみ部の五員環が内側にへこみ，外殻が「引っ張られ」ることで負の圧縮応力，すなわち引張応力を示している．このときつかみ部の五員環は上下共へこんでいたが，その後は上部のへこみのみ進行しており，(d) の応力ピーク後には，中心から第三番目の五員環群 (図 3.5(ii)の (e) の紫丸) が球の内側にへこんでいた．
図 3.6 に応力−ひずみ関係をまとめて再掲する．C540 は ε = 0.17 において σ =2.21[GPa]，C960 は ε = 0.23 において σ = 1.14[GPa]，C1500 は ε = 0.175 において
σ = 1.12[GPa]，C2160 は ε = 0.245 において σ = 0.52[GPa] の最大応力を示しており，その値は高次のフラーレンほど低くなっている．これは原子数が増えたことで bond 数が増加し，変形の自由度が増加したためと考えられる．また図 3.7 に最大応力を示し
た後に応力急減した時の変形部分をまとめて再掲する．図中に記してある点線は，各フラーレンにおける最小構成単位 (図 2.7 参照) のつなぎ目を表している．図 3.7 を見ると，C540 ，C1500 では座屈に寄与する五員環群が点線上に存在するが，C960 ，C2160 では点線よりも外側に，すなわち別の構成単位に属する五員環で変形している．最大応力を生じるひずみは，C540 ，C1500 は C960 ，C2160 より小さいが，応力は高い．構成単位の五角形状に五員環が存在する場合とそうでない場合でこのような傾向を示す可能
性がある．

第3章

単一フラーレンの解析

22

Compressive stress, -σ , GPa

C540
(d)
2

(e)
1

(a)
(f)
0

(b)

(c)

0

0.1

0.2

0.3

Compressive strain, | ε |(i) Stress - strain curve of C540 under compression.

(a) ε =0.005(b) ε =0.03

close-up bottom view(c) ε =0.05zxy(d) ε =0.17(e) ε =0.19

close-up bottom view(f) ε =0.215(ii) Snapshots of C540 under compression.
Fig.3.2 Compression of C540 by holding the ﬁve-membered ring at the top
and bottom.

第3章

単一フラーレンの解析

23

Compressive stress, -σ , GPa

C960
2

(d)
1

(a)
(e)
(f)

0

(b) (c)
0

0.1

0.2

0.3

Compressive strain, | ε |(i) Stress - strain curve of C960 under compression.

(a) ε =0.045(b) ε =0.055

close-up top view(c) ε =0.095zxy(d) ε =0.23(e) ε =0.245

close-up top view(f) ε =0.27(ii) Snapshots of C960 under compression.
Fig.3.3 Compression of C960 by holding the ﬁve-membered ring at the top
and bottom.

第3章

単一フラーレンの解析

24Compressive stress, -σ , GPa

C1500
2

(d)
1

(a)

(e)

0

(f)

(b)

(c)

0

0.1

0.2

0.3

Compressive strain, | ε |(i) Stress - strain curve of C1500 under compression.

(a) ε =0.015(b) ε =0.025

close-up bottom view(c) ε =0.06zxy(d) ε =0.175(e) ε =0.18

close-up bottom view(f) ε =0.295(ii) Snapshots of C1500 under compression.
Fig.3.4 Compression of C1500 by holding the ﬁve-membered ring at the top
and bottom.

第3章

単一フラーレンの解析

25

Compressive stress, -σ , GPa

C2160
2

1

(d)
(a)
(e)

0

(f)

(c)
(b)
0

0.1

0.2

0.3

Compressive strain, | ε |(i) Stress - strain curve of C2160 under compression.

(a) ε =0(b) ε =0.01

close-up top view(c) ε =0.125zxy(d) ε =0.245(e) ε =0.255

close-up top view(f) ε =0.3(ii) Snapshots of C2160 under compression.
Fig.3.5 Compression of C2160 by holding the ﬁve-membered ring at the top
and bottom.

Compressive stress, -σ , GPa

第3章

単一フラーレンの解析

C540
C960
C1500
C2160

2

1

0

0

0.1

0.2

0.3

Compressive strain, | ε |Fig.3.6 Stress - strain curves of fullerenes under compression.

(a) C540(b) C960(c) C1500(d) C2160

Fig.3.7 Top view of fullerenes after stress peak.

26

第3章

3.2.2

単一フラーレンの解析

27

ダイヤモンド壁による圧縮

解析条件
頂点保持による圧縮同様，n=3〜6 のフラーレンを対象に，剛体ダイヤモンド壁を用いて圧縮シミュレーションを行った．まず xy 方向周期境界条件下で，ダイヤモンド壁を各フラーレンの下 0.44[nm] の位置にセットし，さらに 20000[fs] の緩和計算を行う．
その後，図 3.8 に示すように新たなダイヤモンド壁をフラーレンの上 1[nm] の位置に配置し，これを下降させることで圧縮シミュレーションを行う．押し込みは対象の直径 D に対して 3/4D まで行った．ダイヤモンド壁の厚さは 0.7[nm] とし，その下降速
度は 1.0 × 10−4 [nm/fs] とした. セルの周期方向の寸法は 7.1[nm] であり，圧縮してつ
ぶれてもフラーレン同士は相互作用しない大きさであることを確認している．
フラーレンとダイヤモンド壁との間では，共有結合的な反応は考慮せず，Van der

Waals 力のみを考慮した．各原子の初速度は Maxwell Boltzmann 分布に従って乱数で与えており，温度は 10[K]，積分の時間ステップは全て 1.0[fs] である．

1nm

z
x
y
7.1nm
Fig.3.8 Schematic of wall compression (C1500 ).

第3章

単一フラーレンの解析

28

解析結果
図 3.9 に C540 の圧縮シミュレーション中に生じる応力と変位の関係と変形の様子をまとめて示す．図 3.9(ii) のスナップショットは図 3.9(i) 中に矢印で示した各点にそれぞれ対応している．図 3.9(i) の横軸の変位は初期位置を 0 としており，縦軸の応力は上下端の五員環を固定した前節の図 3.6 の 2 倍程度のスケールとなっている．図 3.9(i)の応力を見ると，(a) の変位 0.6[nm] 近傍では負の値を示しているが，これは Van der

Waals 力の引力によるものである．本研究で用いたポテンシャルでは，Van der Waals
力は原子間の距離が 0.37[nm] 以上であれば引力を，以下であれば斥力を生じる．このため変位 0.6[nm] 近傍 (1[nm]-0.37[nm]) で応力が上昇に転じている．(a) から急激に上昇を始めた応力は (d) の変位 1.15[nm] 近傍で一つの応力ピークを示す．この応力ピー
ク前後の構造を図 3.9(ii) の (d)，(e) で確認すると，前項の図 3.2(ii) に示したのと同様の五員環群が上下共に球の内側にへこんで筒状の形状となっている．その後，再び上昇に転じた応力の傾きは，(f) の変位 1.7[nm] 近傍においてさらに変化するが，これは図 3.9(ii) の (f) に示すように，へこんだ上下の面同士が作用し始めるためである．
図 3.10 に C960 の圧縮シミュレーション中の応力−変位関係と変形の様子をまとめて示す．図 3.10(i) の応力−変位関係を見ると，(a)，(b) 間で下降しながらも (c) に向かって応力上昇し，その後は (c) から (e) まであまり上昇せずに進む．(a)，(b) 間の応答を図 3.10(ii) で確認すると，前節の図 3.3(ii) の (b) 同様に，最上部の五員環が内側に向かってへこんでいる．応力は (e) で再上昇し，(f) で小さな応力ピークを示した後に，(h) から再び上昇に転じる．(e) で応力勾配が大きくなるのは，先の C540 と同様に筒状の形状になったフラーレンの両端面が相互作用しはじめたためであるが，(f) で応力がわずかに下降したのは，図 3.10(ii) の (f) で示したように互いにずれたためである．(g)→ (h) で応力が再び上昇するのは，互いにずれることができなくなるまで圧縮されたときである．
図 3.11 に C1500 の圧縮シミュレーションの結果を示す．全体的な傾向は C960 によく似ているが，最後の応力上昇前の小さな応力ピークがほとんどない．図 3.11(ii) の (f)を見ると，原子数が多いために対面の原子同士が流動的に動いて蛇腹を形成しているのがわかる．

第3章

単一フラーレンの解析

29

図 3.12 に C2160 の圧縮シミュレーションの結果を示す．全体的な応答は似ているが，他のフラーレンに比べあまり応力上昇せずに進み，(h) から継続的な上昇を始める．一つ目のピーク前後の構造を図 3.12(ii) の (c)，(d) で確認すると，前節の図 3.5(ii) に示したのと同様の五員環群が球の内側にへこんでいる．また二つ目のピーク前後の構造を図 3.12(ii) の (f)，(g) で確認すると，互いに作用していた両端面の均衡がくずれ，下面が再び外側に折れ返している．終盤 (h) の応力が再上昇し始める点では，図 3.12(ii)の (h) にあるように，折れ返った部分が下端に到達しており，その後 (i) では赤い点線で囲った側面のスペースも，(h) では余裕があったのに対し，Van der Waals 相互作用する距離まで狭まっていた．
図 3.13 に全てのフラーレンの応力−変位関係をまとめて再掲する．C2160 を除き，いずれも前節の図 3.6 と同様，ひずみ 0.2 近傍におけるピークとそれぞれ近い値をとっており，フラーレンの構成単位近傍の五員環群の座屈による応答が表れている．また

C540 は変位 1.85[nm]，C960 は変位 2.5[nm]，C1500 は変位 2.9[nm]，C2160 は変位 3.7[nm]から応力が急上昇する．このときの構造 (C540 は図 3.9(ii) の (f)，C960 は図 3.10(ii) の(h)，C1500 は図 3.11(ii) の (f)，C2160 は図 3.9(ii) の (i)) を図 3.14 に再掲する．C540 は上下面が接近して相互作用し始めるときであるが，その他のフラーレンは折りたたまれた内部の隙間がなくなった点に対応する．特に C1500 ，C2160 は原子数が多いために対面の原子同士が流動的に動いて複雑な形状を示している．こうした柔軟な変形により，いずれのフラーレンでも圧縮後に結合が解けて空孔を生じたり，新たに遷移が発生して結合の切り替わりが起こるような変化はなかった．

第3章

単一フラーレンの解析

Compressive stress, -σ , GPa

C540
4

(d)

(f)

2

(e)
(b)

(c)

0

(a)
0

2

4

Displacement of upper wall, nm(i) Stress - displacement curve of C540 under compression.

(a) 0.6nm indentation(b) 0.85nm indentation(c) 0.9nm indentationzxy(d) 1.15nm indentation(e) 1.4nm indentation(f) 1.85nm indentation(ii) Snapshots of C540 under compression.
Fig.3.9 Compression of C540 by diamond wall.

30

第3章

単一フラーレンの解析

31

Compressive stress, -σ , GPa

C960
4

(f)
(h)

2

(c)

(g)

(a)
(e)
(b)

0
0

(d)

2

4

Displacement of upper wall, nm(i) Stress - displacement curve of C960 under compression.

(a) 0.85nm indentation (b) 0.95nm indentation (c) 1.3nm indentation(d) 1.8nm indentationzxy(e) 2.0nm indentation(f) 2.3nm indentation(g) 2.4nm indentation(ii) Snapshots of C960 under compression.
Fig.3.10 Compression of C960 by diamond wall.

(h) 2.5nm indentation

第3章

単一フラーレンの解析

Compressive stress, -σ , GPa

C1500
4

2

(f)

(c)
(a)
(b)

0
0

(e)

(d)
2

4

Displacement of upper wall, nm(i) Stress - displacement curve of C1500 under compression.

(a) 0.8nm indentation(b) 0.85nm indentation(c) 1.4nm indentationzxy(d) 1.7nm indentation(e) 2.2nm indentation(f) 2.9nm indentation(ii) Snapshots of C1500 under compression.
Fig.3.11 Compression of C1500 by diamond wall.

32

第3章

単一フラーレンの解析

Compressive stress, -σ , GPa

C2160
4

2

(f)
(c)

(i)
(g)

(a)
0

(b)
0

(h)

(d) (e)24Displacement of upper wall, nm(i) Stress - displacement curve of C2160 under compression.

(a) 0.8nm indentation (b) 0.85nm indentation (c) 1.9nm indentationzxy(d) 2.2nm indentation(e) 2.6nm indentation(f) 2.9nm indentation(g) 3.0nm indentation(h) 3.3nm indentation(i) 3.7nm indentation(ii) Snapshots of C2160 under compression.
Fig.3.12 Compression of C2160 by diamond wall.

33

Compressive stress, -σ , GPa

第3章

単一フラーレンの解析

34

C540
C960
C1500
C2160

4

2

0
0

2

4

Displacement of upper wall, nm
Fig.3.13 Stress - displacement curves of fullerenes under compression.

z(a) C540(b) C960xy(c) C1500(d) C2160

Fig.3.14 Snapshots of fullerenes at the point just before the drastic stress increase.

第3章

単一フラーレンの解析

3.3

ダイヤモンド壁による摩擦シミュレーション

3.3.1

解析条件

35

n=3〜6 のフラーレンを対象に，剛体ダイヤモンド壁を用いて摩擦シミュレーションを行った．まずダイヤモンド壁による圧縮時と同様の条件で初期緩和を行う．その後，
図 3.15 に示すように，対象の z 軸方向上部から上に 1[nm] の位置に下端表面がくるよ
うにセットしたダイヤモンド壁を下降させて一定量押し込みを行った後，x 軸方向に移動させることで摩擦シミュレーションを行った．ダイヤモンド壁のサイズは対象の直径の 4 倍程度としており，C540 で 10.7[nm]，C960 で 12.8[nm]，C1500 で 14.9[nm]，C2160
で 17.1[nm] である．押し込み量は 1nm+0.1D，1nm+0.3D，1nm+0.5D，1nm+0.7Dの 4 通り (D はフラーレンの直径) とし，x 軸方向の移動量はいずれも 3D とした．壁の下降・スライド速度はともに 1.0 × 10−4 [nm/fs] とした. 他の条件は前節と同じである．

1nm

z
x
y
14.9nm
Fig.3.15 Schematic of scratch simulation (C1500 ).

第3章

3.3.2

単一フラーレンの解析

36

解析結果

C1500 の場合を例として，各押し込み量における摩擦シミュレーション中の摩擦係数−変位関係を図 3.16 に示す．横軸の変位は壁が x 軸方向に移動を開始する点を 0 としている．また縦軸の摩擦係数は壁が受ける x 軸方向の反力の和を，z 軸方向の反力の和で除して求めた．図中の赤い線はその平均値を示す．フラーレンの微小な動きに反応してマイナス方向にもゆらいではいるが，平均値としてはプラスとなった．またそのゆらぎの幅は，1nm+0.1D，1nm+0.3D 押し込みでは ±0.1 程度，1nm+0.5D 押し込
みでは ±0.05 程度，1nm+0.7D 押し込みでは ±0.02 程度と，押し込み量が大きい場合
はゆらぎの幅が小さい．これは他のモデルについても同様である．

1nm+0.1D 押し込みでの各フラーレンの摩擦係数の平均値を表 3.2 に，摩擦前の押し込み状態での形状を図 3.17 に示す．フラーレンの形状を見ると，(a) の C540 は前節で述べたような応力ピーク直後 (大きな座屈を生じた後) であるが，(b)〜(d) の他のフラーレンは応力ピークに達する前である．このときの摩擦係数 (表 3.2) を見ると，いずれも 10−2 程度のスケールであるが，その値は C960 ，C2160 に比べ，C540 ，C1500 の方が比較的高い値を示した．摩擦過程を図 3.18 に C1500 の場合を例として示す．図中の一部の原子の色が変えてあるのは，フラーレンの回転の有無を判別しやすくするためである．図 3.18 を見ると，圧縮により上下がへこんだ状態で，滑るように移動している．赤色で示した原子の位置は，各スナップショットで多少変化しているものの，その移動は摩擦方向への転がりに対応する y 軸を回転軸としたものではなく，ダイヤモ
ンド壁に垂直な z 軸を回転軸とした回転である．

1nm+0.3D 押し込みでの各フラーレンの平均摩擦係数について表 3.3 にまとめた．図3.19 には摩擦前の形状を示している．図 3.19 に示したフラーレンは，いずれも前節で述べたように五員環群が球の内側に向かってへこみ，大きな座屈を生じた後であり筒状の形状をしている．このとき表 3.3 の摩擦係数は，表 3.2 と同程度のオーダーであり，

C540 ，C1500 の方が比較的高い値を示すという傾向も同じであった．図 3.20 に示した
C1500 の摩擦過程を見ると，やはり筒状の形状のまま z 軸を回転軸とした回転を見せながらも，摩擦方向には回転していない．

第3章

単一フラーレンの解析

37

1nm+0.5D 押し込みの場合について，表 3.4，図 3.21，図 3.22 にこれまで同様にまとめて示す．図 3.21 の摩擦前の形状を見ると，いずれも筒状になったフラーレンの上下端面が相互作用するほど圧縮された状態であり，中でも C540 は前節の図 3.13 に示したように応力が著しく上昇を始めた後の形状である．こうした状態での摩擦係数を表3.4 で見ると，その値は各フラーレンともほぼ同じであり，これまでのような差はあまり見られない．図 3.22 に示した摩擦中のフラーレンの挙動も，摩擦方向には回転していない．

1nm+0.7D 押し込みした場合の結果について，表 3.5，図 3.23，図 3.24 に示す．図3.23 に示した摩擦前のフラーレンは，内部の隙間が Van der Waals 力によって斥力を生じる距離まで狭まっており，いずれも著しい応力上昇を示した後の形状である．摩擦係数は径の大きなフラーレンほど低い値を示しており (表 3.5 参照)，これまでとは違った傾向を示している．しかしながら，図 3.24 に示すように摩擦中のフラーレンの挙動はこれまでと同じく回転せずに移動しており，1nm+0.7D もの押し込みを行った
場合でも，結合の解消や切り替わりなどの変化は見られなかった．
以上のように，摩擦係数はいずれも 10−2 程度のオーダーを示したが，押し込み量によってその大小関係はさまざまである．押し込み量ごとに各フラーレンの摩擦力，圧子と作用する原子中の五員環原子の割合などを調べたが，その傾向は様々であり，摩擦係数の大小関係を説明できる要因は見つからなかった．次の 5 章では大きな凹凸のある圧子を用いた解析を行うが，そこではフラーレンの径の違いによる摩擦係数の傾
向がはっきりと見られている．本章ではフラットな面を持つダイヤモンド壁を用いたために，フラーレンの形状とダイヤモンド壁の表面の原子凹凸に依存した複雑な結果になったものと推測される．

第3章

単一フラーレンの解析

average

0.1

Friction coefficient, µFriction coefficient, µaverage

0

-0.1

0

5

0.1

0

-0.1

10

0

5

10

Displacement of upper wall, nm

Displacement of upper wall, nm(a) 1nm+0.1D indentation(b) 1nm+0.3D indentation

averageaverage

0.1

Friction coefficient, µFriction coefficient, µ

38

0

-0.1

0

5

10

0.1

0

-0.1

0

5

10

Displacement of upper wall, nm

Displacement of upper wall, nm(c) 1nm+0.5D indentation(d) 1nm+0.7D indentation

Fig.3.16 Friction coeﬃcient - displacement curves of C1500 under scratch.

第3章

単一フラーレンの解析

Table 3.2 Friction coeﬃcient of fullerenes (1nm+0.1D indentation).
fullerene

friction coeﬃcient

C540
C960
C1500
C2160

0.97 × 10−20.67 × 10−2
1.38 × 10−20.42 × 10−2zxy(a) C540(c) C1500(b) C960(d) C2160

Fig.3.17 Snapshots of fullerenes after indentation (1nm+0.1D indentation).

z

14.9nm(a) 0nm scratchx(b) D (3.61nm) scratchy(c) 2D (7.22nm) scratch(d) 3D (10.83nm) scratch

Fig.3.18 Snapshots of C1500 under scratch (1nm+0.1D indentation).

39

第3章

単一フラーレンの解析

Table 3.3 Friction coeﬃcient of fullerenes (1nm+0.3D indentation).
fullerene

friction coeﬃcient

C540
C960
C1500
C2160

1.45 × 10−20.63 × 10−2
1.75 × 10−20.75 × 10−2zxy(a) C540(c) C1500(b) C960(d) C2160

Fig.3.19 Snapshots of fullerenes after indentation (1nm+0.3D indentation).

z

14.9nm(a) 0nm scratchx(b) D (3.61nm) scratchy(c) 2D (7.22nm) scratch(d) 3D (10.83nm) scratch

Fig.3.20 Snapshots of C1500 under scratch (1nm+0.3D indentation).

40

第3章

単一フラーレンの解析

Table 3.4 Friction coeﬃcient of fullerenes (1nm+0.5D indentation).
fullerene

friction coeﬃcient

C540
C960
C1500
C2160

1.10 × 10−2
1.19 × 10−2
1.02 × 10−2
1.18 × 10−2zxy(a) C540(c) C1500(b) C960(d) C2160

Fig.3.21 Snapshots of fullerenes after indentation (1nm+0.5D indentation).

14.9nmz(a) 0nm scratchx(b) D (3.61nm) scratchy(c) 2D (7.22nm) scratch(d) 3D (10.83nm) scratch

Fig.3.22 Snapshots of C1500 under scratch (1nm+0.5D indentation).

41

第3章

単一フラーレンの解析

Table 3.5 Friction coeﬃcient of fullerenes (1nm+0.7D indentation).
fullerene

friction coeﬃcient

C540
C960
C1500
C2160

1.10 × 10−20.93 × 10−20.62 × 10−20.52 × 10−2zxy(a) C540(c) C1500(b) C960(d) C2160

Fig.3.23 Snapshots of fullerenes after indentation (1nm+0.7D indentation).

14.9nmz(a) 0nm scratchx(b) D (3.61nm) scratchy(c) 2D (7.22nm) scratch(d) 3D (10.83nm) scratch

Fig.3.24 Snapshots of C1500 under scratch (1nm+0.7D indentation).

42

第4章
単一 OLC の解析
4.1

初期構造緩和シミュレーション

4.1.1

解析条件

前章で説明したフラーレンを入れ子状にすることで作成した OLC を解析対象とし，
全方向自由境界条件のもとで 30000[fs] の緩和計算を行った．他の条件は前章と同じである．以後では，OLC について ＠C1500 のように表記する．ここで，＠ はフラーレン
でなく OLC であることを，また C1500 は最外殻のフラーレンを表している．今回は ＠

C540 ，＠C960 ，＠C1500 ，＠C2160 の 4 つを解析対象としている．

4.1.2

解析結果

例として，図 4.1 に ＠C1500 の緩和計算後の構造を示す．図 4.1 を見ると，フラーレンと同様，多少の凹凸はあるものの球殻に近い形状となっているのがわかる．また，その各層間の距離は，前章のフラーレン単体で比べたときの半径の差と変りなかった．他の OLC でも同様の傾向が見られた．

43

第4章

単一 OLC の解析

z

x

1

y

2

(a) top view(b) side view ( 1 )(c) side view ( 2 )Fig.4.1 Snapshots of ＠C1500 after relaxation.

44

第4章

4.2

圧縮シミュレーション

4.2.1

頂点保持による圧縮

単一 OLC の解析

45

解析条件
前章のフラーレンと同様，各 OLC を対象に，上下頂点の五員環をつかみ部として固定しながら圧縮するシミュレーションを行った．ひずみ制御で圧縮しており，ひずみ
速度は 5.0 × 10−6 [/fs] である．他の条件も前章と同じである．
解析結果
図 4.2 に ＠C540 の圧縮シミュレーション中の応力−ひずみ関係と変形の様子をまとめて示す．図 4.2(i) の ε < 0.05 の圧縮初期を見ると，(a)ε = 0.01 でわずかに応力上昇しているがすぐに低下し，(b)ε = 0.02 から (c)ε = 0.05 まではあまり上昇せずに進む．
このときの様子を図 4.2(ii) で確認すると，やはり (b) で下端のつかみ部である五員環が内側に向かってへこんでいる．その後，図 4.2(i) の応力は (c)ε = 0.05 から急激に上昇を続けるものの，(d)ε = 0.215 において一つの応力ピークを示す．この応力ピーク前後の構造を図 4.2(ii) で確認すると，(e) で最外殻の C540 の下部の五員環群 (前章の図3.7(a) に青丸で示した五員環群) が球の内側にへこんでおり，C540 の座屈応答が表れたものと考えられる．
図 4.3 に ＠C960 の圧縮シミュレーション中の応力とひずみの関係と変形の様子をまとめて示す．図 4.3(i) の応力−ひずみ関係を見ると，圧縮初期からあまり下降することなく上昇を続けている．圧縮初期の形状を図 4.3(ii) の (a) で確認すると，先の ＠C540のように最外殻のつかみ部のみへこむという様子は見られなかった．その後多少下降する様子が見られた (b)，(e) での構造を確認すると，(b) では内部の C540 が，(e) では
最外殻の C960 が座屈していた．なお座屈の判別は，前章のフラーレン単体での頂点保持による圧縮において応力ピークを示した後の形状，すなわち C540 は図 3.7(a) に青丸で示した五員環群が，C960 は図 3.7(b) に赤丸で示した五員環群が球の内側にへこんでいるかどうかで判断している．
図 4.4 に ＠C1500 の圧縮シミュレーションの結果を示す．図 4.4(i) の ε < 0.05 の圧縮

第4章

単一 OLC の解析

46

初期を見ると，応力はあまり上昇せずひずみのみ増加している．図 4.4(ii) の (b) で構造を確認すると，＠C540 と同様に，最外殻のつかみ部である五員環が内側に向かってへこんでいた．その後は応力は単調に上昇するが，(c)〜(e) で一時停滞し踊り場を示す．
図 4.4(ii) に示すように，(c) 点では内部の C540 が，(d) 点では最外殻の C1500 が，(e) 点
では内部の C960 がそれぞれ座屈していた．
図 4.5 に ＠C2160 の圧縮シミュレーションの結果を示す．応力の単調増加の傾向はさらに強まり，中実球の応答に近づいているものと推測される．また応力のゆらぎも小さくなる．変形の様子を図 4.5(ii) で確認すると，(b) で C540 ，(c) で C1500 ，(d) で C960 ，そして (e) で最外殻の C2160 が座屈するという段階的な変化が見られた．
図 4.6 に全ての OLC の応力−ひずみ関係をまとめて示す．また図 4.7 には比較のため，前章のフラーレンに対して行ったひずみ制御による圧縮の図 3.6 を再掲する．ただし，図 4.7 の応力のスケールは，図 4.6 のスケールに合わせているため図 3.6 のそれとは異なっていることに注意されたい．図より，OLC はフラーレンよりも高い抵抗力を示すことがわかる．前述のように，OLC は各層のフラーレンの座屈応答を受けながらも，全体的に大きく下降することなく上昇を続けている．これは多層構造を有することで，一つのフラーレンが座屈しても他のフラーレンが抵抗力を示すためであると考えられる．また図 4.5 の ＠C2160 において，各フラーレンが座屈する順序は C540 →C1500 → C960 → C2160 となっており，図 4.7 の応力ピークを示すひずみの大小関係と一致している．各フラーレンが座屈するひずみの値は，図 4.5 の C960 ，C1500 ，C2160 は図4.7 の応力ピークと同程度であるが，C540 はより早い段階で座屈しており (図 4.5 では
ひずみ 0.16，図 4.7 ではひずみ 0.19)，OLC の内部でより強い力を受けたものと推測さ
れる．

第4章

単一 OLC の解析

Compressive stress, -σ , GPa

@C540
(d)
4

(e) (f)
2

(a)
0

(b)
0

(c)
0.1

0.2

0.3

Compressive strain, | ε |(i) Stress - strain curve of ＠C540 under compression.

(a) ε =0.01(b) ε =0.02(c) ε =0.05zxy(d) ε =0.215(e) ε =0.225(f) ε=0.245(ii) Snapshots of ＠C540 under compression.
Fig.4.2 Compression of ＠C540 by holding the ﬁve-membered ring at the top
and bottom.

47

第4章

単一 OLC の解析

Compressive stress, -σ , GPa

@C960
4

(c)
(d)

(b)
2

0

(e)

(a)
0

0.1

0.2

0.3

Compressive strain, | ε |(i) Stress - strain curve of ＠C960 under compression.

(a) ε =0.065(b) ε =0.245

close-up C540zxy(c) ε =0.255(d) ε =0.265(e) ε =0.275(ii) Snapshots of ＠C960 under compression.
Fig.4.3 Compression of ＠C960 by holding the ﬁve-membered ring at the top
and bottom.

48

第4章

単一 OLC の解析

49

Compressive stress, -σ , GPa@C15004(c) (d)(e)
2

(a)
0

(b)
0

0.1

0.2

0.3

Compressive strain, | ε |(i) Stress - strain curve of ＠C1500 under compression.

(a) ε =0.025(b) ε =0.04(c) ε =0.2(d) ε =0.22(e) ε =0.245

close-up C540zxyclose-up C960

(ii) Snapshots of ＠C1500 under compression.
Fig.4.4 Compression of ＠C1500 by holding the ﬁve-membered ring at the top
and bottom.

第4章

単一 OLC の解析

Compressive stress, -σ , GPa@C21604(d) (e)2(b) (c)

(a)

0
0

0.1

0.2

0.3

Compressive strain, | ε |(i) Stress - strain curve of ＠C2160 under compression.

(a) ε =0.06(b) ε =0.16

close-up C540zxy(c) ε =0.185(d) ε =0.23

close-up C960

close-up C1500(e) ε =0.26(ii) Snapshots of ＠C2160 under compression.
Fig.4.5 Compression of ＠C2160 by holding the ﬁve-membered ring at the top
and bottom.

50

Compressive stress, -σ , GPa

第4章

単一 OLC の解析

@C540
@C960
@C1500@C2160

4

2

0
0

0.1

0.2

0.3

Compressive strain, | ε |Compressive stress, -σ , GPa

Fig.4.6 Stress - strain curves of OLCs under compression.

C540
C960
C1500
C2160

4

2

0
0

0.1

0.2

0.3

Compressive strain, | ε |Fig.4.7 Stress - strain curves of fullerenes under compression.

51

第4章

4.2.2

単一 OLC の解析

52

ダイヤモンド壁による圧縮

解析条件
各 OLC を対象に，剛体ダイヤモンド壁を用いて圧縮シミュレーションを行った．シ
ミュレーション条件はフラーレンの場合と同様であり，圧縮に用いるダイヤモンド壁のサイズは 7.1[nm]×7.1[nm]×0.7[nm]，その下降速度は 1.0 × 10−4 [nm/fs] である.
解析結果
圧縮シミュレーション中に生じる応力と変位の関係と変形の様子を図 4.8 に ＠C540 の場合を例として示す．図 4.8(i) の縦軸の応力は前節の図 4.6 の 2.5 倍程度の範囲になっていることに注意されたい．図 4.8(i) の応力−変位関係を見ると，Van der Waals 力の影響を受けて (a) の変位 0.7[nm] 近傍から応力上昇し，その後は (c) の変位 1.6[nm]で多少減少しながらも，最終的に (e) の変位 1.9[nm] 近傍から急激に立ち上がるという様子が確認できる．変形の様子を図 4.8(ii) で見ると，頂点保持では上部がへこみながら潰れていたのに対し (図 4.2(ii) 参照)，ダイヤモンド壁に沿うように変形し，へこむ
ことなく上下共平面となっている．応力が急激に立ち上がる直前の (c)〜(e) の挙動を細かく観察すると，(c) → (d) では内部のフラーレンの偏りの影響を受けて，最外殻のC540 の一部が内側に向かってへこんでおり，(e) では中心部の C60 が崩壊していた．
図 4.9 に全ての OLC の応力−変位関係をまとめて示す．変位 1.0[nm] 以下の圧縮初期を見ると，高次の OLC も ＠C540 同様に，壁面と Van der Waals 作用する変位 0.7[nm]近傍で応力上昇している．＠C1500 のみ立ち上がりが少し遅れているが，このときの構造を確認すると，前節の図 4.4(ii) の (b) 同様に C1500 の五員環のみ内側に向かってへこんでおり，最外殻の応答が表れたものと考えられる．応力上昇を始めてから，最外層のフラーレンの応答を受けながらも大きく下降することなく上昇を続け，最終的に急
激な立ち上がりを見せる，という全体的な傾向はいずれも同じである．そうした急激な立ち上がりを見せる点では，いずれの OLC も中心部の C60 が崩壊している様子が確認できた．また図 4.9 の応力−変位関係は，前節の図 4.6 に比べて滑らかとなっている
が，これは先の図 4.2(ii) に示したように，上下共平面となって変形するために，OLCを形作る各フラーレンが座屈せずに圧縮されたためであると考えられる．径の大きな第4章

単一 OLC の解析

53

＠C1500 ，＠C2160 は応力が下降することもなく，特に滑らかとなっているが，これは層数が多いことで，内部の偏りの影響が最外殻に表れにくかったことが原因であると考えられる．
比較のため，前章のフラーレンに対して行ったダイヤモンド平板による圧縮応答の
図 3.13 を図 4.10 に再掲する．ここで，図 4.10 の応力のスケールは，図 4.9 のスケールに合わせており，図 3.13 のそれとは異なっている．図を見比べると，図 4.9 の OLC の応力は，最外層しか作用しないような圧縮初期 (変位 0.9[nm] 以前) は同径のフラーレンと大差ないが，それ以降は内部のフラーレンが抵抗となることで，早い段階から格段に高い応力を示している．このように同じ変位においては，フラーレンより OLC の方が，また径の大きなものより小さなものの方が高い応力を示した．

第4章

単一 OLC の解析

Compressive stress, -σ , GPa

@C540
10

(c)
(e)
(d)

5

(b)
0

(a)
0

2

4

Displacement of upper wall, nm(i) Stress - displacement curve of ＠C540 under compression.

(a) 0.7nm indentation(b) 0.85nm indentation(c) 1.6nm indentationzxy(d) 1.7nm indentation(e) 1.9nm indentation

close-up C60(ii) Snapshots of ＠C540 under compression.
Fig.4.8 Compression of ＠C540 by diamond wall.

54

Compressive stress, -σ , GPa

第4章

単一 OLC の解析

@C540
@C960
@C1500@C2160

10

5

0
0

2

4

Displacement of upper wall, nm

Compressive stress, -σ , GPa

Fig.4.9 Stress - displacement curves of OLCs under compression.

C540
C960
C1500
C2160

10

5

0
0

2

4

Displacement of upper wall, nm
Fig.4.10 Stress - displacement curves of fullerenes under compression.

55

第4章

単一 OLC の解析

4.3

ダイヤモンド壁による摩擦シミュレーション

4.3.1

解析条件

56

各 OLC を対象に，剛体ダイヤモンド壁を用いて摩擦シミュレーションを行った．今回用いたダイヤモンド壁のサイズは OLC の直径の 4 倍程度であり，そのサイズは ＠

C540 で 10.7[nm]，＠C960 で 12.8[nm]，＠C1500 で 14.9[nm]，＠C2160 で 17.1[nm] とした．
押し込み量は 1nm+0.1D，1nm+0.3D，1nm+0.5D，1nm+0.7D の 4 通りであり，摩擦量や壁の移動速度など，その他の条件も前節フラーレンと同様としている．

4.3.2

解析結果

1nm+0.1D 押し込みしたときの摩擦シミュレーション中の摩擦係数と変位の関係を，
＠C1500 の場合を例として図 4.11 に示す．図中の赤い線は平均値を表しており，右隣の
表 4.1 には各 OLC のダイヤモンド壁との接触面積を円で近似した場合の半径と摩擦係数の平均値をまとめている．図 4.11 を見ると，前章のフラーレン同様，マイナス方向にもゆらいではいるが，平均値としてはプラスとなっている．またその縦軸のスケー
ルはフラーレンの図 3.16(a) の 5 分の 1 程度の範囲となっており，同じ押し込み量でもOLC のゆらぎはフラーレンより小さくなっている．表 4.1 の摩擦係数を見ると，いずれもフラーレン同様 10−2 程度のスケールとなっており，その値は径の大きな OLC の方が小さくなっている．図 4.12 に摩擦前の押し込み状態での形状を示す．各 OLC とも内部に大きな偏りもなく，最外殻は壁面に沿って平面となっている．＠C1500 の摩擦中の変形の様子を図 4.13 に示す．摩擦中も圧縮された状態で，やはり回転することなく滑るように移動しているのがわかる．壁面と接触する平面の面積は，径の大きな OLC
ほど大きくなっており，摩擦係数は接触面積が大きいほど小さくなっている．
＠C1500 の 1nm+0.3D 押し込み摩擦シミュレーションにより得られた摩擦係数−変位関係を例として図 4.14 に示す．また各 OLC の接触面積の半径，摩擦係数の平均値を表4.2 に，摩擦前の押し込み状態での形状を図 4.15 に示している．図 4.14 の摩擦係数−
変位関係において，マイナス方向にもゆらぎながら推移する傾向は先の 1nm+0.1D の第4章

単一 OLC の解析

57

場合と変わりないが，そのゆらぎの幅は ±0.01 程度と小さくなっている．図 4.15 の摩擦前の形状を見ると，いずれも内部の偏りの影響を受けて一部内側にへこんでいるものの，外殻は壁面に沿った平面となり円筒状となっている．表 4.2 の摩擦係数も，やはり径の大きな OLC ほど低い値を示しているが，いずれも値は先の 1nm+0.1D のそれより小さくなっている．また径が小さい方がその減少割合は大きい．図 4.16 に ＠C1500の摩擦中の様子を例として示す．フラーレンの場合と同じくダイヤモンド壁に垂直なz 軸を回転軸とした回転を生じているが，やはり摩擦方向に転がるような挙動は生じていない．

1nm+0.5D 押し込みの結果についてもこれまで同様に図 4.17〜図 4.19 と表 4.3 にまとめた．図 4.17 の ＠C1500 の摩擦係数−変位曲線を見ると，ゆらぎの幅は先の 1nm+0.3D
押し込みと大差ないが不均一であり，これまでよりマイナス方向にゆらぐことは少な
くなり全体的に高い数値をとっている．また，摩擦初期は (b) → (c) → (d) のように変位 3.5[nm] 近傍で摩擦係数が低下している．図 4.18 の摩擦前の構造は，いずれの OLC
も中心部の C60 が崩壊した形状である．径が大きい OLC ほど摩擦係数が小さくなる傾向は同じであるが (表 4.3)，先の 1nm+0.1D，1nm+0.3D に比べて大きな値を示している．＠C1500 の摩擦過程を図 4.19 で確認すると，これまでとは異なり OLC が摩擦方向に転がる様子が確認できた．なお図 4.19 の (b)〜(d) のスナップショットは図 4.17 中に矢印で示した各点にそれぞれ対応している．摩擦初期の (a) → (b) では内部の C60 ，

C240 が後方に偏り，外側の層が前方に移動して回転を始めている．摩擦係数が低下する (b) → (c) では中心の偏りのため外側の層が変形しやすくなっており，(c) → (d) で下側の凸部が引き込まれるような回転をしている．図 4.17 の (c) で一時的に低い摩擦係数を示したのはこのような転がりによるものと考えられる．その後は (e) → (f) のように転がらずに滑り，摩擦係数も図 4.17 に示すように平均を中心に振動する．1nm+0.5Dの押し込み摩擦では，＠C960 でも転がるような変形が見られたが，＠C540 ，＠C2160 ではそうした変形は見られず，これまで同様滑るように移動している．＠C960 ，＠C1500
で回転を生じたのは，図 4.18 に示したように中心の高密度部分に比べて外側の層に疎な部分が存在したためと考えられる．

1nm+0.7D 押し込みした場合の結果について，図 4.20〜図 4.22 および表 4.4 に示す．

第4章

単一 OLC の解析

58

図 4.20 の ＠C1500 の摩擦係数の変化を見ると，これまでは見られなかった静摩擦に対応するピークが摩擦初期に表れている．それ以外はこれまで同様ゆらぎながらほぼ一定値を示している．表 4.4 の平均より得た摩擦係数は，これまでと同じく大きな OLC
ほど低い摩擦係数を示す．その値は先の 1nm+0.5D より低く，押し込み量による変化
はいずれの OLC も 1nm+0.5D で最も高い摩擦係数を示している．摩擦前の形状を図4.21 で確認すると，OLC はかろうじて最外層が判別できる程度で，それ以外の層はクラスター状になっている．＠C1500 の摩擦過程の図 4.22 を見ると，1nm+0.5D 押し込
み時のような転がり変形はなく，また移動量は他の押し込み量に比べて小さい．なお，
構造に大きな変化は見られなかったものの，摩擦初期には共有結合の bond 数が減少しており，初期の高い摩擦係数はボンドの切断に起因するものと考えられる．

第4章

単一 OLC の解析

59

Friction coefficient, µaverage
0.02

Table 4.1 Radius of contact area and friction coeﬃcient of OLCs (1nm+0.1D indentation).

0

-0.02

0

5

10

Displacement of upper wall, nm

OLC

Rc [nm]friction coeﬃcient

＠C540
＠C960
＠C1500＠C2160

0.78
1.05
1.25
1.47

1.49 × 10−20.54 × 10−20.47 × 10−20.34 × 10−2

Fig.4.11 Friction coeﬃcient - displacement curveof ＠C1500 under scratch (1nm+0.1D indentation).

zxy(a) @C540(c) @C1500(b) @C960(d) @C2160

Fig.4.12 Snapshots of OLCs after indentation (1nm+0.1D indentation).

z

14.9nm(a) 0nm scratchx(b) D (3.61nm) scratchy(c) 2D (7.22nm) scratch(d) 3D (10.83nm) scratch

Fig.4.13 Snapshots of ＠C1500 under scratch (1nm+0.1D indentation).

第4章

単一 OLC の解析

60

Friction coefficient, µaverage
0.02

Table 4.2 Radius of contact area and friction coeﬃcient of OLCs (1nm+0.3D indentation).

0

-0.02

0

5

10

Displacement of upper wall, nm

OLC

Rc [nm]friction coeﬃcient

＠C540
＠C960
＠C1500＠C2160

0.99
1.33
1.62
1.92

0.82 × 10−20.47 × 10−20.43 × 10−20.30 × 10−2

Fig.4.14 Friction coeﬃcient - displacement curveof ＠C1500 under scratch (1nm+0.3D indentation).

zxy(a) @C540(c) @C1500(b) @C960(d) @C2160

Fig.4.15 Snapshots of OLCs after indentation (1nm+0.3D indentation).

z

14.9nm(a) 0nm scratchx(b) D (3.61nm) scratchy(c) 2D (7.22nm) scratch(d) 3D (10.83nm) scratch

Fig.4.16 Snapshots of ＠C1500 under scratch (1nm+0.3D indentation).

第4章

単一 OLC の解析

61

average
Friction coefficient, µ

(d)
(b)

0.02

Table 4.3 Radius of contact area and friction coeﬃcient of OLCs (1nm+0.5D indentation).

0

(c)
-0.02

0

5

10

Displacement of upper wall, nm

OLC

Rc [nm]friction coeﬃcient

＠C540
＠C960
＠C1500＠C2160

1.17
1.60
1.92
2.26

1.66 × 10−2
1.45 × 10−2
1.22 × 10−2
1.20 × 10−2

Fig.4.17 Friction coeﬃcient - displacement curveof ＠C1500 under scratch (1nm+0.5D indentation).

zx(a) @C540y(c) @C1500(b) @C960(d) @C2160

Fig.4.18 Snapshots of OLCs after indentation (1nm+0.5D indentation).

14.9nm(b) 2.8nm scratch(a) 0nm scratch(c) 3.5nm scratch

surface closeup(d) 5.0nm scratchzx(e) 7.0nm scratchy(f) 3D (10.83nm) scratch

Fig.4.19 Snapshots of ＠C1500 under scratch (1nm+0.5D indentation).

第4章

単一 OLC の解析

62

Friction coefficient, µaverage
0.02

Table 4.4 Radius of contact area and friction coeﬃcient of OLCs (1nm+0.7D indentation).

0

-0.02

0

5

10

Displacement of upper wall, nm

OLC

Rc [nm]friction coeﬃcient

＠C540
＠C960
＠C1500＠C2160

1.43
1.98
2.46
2.92

1.24 × 10−2
1.05 × 10−20.55 × 10−2
0.48 × 10−2

Fig.4.20 Friction coeﬃcient - displacement curveof ＠C1500 under scratch (1nm+0.7D indentation).

zxy(a) @C540(c) @C1500(b) @C960(d) @C2160

Fig.4.21 Snapshots of OLCs after indentation (1nm+0.7D indentation).

14.9nmz(a) 0nm scratchx(b) D (3.61nm) scratchy(c) 2D (7.22nm) scratch(d) 3D (10.83nm) scratch

Fig.4.22 Snapshots of ＠C1500 under scratch (1nm+0.7D indentation).

第5章
薄膜構造での摩擦シミュレーション
5.1

圧子の表面起伏が摩擦係数に与える影響

5.1.1

解析条件

フラーレン・OLC を薄膜状に配置したモデルを対象に，表面に起伏を有する剛体
ダイヤモンド壁を用いて摩擦シミュレーションを行った．ここでフラーレンは C540 とC2160 を，OLC は ＠C540 と ＠C2160 を対象としている．ダイヤモンド壁の表面起伏は，
図 5.1 に示すようにのこ歯状の切れ込みとし，切れ込みの大きさによって small serrate
圧子 (S) と large serrate 圧子 (L) の 2 通りとした．まず図 5.2 のように平面上に配置した各モデル (C540 と ＠C540 は 2.5[nm] 間隔，C2160 と ＠C2160 は 5.0[nm] 間隔) について，これまで同様，z 軸方向下端に平滑なダイヤモンド壁を配置して 30000[fs] の緩和を行った．その後，押し込みをしてから摩擦を行うが，押し込み量は圧子が受ける z
軸方向の反力の和が 1000[nN] に達するまで行うこととし，水平方向の摩擦量はいずれも 7.5[nm] とした．壁の下降・スライド速度はともに 1.0 × 10−4 [nm/fs] である.
各モデルについて，フラーレン・OLC 個体の直径，ダイヤモンド壁以外の原子数，
そして底面のダイヤモンド壁面積 (14.9×14.9[nm2 ]) をもとに算出した被覆率を表 5.1にまとめて示す．原子数は異なるが，被覆率はいずれも同程度としている．

63

第5章

薄膜構造での摩擦シミュレーション

64

Table 5.1 Diameter, number of atoms and coverage rate in simulation models.
model namediameter [nm]No. of atoms

coverage rate [%]

C540
C2160
＠C540
＠C2160

2.16
4.34
2.16
4.34

19440
19440
30240
49140

59.5
60.0
59.5
60.0

z
y

x
14.9nm14.9nm0.55nm
1.06nm(a) Small serrate indenter
1.06nm14.9nm14.9nm0.55nm1.24nm(b) Large serrate indenter
2.49nm
Fig.5.1 Dimensions of diamond wall indenter.

第5章

z

薄膜構造での摩擦シミュレーションxy(a) C540(b) C2160

2.5nm

5.0nm

(c) @C540(d) @C2160

Fig.5.2 Top view of fullerene/OLC array.

65

第5章

5.1.2

薄膜構造での摩擦シミュレーション

66

解析結果

small serrate 圧子を用いて C540 ，＠C540 に摩擦を行った時の摩擦係数−変位関係を
図 5.3 に示す．フラーレン・OLC 共に摩擦初期に最大静摩擦に相当するピークを示し，その後はゆらぎながらほぼ一定値を示している．3 章，4 章のダイヤモンド平板による摩擦の図 3.16 や図 4.11 のようにマイナス方向の摩擦係数を示すことはほとんどない．
これは C2160 ，＠C2160 についても同様である．初期のピーク後に再び上昇に転じる位
置からの摩擦係数の平均値を，圧子反力が 1000[nN] のときの押し込み量と共に表 5.2に示す．押し込み量はフラーレンより OLC の方が，径の大きなモデルより小さなモデ
ルの方が低い値を示している．これは 3 章，4 章で示してきたように，フラーレンより

OLC の方が，径の大きなモデルより小さなモデルの方が早い段階で高い反力を示すこ. ＠C の順で高く，小さな C ，＠とによる．摩擦係数は ＠C2160 < C2160 < C540 =
.
540
540
C540 の方が，大きな C2160 ，＠C2160 より高い値を示す．
隣接するフラーレン・OLC の中で Van der Waals 相互作用 (VDW) する結合数変化を図 5.4 に示す．フラーレン・OLC 単体中での VDW や圧子や基板との VDW は含まれていないため注意願いたい．隣のフラーレン・OLC と相互作用する数は，C540 ，＠

C540 の方が C2160 ，＠C2160 より格段に多くなっており，高い摩擦係数を示した要因であると考えられる．図 5.5 に各モデルの圧縮後の形状を示す．横からの図のため全てのフラーレン・OLC が合体しているように見えるが，後に示すように上から見た図 (図5.6，図 5.7) では個々のフラーレン・OLC は独立している．(c)，(d) の OLC は内部のフラーレンが抵抗となるため，圧子表面の切れ込みには最外層が少し入り込む程度であるが，(a)，(b) のフラーレンは圧子形状に合わせて自由に変形しており，特に C2160
はその押し込み量の多さも合わさって切れ込み深くまで入り込んでいる．このようにしっかりと圧子と作用していたことで，C540 は ＠C540 と同程度，C2160 は ＠C2160 よりも高い摩擦係数を示したものと推測される．
図 5.6，図 5.7 に C540 と ＠C540 の摩擦中のスナップショットを上から見たものを示す．図中ピンクで着色したマーカー原子からわかるように，図 5.7 の ＠C540 は多少回転しているものの，その回転は摩擦量 (7.5[nm]) からすると微小であり，他のモデルは

第5章薄膜構造での摩擦シミュレーション

67

図 5.6 の C540 と同様に押し込み後の形状のまま回転することなく，切れ込みに挟まった状態で滑るように移動していた．

large serrate 圧子を用いて C540 ，＠C540 に摩擦を行った時の摩擦係数と変位の関係を図 5.8 に示す．縦軸の摩擦係数は small serrate 圧子を用いた場合の図 5.3 の 2 倍程度のスケールとなっている．初期に最大静摩擦を示すという傾向は同様であるが，その値は small serrate 圧子の場合よりも高くなっている．こうした変化は C2160 ，＠C2160
でも見られた．先ほどと同様に最大静摩擦を除いて平均した摩擦係数を表 5.3 にまとめて示した．大小関係は表 5.2 と同様であるが，摩擦係数は小さい C540 ，＠C540 では増加し，大きな C2160 ，＠C2160 では減少している．
図 5.9 に large serrate 圧子での摩擦シミュレーション中の隣接フラーレン・OLC 間の VDW 相互作用の変化を示す．先ほどの図 5.4 同様，C540 ，＠C540 の方が，C2160 ，
＠C2160 より VDW 作用する原子の数が多い．しかし，先の図 5.4 と比べると，小さな凹凸の摩擦よりフラーレンの C540 ，C2160 では増加し，OLC の ＠C540 ，＠C2160 では減
少しており，フラーレンと OLC 間の大小関係が逆転している．図 5.10 に各モデルの圧縮後の形状を横から見たものを示す．切れ込みが大きくなったことで，フラーレン・

OLC 共に深くまで入り込んでいるのがわかる．図 5.10 の (c)，(d) に示す OLC は，初期配置と切れ込みの関係から，一定の間隔を保ってきれいに整列している．このため押し込み量が増えているにもかかわらず，OLC では VDW の数が減少したものと考えられる．また図 5.10 の (a)，(c) を比較すると，C540 は圧子と接する部分がへこんでしまっているが，＠C540 は内部が抵抗となることで，圧子としっかりと面で作用している．このようにしっかりと作用していたことで，＠C540 は原子間非共有結合の数が減
少しても small serrate 圧子の場合同様，C540 より高い摩擦係数を示したものと考えられる．図 5.10 と図 5.5 を比較すると，それぞれの (a)，(c) に示す C540 ，＠C540 は，径が小さいためどちらの場合でも一つの切れ込みとしか作用していないが，(b)，(d) に示す C2160 と ＠C2160 は，small serrate 圧子の場合はその径の大きさから複数の切れ込みと作用していたのに対し，large serrate 圧子の場合は一つの切れ込みとしか作用していない．このように C2160 と ＠C2160 は切れ込み深くまで入り込んでいるものの，C540 ，
＠C540 ほど圧子と作用する原子が増えたわけではないと推測される．ここで図 5.1 に

第5章

薄膜構造での摩擦シミュレーション68

あるように，large serrate 圧子の切れ込みの傾きは，small serrate 圧子の傾きより緩い
ため，VDW により圧子に与える力の x 軸方向成分はより小さくなるものと考えられる．このような切れ込みの傾きによる影響以上に圧子と作用する原子が増加していた
ことが，表 5.3 における C540 ，＠C540 の摩擦係数が表 5.2 より増加していた要因であると考えられる．
図 5.11，図 5.12 に摩擦シミュレーション中の C540 と ＠C540 を上から見たスナップ
ショットを示す．こちらも small serrate 圧子の場合同様，挟まった状態で滑るように移動しており，切れ込みは大きくなったものの，その中で回転するといった様子も見られなかった．

第5章

薄膜構造での摩擦シミュレーション

average=2.32x10-2
Friction coefficient, µFriction coefficient, µaverage=2.29x10-2
0.1

0.05

69

0.1

0.05

0

0
0

3

6

0

3

6

Displacement of indenter, nm

Displacement of indenter, nm(a) C540(b) @C540

Fig.5.3 Friction coeﬃcient - displacement curves under scratch (small serrate
indenter).

Table 5.2 Indentation depth and friction coeﬃcient (small serrate indenter).

model nameindentation depth [nm]friction coeﬃcient

C540 , small serrate indenter (C540 -S)
C2160 , small serrate indenter (C2160 -S)
＠C540 , small serrate indenter (＠C540 -S)
＠C2160 , small serrate indenter (＠C2160 -S)

2.07
4.18
1.26
1.86

2.29 × 10−2
1.73 × 10−2
2.32 × 10−2
1.18 × 10−2

第5章

薄膜構造での摩擦シミュレーション

C540-S
C2160-S@C540-S@C2160-S(x105)

Number of VDW bonds

2.0

1.5

1.0

0.5
0

3

6

Displacement of indenter, nm
Fig.5.4 Number of VDW bonds - displacement curves under scratch (small
serrate indenter).

(a) C540-S(b) C2160-Szxelement closeup

element closeupy(c) @C540-S(d) @C2160-S

Fig.5.5 Snapshots after indentation (small serrate indenter).

70

第5章

z

薄膜構造での摩擦シミュレーションxy(a) 0nm scratch(b) 2.5nm scratch(c) 5.0nm scratch(d) 7.5nm scratch

Fig.5.6 Snapshots of C540 -S under scratch.

z

x

y

(a) 0nm scratch(b) 2.5nm scratch(c) 5.0nm scratch(d) 7.5nm scratch

Fig.5.7 Snapshots of ＠C540 -S under scratch.

71

第5章

薄膜構造での摩擦シミュレーション

average=2.71x10-2

average=2.56x10-2
0.2

Friction coefficient, µ0.2

Friction coefficient, µ

72

0.1

0.1

0

0
0

3

6

0

3

6

Displacement of indenter, nm

Displacement of indenter, nm(a) C540(b) @C540

Fig.5.8 Friction coeﬃcient - displacement curves under scratch (large serrate
indenter).

Table 5.3 Indentation depth and friction coeﬃcient (large serrate indenter).

model nameindentation depth [nm]friction coeﬃcient

C540 , large serrate indenter (C540 -L)
C2160 , large serrate indenter (C2160 -L)
＠C540 , large serrate indenter (＠C540 -L)
＠C2160 , large serrate indenter (＠C2160 -L)

2.48
4.45
1.82
2.16

2.56 × 10−2
1.46 × 10−2
2.71 × 10−2
1.07 × 10−2

第5章

薄膜構造での摩擦シミュレーション

C540-L
C2160-L@C540-L@C2160-L(x105)

Number of VDW bonds

2.0

1.5

1.0

0.5
0

3

6

Displacement of indenter, nm
Fig.5.9 Number of VDW bonds - displacement curves under scratch (large
serrate indenter).

(a) C540-L(b) C2160-Lzxelement closeup

element closeupy(c) @C540-L(d) @C2160-L

Fig.5.10 Snapshots after indentation (large serrate indenter).

73

第5章

z

薄膜構造での摩擦シミュレーションxy(a) 0nm scratch(b) 2.5nm scratch(c) 5.0nm scratch(d) 7.5nm scratch

Fig.5.11 Snapshots of C540 -L under scratch.

z

x

y

(a) 0nm scratch(b) 2.5nm scratch(c) 5.0nm scratch(d) 7.5nm scratch

Fig.5.12 Snapshots of ＠C540 -L under scratch.

74

第5章

薄膜構造での摩擦シミュレーション

5.2

基板の表面起伏が摩擦係数に与える影響

5.2.1

解析条件

75

前節と同様，フラーレン・OLC を薄膜状に配置したモデルを対象に，表面に起伏を有する剛体ダイヤモンド壁を用いて摩擦シミュレーションを行った．対象とするモデ
ルや圧子形状は先ほどと同様であるが，今回は図 5.13 にあるように，下端の基板にものこ歯状の切れ込みがあるものを用いることとした．基板は先の small serrate 圧子と同様の形状としている．基板を配置してからの初期緩和は，フラーレンで 80000[fs]，

OLC で 30000[fs] と先ほどより多くとったが，その他の条件はいずれも同様である．

1nm

z
x
y
14.9nm
Fig.5.13 Schematic of simulation model (＠C540 , small serrate substrate, small
serrate indenter).

5.2.2

解析結果

small serrate 圧子を用いた場合の摩擦シミュレーション中の摩擦係数−変位関係を
図 5.14 に示す．図 5.14 を見ると，(d) の ＠C2160 は小さなゆらぎを生じつつもほぼ一定値を示しているが，(a)〜(c) の他のモデルは大きなのこ歯状の波形を描いている．赤線で示した平均値は初期からの摩擦係数全体に対するものである．この波形の間隔は約1.0[nm] であり，small serrate 圧子の切れ込み幅と対応する．摩擦係数の平均値を，表5.4 に前節と同様にまとめる．押し込み量は前節の表 5.2 よりいずれも高い値を示して

第5章

薄膜構造での摩擦シミュレーション

76

いるものの，その大小関係は先ほどと同様である．また摩擦係数も高い値を示し，その増加量は C540 ，C2160 で 8 倍，＠C540 で 4 倍，＠C2160 で 3 倍程度とフラーレンの方が大きく，大小関係は先ほどとは異なり，＠C2160 < ＠C540 < C2160 < C540 の順で高くなっていた．
図 5.15 に横から見た各モデルの圧縮後の形状を示す．いずれのモデルも圧子だけでなく基板の切れ込みにも入り込んでおり，その入り込みは中空の形状であるフラーレ
ンの方が深い．このようにフラーレンは基板の切れ込みにも深く入り込みやすいために，摩擦時により大きな抵抗力を示したものと考えられる．また ＠C2160 はフラーレンのように中空ではなく，径の大きさに対して凹凸のサイズが小さかったことで，図5.14(d) で他のモデルのようなゆらぎ波形を描かなかったものと推測される．摩擦シ
ミュレーション中の変形の様子を図 5.16〜図 5.19 に示す．それぞれの図の上段は対象全体を z 軸方向上部から見たスナップショットであり，下段はそのうちの一つを y 軸方
向から見たスナップショットである．図 5.16〜図 5.19 の上段のスナップショットを見ると，いずれも基板の切れ込みの抵抗により，前節の図 5.6，図 5.7 よりも移動量が減
少している．また下段のスナップショットを見ると，上下の切れ込みによる挟み込みによって，いずれのモデルも摩擦中に回転する様子が確認できた．

large serrate 圧子を用いて摩擦を行った時の摩擦係数と変位の関係を図 5.20 に示す．
縦軸の摩擦係数の範囲は図 5.14 と多少異なっているので注意されたい．図 5.20(b) のC2160 は，先の図 5.14(a)〜(c) 同様にのこ歯状の波形を描いているが，その波形の間隔は一様ではない．(a) の C540 は変位 3[nm] 以下では平均より高い値をとるものの，それ以降は (c)，(d) の OLC と同様にある程度一定の値をとっている．表 5.5 にまとめた
押し込み量，摩擦係数を前節の表 5.3 と比較すると，押し込み量，摩擦係数共に増加しており，特に摩擦係数の増加量は C540 で 13 倍，C2160 で 18 倍，＠C540 で 8 倍，＠

C2160 で 9 倍程度となっている．また表 5.5 の押し込み量，摩擦係数の大小関係は表 5.4と一致しているが，その値は表 5.5 の方がいずれも大きかった．
図 5.21 の圧縮後の形状を見ると，いずれも基板より切れ込みの大きな圧子側により深く入り込んでいることがわかる．このように large serrate 圧子はその凹凸の大きさから small serrate 圧子よりも多くの原子と作用しており，基板の影響が表れやすかった

第5章

薄膜構造での摩擦シミュレーション

77

ことが，表 5.4 より表 5.5 の摩擦係数の方が大きかった原因であると推測される．large

serrate 圧子を用いた場合の摩擦過程について，図 5.22〜図 5.25 に示す．上段のスナッ
プショットを見ると，図 5.23 の C2160 は small serrate 圧子を用いた場合と同様に移動量が減少しているが，図 5.22 の C540 ，図 5.24 の ＠C540 ，図 5.25 の ＠C2160 にはそのよ
うな移動量の減少は見られない．C540 ，＠C540 ，＠C2160 の移動量が減少していなかったのは，図 5.21 に示したように，それぞれ切れ込みの大きな圧子側にしっかりと入り込んでいたためであり，C2160 の移動量が減少していたのは，径の大きなフラーレンであったことで，基板側の切れ込みにも複数にわたって深く入り込んでいたためである．
また下段のスナップショットを見ると，図 5.22〜図 5.24 の C540 ，C2160 ，＠C540 は先ほど同様に回転しているが，図 5.25 の ＠C2160 は他に比べてあまり回転していない．これは ＠C2160 が圧子の切れ込みに深く入り込んでいたこと，またその径の大きさに対して基板の切れ込みが小さかったことの２つが要因となり，基板側の外層のみ波打ちながら応答したためであると考えられる．

第5章

薄膜構造での摩擦シミュレーション

average=1.43x10-1

average=1.88x10-1
0.3

Friction coefficient, µFriction coefficient, µ

0.3

0.2

0.1

0.2

0.1

0

0
0

3

0

6

3

6

Displacement of indenter, nm

Displacement of indenter, nm(a) C540(b) C2160

average=0.90x10-1

average=0.33x10-1

0.3

0.3

Friction coefficient, µFriction coefficient, µ

78

0.2

0.1

0

0.2

0.1

0
0

3

6

0

3

6

Displacement of indenter, nm

Displacement of indenter, nm(c) @C540(d) @C2160

Fig.5.14 Friction coeﬃcient - displacement curves under scratch (small serrate
substrate, small serrate indenter).

第5章

薄膜構造での摩擦シミュレーション

79

Table 5.4 Indentation depth and friction coeﬃcient (small serrate substrate,
small serrate indenter).

model nameindentation depth [nm]friction coeﬃcient

C540 , small serrate substrate, small serrate indenter (C540 -SS)
C2160 , small serrate substrate, small serrate indenter (C2160 -SS)
＠C540 , small serrate substrate, small serrate indenter (＠C540 -SS)
＠C2160 , small serrate substrate, small serrate indenter (＠C2160 -SS)

2.23
4.35
1.47
2.01

1.88 × 10−1
1.43 × 10−10.90 × 10−1
0.33 × 10−1(a) C540-SS(b) C2160-SSzxelement closeup

element closeupy(c) @C540-SS(d) @C2160-SS

Fig.5.15 Snapshots after indentation (small serrate substrate, small serrate
indenter).

第5章

薄膜構造での摩擦シミュレーション

80

x

z

y

(a) 0nm scratch(b) 2.5nm scratch(c) 5.0nm scratch(d) 7.5nm scratchzxycloseupcloseupcloseupcloseupFig.5.16 Snapshots of C540 -SS under scratch.
xzy(a) 0nm scratch(b) 2.5nm scratch(c) 5.0nm scratch(d) 7.5nm scratchzxycloseupcloseupcloseupFig.5.17 Snapshots of C2160 -SS under scratch.

closeup

第5章

薄膜構造での摩擦シミュレーション

81

x

z

y

(a) 0nm scratch(b) 2.5nm scratch(c) 5.0nm scratch(d) 7.5nm scratchzxycloseupcloseupcloseupcloseupFig.5.18 Snapshots of ＠C540 -SS under scratch.
xzy(a) 0nm scratch(b) 2.5nm scratch(c) 5.0nm scratch(d) 7.5nm scratch

closeupcloseupcloseupzxycloseupFig.5.19 Snapshots of ＠C2160 -SS under scratch.

第5章

薄膜構造での摩擦シミュレーション

average=2.61x10-1

average=3.29x10-1
0.4

Friction coefficient, µFriction coefficient, µ

0.4

0.2

0

0

3

0

3

6

Displacement of indenter, nm

Displacement of indenter, nm(a) C540(b) C2160

average=2.10x10-1

average=1.00x10-1
0.4

Friction coefficient, µFriction coefficient, µ

0.2

0

6

0.4

0.2

0

82

0

3

6

0.2

0

0

3

6

Displacement of indenter, nm

Displacement of indenter, nm(c) @C540(d) @C2160

Fig.5.20 Friction coeﬃcient - displacement curves under scratch (small serrate
substrate, large serrate indenter).

第5章

薄膜構造での摩擦シミュレーション

83

Table 5.5 Indentation depth and friction coeﬃcient (small serrate substrate,
large serrate indenter).

model nameindentation depth [nm]friction coeﬃcient

C540 , small serrate substrate, large serrate indenter (C540 -SL)
C2160 , small serrate substrate, large serrate indenter (C2160 -SL)
＠C540 , small serrate substrate, large serrate indenter (＠C540 -SL)＠C2160 , small serrate substrate, large serrate indenter (＠C2160 -SL)

2.68
4.56
2.04
2.34

3.29 × 10−1
2.61 × 10−1
2.10 × 10−1
1.00 × 10−1(a) C540-SL(b) C2160-SLzxelement closeup

element closeupy(c) @C540-SL(d) @C2160-SL

Fig.5.21 Snapshots after indentation (small serrate substrate, large serrate indenter).

第5章

薄膜構造での摩擦シミュレーション

84

x

z

y

(a) 0nm scratch(b) 2.5nm scratch(c) 5.0nm scratch(d) 7.5nm scratchzxycloseupcloseupcloseupcloseupFig.5.22 Snapshots of C540 -SL under scratch.
xzy(a) 0nm scratch(b) 2.5nm scratch(c) 5.0nm scratch(d) 7.5nm scratchzxycloseupcloseupcloseupFig.5.23 Snapshots of C2160 -SL under scratch.

closeup

第5章

薄膜構造での摩擦シミュレーション

85

x

z

y

(a) 0nm scratch(b) 2.5nm scratch(c) 5.0nm scratch(d) 7.5nm scratch

closeupcloseupzxycloseupcloseupFig.5.24 Snapshots of ＠C540 -SL under scratch.
xzy(a) 0nm scratch(b) 2.5nm scratch(c) 5.0nm scratch(d) 7.5nm scratch

closeupcloseupcloseupzxycloseupFig.5.25 Snapshots of ＠C2160 -SL under scratch.

第6章
結論
本研究では OLC 薄膜の摩擦特性について評価するために，まず単一フラーレン並びに単一 OLC の圧縮特性，摩擦特性を分子動力学シミュレーションにより検討した後,
フラーレン・OLC を平面上に並べた薄膜構造に対して摩擦シミュレーションを行った.
以下に，得られた結果を総括する．
第 2 章では，本研究で用いた解析手法の基礎について述べた．まず，分子動力学法の概要ならびに基礎方程式を示し，本研究で用いた数値積分法について説明した．次に，原子間相互作用の評価に用いられるポテンシャルエネルギーについて述べ，炭素原子に関するポテンシャル関数を具体的に説明した．さらに，大規模シミュレーショ
ンに必要な計算の高速化手法, ならびに Stone-Wales 欠陥を導入したフラーレンの幾何形状について述べた．
第 3 章では径の大きさの異なる単一フラーレンに対して，初期構造緩和，圧縮，摩擦のシミュレーションを行った．初期構造緩和シミュレーションにより，作成した C60の半径は実験値の約 0.35nm とほぼ等しく，原子数を 60n2 の形で表した場合の n の値
が 1 異なるフラーレン同士の半径の差はグラファイトの層間距離に近いことが示された．圧縮シミュレーションでは，頂点保持による圧縮により，各フラーレンの座屈には構成単位近傍の五員環が関係していることがわかった．そうした座屈を生じるときの応力は高次のフラーレンほど低くなっており，bond 数の増加による変形の自由度の増加がその理由として挙げられる. 平滑なダイヤモンド壁による圧縮では，C540 は圧縮された上下面が接近して相互作用し始めるとき，C960 ，C1500 ，C2160 は折りたたま

86

第6章

結論

87

れた内部の隙間がなくなったときに応力が急上昇することが示された．ダイヤモンド壁を用いた摩擦シミュレーションでは，フラーレンの摩擦係数は 10−2 程度のオーダーであること，その大小関係は径の大きさや押し込み量では整理することが難しいこと，を明らかにした．なお摩擦中の挙動は，摩擦方向への回転は生じず滑るように移動している．
第 4 章では径の大きさの異なる単一 OLC に対して，3 章と同様のシミュレーションを行い，フラーレンとの違いについて議論した．初期構造緩和シミュレーションを行った結果，作成した OLC の層間距離はフラーレン単体で比べたときの半径の差と変りないことが示された．頂点保持による圧縮シミュレーションでは，OLC は各層のフラーレンの座屈応答を受けながらも，多層構造を有することであまり応力が減少することなく，フラーレンよりも高い抵抗力を示すこと，また圧縮時は内部のフラーレンほど強い力を受けることがわかった．ダイヤモンド壁による圧縮では，上下共平面となって変形するため各層のフラーレンの座屈による応力減少は見られないこと，中心部のC60 が崩壊するまで圧縮した時に応力が急上昇することがわかった．同じ押し込み量で比較するとフラーレンより OLC の方が，径の大きなものより小さなものの方が高い応力を示していた．摩擦シミュレーションでは，OLC の摩擦係数はやはり 10−2 程度のオーダーであり，フラーレンの摩擦係数と同程度であった．多くはフラーレン同様，
摩擦方向に回転することなく滑るように移動していたが，1nm+0.5D(D:直径) 押し込
みした場合の ＠C960 ，＠C1500 は，中心と外側の層の密度の違いにより回転する様子が見られた．
第 5 章ではフラーレン・OLC を同程度の被覆率になるように平面上に並べた薄膜構造に対し，圧子や基板の形状を変えた摩擦シミュレーションを行った．平滑な基板上に配置した C540 ，C2160 ，＠C540 ，＠C2160 に対して，のこ歯状の切れ込みのある圧子(大小 2 通り) を用いて摩擦を行った場合，圧子の形状に関係なく摩擦係数は 10−2 程度のオーダーであり，その大小関係は ＠C2160 < C2160 < C540 < ＠C540 であった．このとき摩擦中の挙動はいずれも滑るような移動であった．基板にものこ歯状の切れ込みを入れて同様のシミュレーションを行った場合，摩擦係数のオーダーは 10−1 程度と一桁増加し，その大小関係も ＠C2160 < ＠C540 < C2160 < C540 と OLC よりもフラーレ

第6章

結論

88

ンの方が高い値を示すことがわかった．摩擦過程を確認すると，回転量に差はあるものの全てのモデルが摩擦中に回転していた．フラーレンの方が圧子だけでなく基板の切れ込みにも深く入り込みやすいため，回転時大きな変形を受け大きな抵抗力を示す
ことが摩擦係数に表れたものと考えられる．
以上のように，フラーレン・OLC 薄膜の摩擦係数は単体の径の大きさだけでなく，基板形状や圧子形状によってその値を大きく左右される．中でも実現象に近い圧子，基板共に表面起伏のあるダイヤモンド壁を用いた場合では，3 章，4 章の平滑なダイヤモ
ンド壁を用いた場合より摩擦係数は格段に増加しており，対象 (フラーレンか OLC か)の違いや径の大きさの違いによる摩擦係数の傾向や対象が摩擦方向に回転する様子をはっきりと確認することができた．今後はより実現象に近い対象を積層させたモデルに対するシミュレーション，表面に曲面を導入した圧子や基板を用いてのシミュレー
ションといった検討が望まれる．

参考文献[1] Kroto, H.W., Heath,J.R., O’Brien, S.C., Curl, R.F., and Smalley, R.E., Nature ,318, (1985), 162.
[2] Iijima, S., Nature, 354, (1991), 56.
[3] Ugarte, D., Nature, 359, (1992), 707.
[4] Kuznetsov, V.L., Butenko, Yu.V., Chuvilin, A.L., Romanenko, A.I., Okotrub,
A.V., Chem.Phys.Lett., 336, (2001), 397-404.
[5] 垣内, 平田, 精密工学会誌, 67, (2001), 1175.
[6] Kuznetsov, V.L., Chuvilin, A.L., Butenko, Yu.V. et al, Chem.Phys.Lett., 222,
(1994), 343.
[7] Hirata, A., Igarashi, M., Kaito, T., Tribol.Int., 37, (2004), 899-905.
[8] Matsumoto,N.,

Joly-Pottuz,

L.,

Kinoshita,

H.,

Ohmae,

N.,

Dia-

mond.Relat.Mater, 16, (2007), 1227-1230.
[9] Joly-Pottuz, L., Vacher, B., Ohmae, N. et al, Tribol.Lett., 30, (2008), 69-80.
[10] Joly-Pottuz, L., Kinoshita, H. et al, Tribol.Int., 41, (2008), 69-78.
[11] Cabioc’h, T., Riviere, J.P., Delafond, J., J.Mater.Sci., 30, (1995), 4787-4792.
[12] Banhart, F., Fuller, T., Redlich, Ph.D., Ajayan, P.M., Chem.Phys.Lett., 269,
(1997), 349-355.
[13] Cabioc’h, T., Jaouen, M. et al, Appl.Phys.Lett., 73, (1998), 3096-3098.
[14] 山口，丸山，機械学会論文集 (B 編)，63, 611, (1997), 2398.
[15] 塩見，五十嵐，丸山，日本機械学会論文集 (B 編)，76, 764, (2010), 642.
89

参考文献

90

[16] 尾方, 花生, 渋谷, 材料, 55, (2006), 754.
[17] Hirai, Y., et al., Jpn J Appl Phys., 42 6B, (2003), 4120.
[18] Deguchi, H., Yamaguchi, Y. et al, Chem.Phys.Lett., 503, (2011), 272-276.
[19] 西村，チュン，荒井，計算数理工学論文集，10, 13, (2010).
[20] Chico, L., Crespi, V.H., Benedict, L.X., Louie, S.G., Cohen, M.L., Phys.Rev.Lett.,
76, 6, (1996), 971.
[21] Lee, Y.H., Kim, S.G., Tomanek, D., Phys.Rev.Lett., 78, 12, (1997), 2393.
[22] 飯島，遠藤，ナノカーボンハンドブック, (2007), 531, 株式会社エヌ・ティー・エヌ.
[23] 篠原，斎藤, フラーレンの化学と物理, (1997), 79, 名古屋大学出版会.
[24] Saito, R. et al, Chem.Phys.Lett., 195, (1992), 537-542.
[25] Terrones, M. et al, Structural.Chem., 13, (2002), 373-384.
[26] Wang, B.-C., Wang, H.-W., Chang, J.-C., Tso, H.-C., J.Mol.Struct., 540, (2001),
171-176.
[27] Lu, J.P., Yang, W., Phys.Rev.B., 49, 16, (1994), 11421.
[28] Bakowies, D. et al, J.Am.Chem.Soc., 117, (1995), 10113-10118.
[29] Joly-Pottuz, L., Bucholz, E.W., Matsumoto, N. et al, Tribol.Lett., 37, (2010),
75-81.
[30] Bucholz, E.W., Phillpot, S.R. et al, Com.Mater.Sci., 54, (2012), 91-96.
[31] Brenner, D.W., Phys.Rev.B., 42, 15, (1990), 9458.
[32] Ahmadieh, A.A., Raﬁzadeh, H.A., Phys.Rev.B., 7, (1973), 4527.
[33] 上田顯, コンピューターシミュレーション, (1990), 朝倉書店.
[34] 西村英晃 神戸大学工学部機械工学科 卒業論文

(2011).

学術論文・学術講演
学術講演

 OLC の構造安定性並びに変形・破壊挙動に関する分子動力学シミュレーション
西村

英晃，屋代

如月

日本材料学会 第 1 回マルチスケールマテリアルモデリングシンポジウム (第 16回分子動力学シンポジウム), 大阪大学, (2011.5)

 分子動力学法による球殻形状の OLC に対する変形・破壊シミュレーション
西村

英晃，屋代

如月

日本機械学会 第 24 回計算力学講演会 CMD2011, 岡山大学, (2011.10)

 分子動力学法を用いたオニオンライクカーボンの力学的挙動に関する研究
西村

英晃，屋代

如月

日本機械学会 第 25 回計算力学講演会 CMD2012, 神戸ポートアイランド南地区,(2012.10)91学術論文・学術講演
論文貼り付け用ページ

92学術論文・学術講演
論文貼り付け用ページ

93学術論文・学術講演
論文貼り付け用ページ

94学術論文・学術講演
論文貼り付け用ページ

95学術論文・学術講演
論文貼り付け用ページ

96学術論文・学術講演
論文貼り付け用ページ97学術論文・学術講演
論文貼り付け用ページ98

謝辞
本研究を遂行するにあたり，仔細に渡る懇切丁寧な御指導を賜りました屋代如月准教授に，心より感謝の意を表します．浅学非才の著者に対し，三年間に渡り親切に，時に厳しく御指導頂いたことは人生の大きな財産となりました．また，本論文を完成するにあたり，研究報告会などにおいて様々な御助言を賜りました田中克志教授，学会などにおいて丁寧に御指導いただいた西村正臣助教 (現信州大学)，日々の研究活動が円滑に行えるよう便宜を図っていただきました古宇田由夫技術職員に深く感謝いたします．
若輩だった私に貴重な御助言を頂き，有益な議論を交わした坂田了介氏 (現三菱電気株式会社) に心より御礼申し上げます．日々の学生生活や大学院一年次の就職活動において，多くの有益なアドバイスを与えてくださった表面・界面工学研究室の諸先輩方に厚く御礼申し上げます．また，同じ研究グループとして活発な議論を交わした後輩の毛利友宙君に感謝いたします．
同じ研究室に配属されて以来共に研究に取り組み，互いに切磋琢磨し，時には励ましあった池宮一繁君，住谷昂大君，竹下和也君，中田伸哉君，藤原正大君に御礼申し上げます．また，同じ材料系として議論を交わし，多くの時間を共有した固体力学研
究室の同輩諸子に御礼申し上げます．
最後に，六年間の学生生活を暖かく見守り，経済的にも精神的にも支えて頂いた家族に最大の敬意と感謝の意を表します．ありがとうございました.
平成 25 年 2 月
西村99

英晃