— 1 —SENAC Vol. 48, No. 3（2015. 7）［共同研究成果］

X 線自由電子レーザーパルスによるフラーレン
超多価カチオン C60q+の爆発解離の動力学シミュレーション―

2 段階爆発機構の提案―

山崎 馨 1,2・上田 潔 3・河野 裕彦 1
1
東北大学大学院理学研究科化学専攻，
2
北海道大学大学院理学研究科化学部門，
3
東北大学多元物質科学研究所
X 線自由電子レーザー照射によって生成する超多価カチオン C60q+ (q = 20‒60) の解離機構を
Self-consistent charge density-functional based tight-binding (SCC-DFTB) 法に基づく分子動力学計算を用いて研究した．その結果，電荷 q によらず(1) 正電荷間のクーロン反発によって多価原子カチオンが 10‒20 fs 程度で非統計的に放出されるクーロン爆発と(2)１価・中性分子フラグメントの統計的脱離(100 fs から 1 ps 程度)の 2 段階機構で解離することを明らかにした．

1. 序論：X 線自由電子レーザーとクーロン爆発
我々が利用する医薬品は，有効成分の分子がタンパク質などの生体分子の特定の部位に作用して効果を示す物が多い．このため，タンパク質の構造を決めることは病気の発症機構を解明し有効な医薬品を開発するための第 1 歩である．我々が目に見えない小さなものを観測する際に使用する顕微鏡は使用する光の波長の半分程度の大きさの構造を見分けることができる．例えば，中学・高校で使用する普通の光学顕微鏡は，波長数百 nm の可視光を使用するため，μm 程度の大きさの細胞の内部構造を観測することができるが，細胞の中にある数十から数百 nm 程度のタンパク質の構造を観測することはできない．このため，ナノ分子の構造や化学結合の組み替えを観測できる新たな顕微鏡の開発が求められていた．この様な目的のために作られた X 線レーザーパル
ス顕微鏡が日本の兵庫県にある SACLA や米国カリフォルニア州にある LCLS などの X 線自由電子レーザー(XFEL)施設である．XFEL 施設では化学結合の長さと同程度(0.1 nm)かそれ以下の波長を持つ X 線を用いて今まで測定が困難であったタンパク質などのナノ分子 1 分子の構造や化学反応の様子を直接観察することを究極的な目標の一つとしている[1]．しかし，このような測定を可能にする高強度の X 線レーザーパルス(XFEL)1では，タンパク質分子が容易にイオン化され，分子内の強いクーロン反発によって爆発してしまうクーロン爆発という現象がおきる[1-3]．我々の体の中で有用な働きをしているタンパク質は，構造や反応の仕組みを決めるために質・量共に十分な試料を用意することが難しいものも多く，XFEL パルスによる試料分子のクーロン爆発は可能な限り回避しなければならない．よって，この爆発の影響を最小限に抑えるために，ナノ分子が X 線レーザーパルスによってどの程度の時間でどの様に爆発していくかを明らかにすることが喫緊の課題となっている．そこで筆者らはナノ分子としてフラーレン C60 を取り上げ，XFEL パルスの照射によって数十個の電子がはぎ取られてできた多価のフラーレン正イオン C60q+[2] 2がどの様に壊れていくかを，量子化学計算と分子動力学シミュレーションを組み合わせて調べた[4]．

1

典型的な反応動力学イメージング実験の場合，パルスの半値幅は 10 fs 前後，パルス 1 つに含まれる光子の数は 1 μm2 あたり 1011 個程度である[3]．
2
実験的には 100 個ぐらいまでの電子がはぎ取られる可能性が示唆されている．詳しくは[2]を参照．

— 2 —

SENAC Vol. 48, No. 3（2015. 7）

2. 計算手法
2.1 内殻イオン化とオージェイオン化の垂直イオン化モデル
XFEL パルスが C60 に照射されると，図 1 に示すように，強いクーロン引力で原子核近くに束縛された内殻 1s 軌道から電子が放出され，1s 軌道に非常に不安定な空孔ができる(内殻イオン化)．
この空孔を埋めるために，より外側にゆるく束縛された 2p 軌道などの価電子軌道の電子が 1s 軌
道へ移ると同時に価電子軌道の電子がもう 1 電子イオン化されるオージェイオン化がおこる．
XFEL パルスは非常に多くの光子を含んでいるために，この内殻イオン化・オージェイオン化は連続して複数回起き，C60q+が生成する[2]．1 回の内殻イオン化とオージェイオン化で計 2 個の電子が分子から放出されるので，q 価の多価カチオンを作るためには，q/2 回の内殻イオン化・オージェイオン化過程が必要である．このイオン化過程では，電子と原子核との相互作用によって，1回当り 10 eV，q/2 回で 5q eV の余剰振動エネルギーEin が原子核の自由度に注入される[5]．本研
究では， C60q+の解離機構が電荷 q に応じてどの様に変化するかを調べるため，まずこのイオン化過程を大幅に簡略化した垂直イオン化モデルを採用した．すなわち，C60 がサッカーボール状の平衡構造を保ったままシミュレーション開始時に瞬間的に q 価までイオン化され，それと同時にEin = 5q eV が注入されると仮定した．

図 1：X 線自由電子レーザーによる C60 の多重イオン化過程．1s 軌道の内殻イオン化とそれに続いて起きるオージェ緩和が q/2 回繰り返されて C60q+が生成する．

2.2 Self-consistent charge density-functional based tight binding (SCC-DFTB)法
本研究では，半経験的な分子軌道法の一つである Self-consistent charge density-functional based
tight binding (SCC-DFTB)法[6,7]を用い, C60q+に対して全エネルギー一定の条件で分子動力学シミュレーションを行った．この動力学計算には量子化学計算パッケージ DFTB+ 1.2.2[8]を用いた．
SCC-DFTB 法は密度汎関数(DFT)法を使って決められたパラメータを用いると共に分子の構造変化による電荷分布の揺らぎを考慮に入れることによって，DFT 法に匹敵する計算精度の計算を1/100 以下の計算コストで可能にしている．例えば，XFEL の光エネルギーが C6060+を生成するために全て使われた場合に起きると考えられる等方的な爆発過程では，生成する 60 個の C+カチオンの平均運動エネルギーは SCC-DFTB 法では 94 eV と計算され，代表的な DFT の汎関数である
B3LYP を用いた場合(90 eV)とよく一致する3．この際の計算時間は 16 コアをノード内並列で用いた場合，B3LYP 法で 500 fs あたり約 1 ヶ月，SCC-DFTB 法では 1 時間以下(設定にもよるが，15-303解離極限では，C+の電子配置は 2p 軌道に 1 個不対電子を持つラジカルになると考えられる．この様な開殻電子配置をより適切に記述するために, α スピンと β スピンの電子が入る分子軌道の形状を別々に最適化する非制限 DFT 法および spin-polarized DFTB 法[7]を用いた．また，多価カチオンの最高占有軌道と最低非占有軌道間のエネルギー差は非常に小さくなる．このため，占有軌道の電子が非占有軌道に遷移している励起電子配置が無視できなくなり，エネルギー計算の収束性が非常に悪くなる．この様な状況を改善するために，SCC-DFTB 法における計算では電子を温度 5000 K(~0.43 eV)で被占有軌道にも Fermi-Dirac 分布させることによって（励起電子配置の効果を熱平均する），収束性を改善すると共に低エネルギー励起電子配置の効果を現象論的に取り込んだ．

— 3 —

X 線自由電子レーザーパルスによるフラーレン超多価カチオン C60q + の爆発解離の動力学シミュレーション分程度)となった．そこで本研究では q = 20, 28, 60 において，Ein = 5q eV を注入した場合の全長10 ps のトラジェクトリをそれぞれ 31 回計算し，C60q+の解離機構を統計的に解析した．

2.3 フラグメントの解離判別法および形式電荷 zn の割り当て
トラジェクトリごとに C–C 間の結合距離 rC–C を逐次計算してフラグメントの解離判別を行った.今回は，rc–c ≤ 3.2 Å となる炭素原子を一つのフラグメントとしてまとめた4．なお，rc–c = 3.2 Åという値は， 直線 Cn+ (n = 2–5)フラグメントの最低 2 重項状態において，C–C 結合の開裂が完了して SCC-DFTB 法によるポテンシャルエネルギー曲線がほぼ平坦になる結合長である.
SCC-DFTB 法によって計算されたフラグメントのマリケン電荷は非整数の値をとる．このため，フラグメントの分岐比を計算するには この非整数電荷から整数の形式電荷 zn へと何らかのモデルを用いて変換してやる必要がある．我々はある時間において n 個の炭素原子からなるフラグメントの形式電荷が zn と zn + 1 の二つしかとらないと仮定した．さらに，zn と zn + 1 は SCC-DFTB法によって求められたサイズ n のフラグメントのマリケン電荷の平均値 ⟨zn⟩に対して zn < ⟨zn⟩
< zn + 1 の関係を持つと仮定した．以上 2 つの条件の下で，サイズ n のフラグメントのマリケン電荷の総計と形式電荷の総計の差の絶対値が最小になるように，zn と zn + 1 の電荷を持つフラグメントの数を決定した．詳細については[4]の Supplemental Material を参照されたい．

3．C60q+ (q = 20-60) 多価カチオンの 2 段階爆発機構[4]本研究で， q ≥ 20 および Ein ≥ 100 eV の条件において，C60q+ 多価カチオンの爆発的解離が 2段階で起きることを見いだした．この 2 段階機構では，図 2 に示す様に多価の原子フラグメントCz+ (z ≥ 2)が 10 fs 程度の時間スケールで放出される．この Cz+ の放出に伴って多くの正電荷が放出され，解離せずに残っている core cluster と呼ばれる構造が静電的に安定化される．この corecluster は，内部に残っている振動(熱)エネルギーによって 100 fs から 1 ps 程度の時間スケールで1 価の原子フラグメントや C2+，C3+などの小さな分子フラグメントに分解する．ここでは具体例
として， C6060+ に Ein = 300 (= 5q) eV を注入した場合について詳しく見ていくことにする.
クーロン 爆発(t ≥ 10 fs)

統計的 脱離C4+

C2+Core cluster

t =0 fs~12 fs~25 fs

図 2 フラーレン超多価カチオン C60q+ (q = 20‒60)における 2 段階解離機構の q = 60 の例．多価の原子カチオン Cz+( z = 2−4)がクーロン反発によって 10‒20 fs 程度で放出された後に，解離せずに残っている炭素クラスターから Cnz+(n ≥ 1, z = 1, 2)が統計的に脱離する5．

4

このグループ化は，階層クラスター解析法の一つである Minimum linkage アルゴリズム（例えば文献[9,10]を参照）を応用して実装した．
5
なお，この爆発過程のトラジェクトリはサイバーサイエンスセンターの協力により 3 次元動画として可視化されており，本センターにおいて閲覧することが可能である．

— 4 —

SENAC Vol. 48, No. 3（2015. 7）

図 3：C6060+ に Ein = 300 (= 5q) eV を注入した場合の(a)フラグメントの種類と数の時間変化および(b)シミュレーション開始から 10 ps 後におけるフラグメントの並進運動エネルギー分布．
図 2 および図 3(a)に示す様に，C6060+ の爆発過程は t ~ 12 fs からの C4+ 放出によって始まる．
この原子性フラグメントの放出過程は正電荷を持つ炭素原子間の強いクーロン反発によるものである. このクーロン爆発過程に引き続いて core cluster Cnz+ (n ~ 25 and z ~ 5 ) が Cz+ (z = 1‒3) に分解する過程 t ~ 25 fs 付近から始まる． C–C 間のクーロン反発エネルギーを計算すると，シミュレーションを開始した t = 0 fs における C6060+ では結合一本当り 69.8 eV であった．これは，典型的な C=C 2 重結合を回裂させるために必要な解離エネルギー(6.3‒6.5 eV)の約 10 倍に相当する．一方クーロン爆発が終了して分子フラグメントが放出されるようになる t = 25 fs では，core cluster(~C255+)におけるクーロンエネルギーは結合 1 本当り約 0.5 eV まで減少していた. なお，0.5 eV という値は C=C 2 重結合の解離エネルギーのおよそ 10 % にしか相当しない. その代わり corecluster の中には 24 本の C=C 2 重結合を切断できるエネルギーに相当する 152.3 eV(4.9 万℃に相当)もの振動・回転エネルギーが熱として蓄えられており，この熱エネルギーを使って分子フラグメントの生成が進行する．以上より，分子フラグメントの脱離は，core cluster 内の振動・回転エネルギーがその内部で統計的に分配されて起きる過程だと考えられる．
この 2 段階爆発機構をより詳細に解析するために，筆者らは解離フラグメントの並進運動エネルギー分布を解析した．一般に，解離する親分子(この場合は C6060+)が高いエネルギーを持っている状態から放出されたフラグメントは高い並進運動エネルギーを持ち，系からのエネルギー放出が進んで系内のエネルギーが少なくなってから放出されるフラグメントの並進運動エネルギーは低くなる傾向がある．我々の統計解析の結果，各フラグメントの並進運動エネルギーの分布関数 fは式(1)にフラグメントの並進運動エネルギー KT の関数である Maxwell-Boltzmann エネルギー分
布[11] 1 つまたは 2 つでフィッティングできることが明らかになった ．
����� � � � �����

2 ���√π��� �� ���� �����

���,�� ��(1)

ここで Ti と Ni は i 番目の運動エネルギー分布関数における実効的な並進運動温度と対応するフラグメントの相対強度である．また，kB は Boltzmann 定数を表す. Ti と Ni の値は，シミュレーションで得られた並進運動エネルギー分布に f (KT)を最小自乗法でフィッティングして求めた.
C6060+ に Ein = 300 eV を注入した場合におけるサイズ n (Cnz+, z = 0, 1, …) のフラグメントの並進運動エネルギー分布を図 3(b)に示す．C+ における並進運動温度の高速成分 kBT1 と低速成分 kBT2の値はそれぞれ 150.9, 14.3 eV と求められた．C+ の kBT1 の値は，電荷 q が小さくなるにつれて急激に減少する．これは， Xe25+衝突による C60q+のクーロン爆発実験において，C+ の運動エネルギー分布の高速成分の値が，q の減少と共に小さくなるという Tomita らの実験結果と一致す

X 線自由電子レーザーパルスによるフラーレン超多価カチオン C60q + の爆発解離の動力学シミュレーション— 5 —る[12]．また，n = 2, 3 の分子フラグメントにおける kBT1 と kBT2 はどちらもおよそ 10 eV であり，系全体の電荷 q および余剰振動エネルギー Ein によらず定性的には一定であると見なせることが分かった．このことは，Tomita らの実験で，分子フラグメントの平均並進運動エネルギーがq によらず 10‒20 eV とほぼ一定であったこと[12]に一致する. 以上より，C60q+ から生成したフラグメントの高速成分と低速成分はそれぞれクーロン爆発とそれに引き続いて起きる統計的な解離に起因することが明らかになった. この 2 段階爆発機構を実験的に証明するためには，並進運動エネルギーの大きな C+ と低エネルギーの Cn+ (n = 1‒3)を同時に質量分析装置で観測するか，2つ以上の種類のフラグメントの運動量または並進運動エネルギーの相関を解析すればよい．また，
図 3(a)に示したフラグメント数の時間変化については XFEL パルスでクーロン爆発を引き起こし，紫外線パルスでフラグメント数の変化を追跡する XFEL ポンプ―紫外プローブ分光法で実験的に観測できると期待される.

4．結論
本研究では，X 線自由電子レーザー照射によって生成する超多価カチオン C60q+ (q = 20‒60) の解離機構を Self-consistent charge density-functional based tight-binding (SCC-DFTB) 法に基づく分子動力学計算を用いて研究した．その結果，電荷 q によらず(1) 正電荷間のクーロン反発によって多価原子カチオンが 10‒20 fs 程度で非統計的に放出されるクーロン爆発と(2)1 価・中性分子フラグメントの統計的脱離(100 fs から 1 ps 程度)の 2 段階機構で解離することを明らかにした．本機
構は幅広い電荷において成り立つため，タンパク質などの他のナノ・バイオ分子に対しても一般的に成り立つことが期待される．このため，この 2 段階爆発機構の適応範囲を理論と実験の緊密なインタープレーによって明らかにすることで，ナノ・バイオ分子の XFEL による反応ダイナミクスイメージングの実用化に繋がることが期待される．

謝辞
本研究の計算の一部は東北大学サイバーサイエンスセンターの並列コンピュータを利用することで行った．また，本センターの山下毅氏には，DFTB+1.2.2 の並列コンピュータへの最適化および計算結果の 3 次元可視化についてご協力いただいた．なお，本研究は日本学術振興会特別研究員研究奨励費（山崎：No. 251672），科研費（河野：No.24245001），文部科学省 X 線自由電子レーザー重点戦略研究課題（上田）の支援によって行われた．

参考文献[1] K. J. Gaffney & H. N. Chapman, “Imaging Atomic Structure and Dynamics with Ultrafast X-rayScattering,” Science 316, 1444-1448 (2007).
[2] B. F. Murphy, K. Ueda et al., “Femtosecond X-ray-induced explosion of C60 at extreme intensity,”
Nature Commun. 5, 4281 (2014).
[3] M. Hoener et al., Phys. Rev. Lett. 104, 253002 (2010).
[4] K. Yamazaki, T. Nakamura, N. Niitsu, M. Kanno, K. Ueda, H. Kono, “Communication: Two-stepexplosion processes of highly charged fullerene cations C60 q + (q = 20–60) ,” J. Chem. Phys. 141,121105 (2014).
[5] K. Endo et al.“Analysis of Electron Spectra of Carbon Allotropes (Diamond, Graphite, Fullerene) by
Density Functional Theory Calculations Using the Model Molecules,” J. Phys. Chem. A 107,9403-9408 (2003).
[6] M. Elstner, D. Porezag, G. Jungnickel, J. Elsner, M. Haugk, Th. Frauenheim, S. Suhai, and G. Seifert“Self-consistent-charge density-functional tight-binding method for simulations of complex materials

— 6 —

SENAC Vol. 48, No. 3（2015. 7）

properties,” Phys. Rev. B 58, 7260-7268 (1998).
[7] C. Köhler, G. Seifert, and T. Frauenheim, “Density-Functional based calculations for Fen (n≤32),”
Chem. Phys. 309, 23 (2005).
[8] B. Aradi , B. Hourahine , Th. Frauenheim, “DFTB+, a Sparse Matrix-Based Implementation of theDFTB Method,” J. Phys. Chem. A 111, 5678–5684 (2007).
[9] W.H. Press, S. A. Teukolsky, W. T. Vetterling, B. P. Flannery, Numerical Recipes (3rd ed.), Cambridge
Univ. Press (2007), Chap. 16.4
[10] 小西 貞則「多変量解析入門―線系から非線形へ―」岩波書店 (2010) 第 10 章．
[11] C. Nordling, & J. Österman, Physics handbook for science and engineering 8th ed. (Studentlitteratur:Lund, Sweden, 2006) pp 346[12] S. Tomita, H. Lebius, A. Brenac, F. Chandezon, B. A. Huber “Kinetic-energy release and fragment
distribution of exploding highly charged C60 molecules.” Phys. Rev. A 65, 053201 (2002).