2E11
分子動力学シミュレーションによるシクロデキストリンの水和構造の理論的研究
†

‡

§

†‡

†‡

（広島大院・理 , 広島大 QuLiS , PNNL ） 赤瀬 大 , , 相田 美砂子 , , Sotiris S. Xantheas§

A theoretical study on solvation structure of cyclodextrins
using MD simulation‡†( Center for Quantum Life Sciences, Graduate School of Science, Hiroshima Univ.,§

Pacific Northwest National Laboratory)†‡†‡

Dai Akase , , Misako Aida , , Sotiris S. Xantheas【序】 シクロデキストリン(CD)は α-D-グルコースが 1,4-グリコシド結§OH合によって結合した環状オリゴ糖であり、天然ではでんぷんに酵素が

O

O

OH

HO

作用することで生成される。最も一般的なシクロデキストリンは、グル

O
OH

HO
O

OH
O

OH
O

HO

O

らのシクロデキストリンはドーナツ型の分子構造をもち、その中心の

OH

OH

HO

O

O

OH

HO

α-CD

空洞に様々なゲスト分子を包接する。その包接機能を利用して、シクOH

ロデキストリンは食品、化粧品、医薬品、家庭用品など様々な分野にHO広く応用されている。

O

HO

HO
O

O
OH

O

OH

HO

HO

O

OH

OH

O

の外側は親水的であるが、空洞内部は幾分疎水的である。水溶液

OH

O
HO

O

シクロデキストリン分子は OH 基を残基あたり 3 つもつため分子HO中ではこの疎水的な空洞に溶媒水分子が取り込まれており、ゲスト

O

HO
O

OH

HO

O

HO

OH

O
O

O

HO
OH

HO

分子が空洞内の不安定な水分子を置換することが包接体生成の推β-CD進力と考えられている。したがって、シクロデキストリン分子の水和構OHHO造を明らかにすることは包接機構の理解に重要である。

OH

我々はこれまでに量子化学計算を用いて気相中におけるシクロデ

HO

キストリン分子および 2 量体の様々なコンフォマーの安定構造を計

HO

では、分子動力学(MD)法を用いて α-, β-, γ-シクロデキストリン水

OH

O

HO

O
HO

コース 6,7,8 量体の α-, β-, γ-シクロデキストリンである(図 1)。これ算し、気相中におけるシクロデキストリンの安定構造を得た。本研究

HO
O

HO

O

O
O

HO
O

OH

O

O

HO

OH

HO

HO
OH

OH

O
OH
O

HO

OH

O
HO

OH
O

OH

O
HO

O

O

O
OH

O

O

HO

OH

HO

溶液のシミュレーションをおこない、水溶液中におけるシクロデキスト

γ-CD

リン分子の構造と気相中における安定構造との比較およびシクロデ

図 1 シクロデキストリンの化学構造

キストリン分子の水和構造について議論する。
【計算方法】 まず、シクロデキストリンの力場を評価するために気相中でシクロデキストリン分子の様々なコンフォマーの相対エネルギーを計算し、量子化学計算 (B3LYP, ωB97X-D, MP2)の結果と比較した。
シクロデキストリンの力場には CHARMM[1]と GLYCAM[2]を用いた。計算プログラムは NAMD を使用した。水溶液の MD シミュレーションはシクロデキストリン分子と水 5000 分子の系で立方体のユニットセルでおこなった。水分子の力場は TIP3P を使用した。周期的境界条件を用い、静電相互作用はParticle Mesh Ewald 法で、van der Waals 相互作用はスイッチング距離 10 Å、カットオフ距離 12
Å で取り扱った。また、水素原子の共有結合は SHAKE 法により固定し、タイムステップ 1.0 fs で時間発展させた。温度制御には Langevin dynamics、圧力制御には Nosé-Hoover Langevin piston 法を使用した。シクロデキストリンの初期構造は気相中における量子化学計算で安定であった 4 つのコンフォマーの構造を用いた。5000 ステップのエネルギー最小化の後に、NpT 一定(p = 1.01325 bar, T= 300 K)の MD シミュレーションをおこなった。さらに、NpT 一定の MD シミュレーションから得られた

V を用いて、NVT 一定の MD シミュレーションをおこなった。
【結果と考察】 表 1 にシクロデキストリン分子の様々なコン
フォマーの相対エネルギーの量子化学計算の結果に対する

表 1 気相中におけるコンフォマーの相対エネル
ギーの平均絶対誤差 (kJ/mol)na) CHARMM

systemGLYCAM

α-CD B3LYP

32

9.2

19.4

α-CD ωB97X-D

29

14.7

24.6

α-CD MP2

13

10.8

23.9

β-CD B3LYP

15

8.3

18.9

β-CD ωB97X-D

15

19.3

27.5

表 2 に NVT 一定の MD シミュレーションから得られたシβ-CD MP2

6

4.9

22.3

クロデキストリン分子の水素結合の数の平均を示す。量子化

γ-CD B3LYP

22

12.8

21.0

学計算の結果から α-, β-, γ-シクロデキストリン分子の気相

γ-CD ωB97X-D

22

21.4

29.6

中の最安定構造は homodromic な水素結合ネットワーク

γ-CD MP2

6

7.5

24.8

力場の平均絶対誤差を示す。すべての場合で GLYCAM より CHARMM のほうが小さい平均絶対誤差を示した。この結果より、水溶液の MD シミュレーションには CHARMMを用いた。

をもちそれぞれ 12, 14, 16 本の分子内水素結合を形成していること分かっているが、水溶液中での MD シミュレーa) number of conformer

表 2 水素結合の数の平均ションではそのような水素結合ネットワークは形成されず、分子内水素結合の数の平均はかなり少なくなっている。溶媒
水分子との分子間水素結合及び合計の水素結合の数の平均は分子が大きくなるにつれて多くなるが、分子内水素結合

分子内

分子間

合計

α-CD

3.9

34.3

38.2

β-CD

5.8

38.7

44.5

γ-CD

5.3

46.2

51.6

は β-シクロデキストリンで最大となった。これは、γ-シクロデキストリン分子の柔軟性が高いため分子間水
素結合を多く形成できるコンフォメーションをとれるためであると考えられる。
あるトラジェクトリーにおける α-シクロデキストリンの分子内及び分子間水素結合の本数の時間変化を
図 2(a)に示す。約 t=3 ns から 9 ns の間で、分子内水素結合の本数が少なくなり分子間水素結合の本数が増えているのがわかる。このとき α-シクロデキストリン分子はドーナツ型の大環状構造が部分的に崩れ図 2(c)に示すような 1 つのグルコース残基の CH2OH 基が空洞内部に貫入する構造になっていることがわかった。この構造は約 t=9 ns に元の構造に戻った。同じトラジェクトリーにおける α-シクロデキス
トリン分子の空洞内の水分子の数の時間変化を図 2(b)に示す。空洞内には最大 4 個の水分子が存在しており、約 t=3 ns から 9 ns の分子内水素結合が減少している間は、空洞内の水分子の数も減少している。これは、CH2OH 基が貫入し空洞内に存在した水分子を排出したためである。このように、 α-シクロ
デキストリン分子は水溶液中で少なくとも 2 つの状態で存在することが示唆された。
(a)

(b)

(c)

図 2 (a) α-CD の分子内および分子間水素結合の本数, (b) α-CD の空洞内の水分子の数, (c) グルコース残基が空洞に貫入した α-CD の構造[1] O. Guvench, S. N. Greene, G. Kamath, J. W. Brady, R. M. Venable, R. W. Pastor
and J. Mackerell, Alexander D., J. Comput. Chem., 29, 2543 (2008).
[2] K. N. Kirschner, A. B. Yongye, S. M. Tschampel, J. Gonzalez-Outeirino, C. R.
Daniels, B. L. Foley and R. J. Woods, J. Comput. Chem., 29, 622 (2007).

	
  2E13

化学反応分子動力学シミュレーションの高速化について(東京大学) ○山下雄史

On fast algorithms for the reactive molecular
dynamics simulation(Univ. Tokyo) ○Takefumi Yamashita
近年、計算機の演算能力が高まり、より長い時間スケールの分子動力学(MD)シミュレ
ーションが可能になってきた。特に、MD 計算でタンパク質の本質的な構造変化の時間スケールに迫れたことで、タンパク質の機能解明が期待されるようになってきた。しかし、興味のある系が化学反応（化学結合の組み替え）を起こす可能性がある場合は、通常の古典 MD 計算で仮定される経験的力場が利用できないため、シミュレーションできる時間がより制限されている。例えば、化学反応を起こす系を記述する力場を求めるた
めに QM/MM 法を利用するとすれば、QM 計算に大きな計算リソースを割かれて長時間シ
ミュレーションは難しくなる。これは、溶液中のプロトン移動反応のように、連続的に過剰なプロトンがリレーされていくような現象を追跡したい場合、QM 領域が広くなりさらに難しくなるだろう。	
 これを解決する処方箋の１つは、化学反応を記述する経験的力場を構築することであろう。こうした経験的力場構築手法の１つが多状態経験的原子価結合(MS-EVB)法[1]であり、実際、我々は MS-EVB 法を用いて生体系や水溶液系におけるプロトン移動の分子
動力学シミュレーションを行い、溶液系や生命現象の分子論的理解を深めてきた。[2-4]こうした研究が可能になったのは、MS-EVB 法が QM/MM 法などのような直接的電子状態計算を含む方法に比べ計算コストが低いためである。しかしながら、MS-EVB 法を用いた MD 計算も通常の古典 MD 計算に比べると、計算コストは高くなる。これは、MS-EVB法において、スムーズな結合様式の切り替えを記述するために、量子力学計算をまねた行列形式を採用しているため、その要素の数の分だけ計算コストが増加するのである。	
 まず、我々は、行列要素が増加した際にどこが計算のボトルネックになっているのかを詳細に調べた。そうすると、現在、
周期境界条件を持つ系のクーロン力を計算する上で最も良く利用されている
Particle	
 Mesh	
 Ewald 法においては、クーロン相互作用の長距離項の計算であることが分かった。支配的な相互作用は近距離項であるので、長距離項に多くの計算資源

図１：ELRI 近似の精度をつぎ込むのは「もったいない」と感じられる。そこで、我々は、長距離項の計算コストを下げるために、計算結果に大きな影響を与えない程度の良い近似法（図１）を開発することにした。	
 本発表で紹介する近似方法は２つのアイデアとその組み合わせである。第 1 のアイデ
アは、有効電荷を導入するもので、Effective	 Long-range	 Interaction(ELRI)法と呼ぶ。
第２のアイデアは、長距離力の効果を近距離項の中に繰り込むもので、Coarse-Graining	
 in	
 Interaction	
 Space(CGIS)法と呼ぶ。ELRI 法では、MS-EVB 行列の各要素に対して遠距離力計算をおこなわず、有効電荷に対する遠距離力の計算で置き換えることで計算が速くなっている。CGIS 法では、そもそも近距離力計算しか使わないので計算が速い。
詳細な解説は当日の発表でおこなうが、最後に本近似を使うことで、MS-EVB 法を使った化学反応の MD 計算は、標準的な古典 MD 計算の数倍程度にしか計算コストがかからな
くなったことを強調しておきたい。[5]	
 References:[1] G. A. Voth, Acc. Chem. Res. 39 143 (2006), and more references therein[2] T. Yamashita and G. A. Voth, J. Phys. Chem. B 114 592 (2010)[3] T. Yamashita and G. A. Voth, J. Am. Chem. Soc. 134 1147 (2012)[4] J. Xu, T. Yamashita, N. Agmon, and G. A. Voth, J. Phys. Chem. B. (in press)[5] T. Yamashita, Y. Peng, C. Knight, G. A. Voth, J. Chem. Theory Comp. 8 4863 (2012)2E14

交換対称性を利用した QM/MM 法の境界条件（原子力機構 1、サッサリ大 2）○志賀 基之 1、マシア マルコ 2

QM/MM (quantum mechanics/molecular mechanics) 法は、着目する部分系(QM 領域)を第一原
理的な電子状態計算で扱い、その周囲の部分系(MM 領域)を経験的な古典力場計算で扱うという、ハイブリッド法である[1]。全電子計算と比べて計算量が格段に減るので、溶液や蛋白質など複雑多体系を扱うのに適している[2-5]。ただし、境界条件の決め方など、計算法としては未熟で発展途上の面もある。特に最近では、分子動力学法など、分子シミュレーションとの組み合わせにおいて、計算の信頼性や効率化が議論の対象になっている[6-10]。
本研究で取り上げる問題は、QM 領域と MM 領域の間を行き交う分子が存在するときの不都合である。例として、希薄なイオン水溶液系を考えよう。いま、イオンと第一水和圏の水分子に間に働く複雑な相互作用を第一原理計算で扱いたいので、イオンとイオンに最も近い
n 個の水分子を QM 領域、残りの m 個の水分子を MM 領域に設定し、m>>n とする。この条
件で分子動力学シミュレーションを実行すると、水和圏内の水とバルクの水が自然に入れ替わってしまう。もしこれを放置すれば、やがてイオンの周りには QM の水分子はいなくなり、MM の水分子により取り囲まれるだろう。m>>n の場合、交換したほうがエントロピー的に有利なので、イオンと水の相互作用の強さとは無関係に、時間が経てば交換はいつか必ず起こる。しかし、これでは第一水和圏をきちんと扱おうとした目的にそぐわない。
最も安易な解決法は、水和圏内の水とバルクの水の入れ替わる瞬間に、QM と MM のラベルを入れ替えることである。しかし、ふつう QM/MM ポテンシャルは QM 原子と MM 原子の座標交換に対して対称性が崩れているので、そこでポテンシャルが不連続になり、分子動力学シミュレーションに支障をきたす。これを防ぐのに、ポテンシャルの段差を和らげる方法がいくつか提案されている[11-16]。
本研究では、もう一つの解決法として、外部から荷重ポテンシャルを加えて交換を抑制する方法を提案する[17]。荷重ポテンシャルは人工的なものだから、いかなる場合も動的な情報は損なわれる。しかし、時間に依存しない熱力学量（イオン水溶液では、動径分布関数、水和エネルギーなど）を全く変えないようにすることが可能である。これは、厳密な定式化によれば、QM 原子と MM 原子の座標交換で得られるあらゆる配置に対して存在確率を適切に決め、それに応じて荷重ポテンシャルを与えればよい[17]。しかし、実際の QM/MM 分子動
力学計算では、次のように最も簡略化された“単交換近似”でも、かなりうまくいくことがわかった。これには、まず、原点（イオン水溶液ではイオンの位置とする）から最も遠い QM
原子 i への距離 ri と、原点から最も近い MM 原子 j への距離 r j の差 rij  rj  ri を計算する（図１）。次に、 rij を変数とした Fermi 関数の対数を荷重ポテンシャルとしてV bias  k BT log f ij ,f ij 

1
1  exp rij と与える。ここで、k B は Boltzmann 定数、T は温度であり、 f ij と f ji はそれぞれ配置(i,j)と(j,i)の現れる確率因子を表していて、和が 1 である。最後に、 V bias と、その力   R V bias を
QM/MM 計算において加えるだけである。これは実装が簡単で、計算負荷も極めて少ない。

図１：原子の座標交換。●は QM 原子、○は MM 原子を表す。
ここで、指数  は QM 原子と MM 原子の分離しやすさを制御する入力パラメータである。

 を無限大にすれば rij が常に負となり、完全に分離される。これは都合がよい反面、荷重ポテンシャルが急峻になってしまう問題がある[18]。したがって、実際の分子動力学計算では

 を大きくしすぎないようにすべきである。
以上の方法は、荷重ポテンシャルの原点の定義を変えれば、界面など、他の境界にも応用できる。また、QM/MM 法に限らず、ONIOM 法などハイブリッド法全般に使えるものである。
本発表では、多重時間ステップ法[6,7]と組み合わせた効率のよい QM/MM 分子動力学計算の具体例を紹介し、その有用性を示す予定である。
[1] A. Warshel, M. Levitt, J. Mol. Biol. 103, 227 (1976). [2] H. Lin, D. G. Truhlar, Theor. Chem. Acc.
117, 185 (2007). [3] H. M. Senn, W. Thiel, Angew. Chemie Intl. Ed. 48, 1198 (2009). [4] H. Hu, W.
Yang, Ann. Rev. Phys. Chem. 59 (2008). [5] M. Boeckmann, D. Marx, C. Peter, L. D. Site, K. Kremer,
N. L. Doltsnis, Phys. Chem. Chem. Phys. 13, 7604 (2011). [6] M. Shiga, M. Tachikawa, Mol. Simul.
33, 171 (2007). [7] T. K. Woo, P. Margl, P. R. Blochl, T. Ziegler, J. Phys. Chem. A 106, 1173 (2002).
[8] H. Takahashi, Y. Kawashima, T. Nitta, N. Matsubayasi, J. Chem. Phys. 123, 124504 (2005). [9] T.
Yamamoto, J. Chem. Phys. 129, 244104 (2008). [10] T. Janowski, K. Wolinski, P. Pulay, Chem. Phys.
Lett. 530, 1 (2012). [11] T. Kerdcharoen, K. Morokuma, Chem. Phys. Lett. 355, 257 (2002). [12] B. M.
Rode, T. S. Hofer, B. R. Randolf, C. F. Schwenk, D. Xenides, V. Vchirawongkwin, Theor. Chem. Acc.
115, 77 (2006). [13] A. Heyden, H. Lin, D. G. Truhlar, J. Phys. Chem. B 111, 2231 (2007). [14] R. E.
Bulo, B. Ensing, J. Sikkema, L. Visscher, J. Chem. Theory Comput. 5, 2212 (2009). [15] N. Takenaka,
Y. Kitamura, Y. Koyano, M. Nagaoka, J. Chem. Phys. 137, 024501 (2012). [16] N. Bernstein et al.
Phys. Chem. Chem. Phys. 14, 646 (2012). [18] M. Shiga and M. Masia, J. Chem. Phys. 139, 144120(2013). [17] C. N. Rowley and B. Roux, J. Chem. Theor. Comput. 8, 3526 (2012).

2E15
超高速光化学的電子環状反応機構に関する新しい観点（三菱化学科学技術研究センター ）・理研 ・原研 ）○小林高雄 ，中村振一郎 ，志賀基之
a

b)

c)

a)

b)

c)

A New Perspective of Ultrafast Photochemical Electrocyclic Reaction Mechanism(Mitsubishi Chemical Group Science and Technology Reseach Center, Inc.a）, RIKENb）, JAEAc）)○Takao Kobayashi, Shinichiro Nakamura, Motoyuki Shiga【序】超高速光化学反応機構において、非断熱遷移が重要な役割を果たすことが知られている。これを理論的に理解するにあたり、“静的”なポテンシャルエネルギー曲面(PES)（ここではポテンシャルは
断熱ポテンシャルを意味することとする）の情報を使った解析が第一歩である。そこでは、最小エネ
ルギー経路の概念が導入され、それを繋ぐ（特に最小エネルギー）円錐交差(CI)点（エネルギーギャップ(∆E)=0）を求めることが常套手段になっている。しかし、非断熱遷移は、本来、“動的”な情報に基づいて理解されなければならない。なぜなら、系が大きな運動量を持っている場合には、遷移は円錐交差点近傍のみで生じるとは限らず、より広範に起こる可能性もあるからである。これは、運動量が PES
上の軌跡を決める運動方程式に関係するほかに、PES 間の遷移確率にも関係しているからである。さ
らに、ある地点での遷移確率はその地点に滞在する時間に比例する。例えば、２準位モデルにおいてTully の fewest switches (TFS) アルゴリズムでは時間ステップ∆t の間の状態 2 から状態 1 への遷移（スイッチ）確率 p （ここではスイッチを遷移と呼ぶ）は21

p 21 = J 12 ∆t

J 12 =(

a& 11− 2 Re c 1* c 2 v ⋅ d 12
=
2
a 22
c2

)

（ c , a はそれぞれ状態 i の振幅と分布、 v は核の速度、 d は状態 1,2 間の非断熱結合ベクトル）で与えられるが、軌跡におけるの J （単位時間当たりの遷移確率）の出現に関しては分布( ρ ( J ) )が存在し、軌跡が状態 2 にある場合の J 地点での遷移確率密度( P ( J ) )は、 τ を状態 2 の滞在時間とした
とき、 ρ (J )τ を J 地点における滞在時間と見なすことができ、概ね P ( J ) = J ρ ( J )τ と見積もることができる。この ρ ( J ) により、Massey パラメーター ζ = ∆E /(hv ⋅ d ) が 1 より大きい地
点でも遷移確率は高くなる可能性があり、特に v ⋅ h （１次の状態間結合 h = 1 ∇ H 2 = ∆E d ）大きい場合は、 ∆E がかなりの大きさにまで遷移が生じることになる。また、遷移が起こるときの運動量ベクトルは、その後の軌跡の分岐の仕方に関係するので、最 2A
1B
終的な反応の分岐比にまで影響が及ぶ可能性もある。
本研究では、代表的な光化学的電子環状反応の１つである 1,3シクロヘキサジエン（CHD）から 1,3,5-ヘキサトリエン（HT）への光開環反応を実例として、非断熱遷移を考慮した ab initio 分子 hν
動力学計算を行い、以上の動的効果を検証した。その結果、(1) 1A
滞在時間分布と v ⋅ h が大きいことに起因して、 ∆E がかなりの大きさに於いても非断熱遷移確率が生じること (2)非断熱遷移する際の軌跡の分岐は、遷移する際の分子構造そのものとは相関
CHD
cZccZc-HTが小さく、むしろ構造の変化速度と相関していることがわかった。
i

ii

12

12

12

12

12

2

21

2

12

12

21

12

21

12

21

12

12

2

12

12

R

12

21

12

21

S2O2O-FC

S2C2C-FC

S1O1O-FC

S1C1C-FC

S2/S1 CIO

S2/S1 CIC
S1P1P-C2

S1C1C-MIN

S1/S0 CIP-MIN
S1P1P-MIN
S0-TS

12

21

S0C0C-MIN

S0O0O-MIN

図１ CHD の光開環反応経路に沿った PES以下、その詳細を述べる。
【計算】CHD は 1B(S )状態へ光励起すると、S から S へ超高速で緩和して 200fs 程度で cZc-HT が生成することが知られている。ここで、S /S /S（C2 対象制限下では 1A(covalent) / 1B(ionic) / 2A(covalent)）の３状態が関与する（図１参照）が、反応経路全域に渡って、covalent と ionic２つの性格の状態が混合した３状態の PES の相対エネルギー関係を正確に記述するには、本来であれば動的電子相関を取り込む必要があるので、CASPT2 レベルで PES を求める必要がある。しかし、CAS/RAS の様々な活性空間を調査した結果、RASSCF(4,4)レベルの PES が MS-CAS(6,6)PT2 レベルの PES との間に、そのトポ
ロジー、ペリ環状反応において重要な点である S ペリ環状最小点(S )や最小エネルギーCI 点の構
造や性格、エネルギー関係において定性的な類似性があることを見出した。RASSCF(4,4)レベルの PES
上であれば、CHD の 1B Franck-Condon 状態から on-the-fly 非断熱ダイナミクス計算を行うことが可能であるので、それを実行した。半古典非断熱ダイナミクス法としては Zhu-Truhlar の SCDM(self-consistent
decay-of-mixing)法を採用し、初期座標・速度は CHD の S 安定点における基準振動モードランダムサ
ンプリングにより発生させ、初期エネルギーとして各モードにその零点振動エネルギーの半分を与えた。基底関数は 6-31G*、電子状態計算には Molpro2006 を用いた。
【結果】図２は S →S 遷移（最後の S →S スイッチ）における ∆Eと hv⋅d （状態 1,2 をそれぞれ S , S とする）の値をプロットしたものであるが、概ね反比例関係になっている。これは非対角
Hellmann-Feynman 式から導かれる関係式 hv ⋅ d = hv ⋅ h / ∆E が成立することによるものであり、 hv ⋅h = 25 − 250 (kcal/mol) である。この分布は予想に反して、遷移の場所には大きな分散があり、 図２ S →S 遷移における ∆E と hv ⋅ dMassey パラメーター ζ ≤ 1 の領域だけでなく、ζ > 1 の領域でも遷移が多い。これ原因を調べるべく、図３に示したように J とMasssey パラメーターの逆数 ζ はほぼ線型に相関していることから、 J の代わりに ζ の関数として、S S 遷移した軌跡数と全軌跡の S における滞在時間比率（ ρ (ζ ) ）を調べた。その結果が図４である。この結果から、この系に於いては S の滞在時間は ζが小さなところで極めて長く、 ζ が大きくなると急速に短くなっており、このことが ζ が小さい（すなわち ζ > 1 ）地点でも遷移を可能にしていることが明らかになった。さらに、全体的に 図３ S →S 遷移における J と ζv ⋅ h の絶対値が大きなことが、比較的大きな ∆E での遷移を可能にしている。
また、ここでは紙幅の都合上、結果は示さないが、Sへ非断熱遷移する際の「分岐と構造及び構造変化速度との相関」に関しては、非断熱遷移する際の分岐と構造との間に相関はほとんど無いが、分岐と構造変化速度との間には極めて明確な相関があることが見出された。
1

1

0

1

0

2

1

1-PMIN

0

0

1

12

0

0

21

1

12

12

hbar∗|v∗d12| at last switch (kcal/mol)

1

100
95
90
85
80
75
70
65
60
55
50
45
40
35
30
25
20
15
10
5
0

21

2

12

Cl-trajectory
Op-trajectory
y=250/xy=25/x
ζ21=1

0

5

10

15

20

25

30

35

40

45

50

55

60

∆E21 at S1−>S0 last switch (kcal/mol)

1

0

21

12

10

21

21

Cl-trajectory
Op-trajectory

12

−1
21

12

1

J12(1/fs)

1

−1
21

0

−1
21

1

0.1

−1
21

1

−1
21

0.01

21

0.01

0.1
1
10
1/ζ21 (hbar*∗|v∗d12|/∆E21)

1

12

−1
21

14

21

0

Number of last switch events

12

0

100

1e+0

13

nlast-switch(cl-traj)

12

nlast-switch(op-traj)11

ρ(1/ζ21) on S1

10

1e-1

1e-2

9
8
7

1e-3

6
5

1e-4

4
3

1e-5

2
1
0

1e-6

0.0

0.5

1.0

1.5

2.0

2.5

3.0

3.5

4.0

4.5

1/ζ21(=hbar∗|v∗d12|/∆E21)

5.0

図４ S →S 遷移した軌跡数と滞在時間比率の ζ 依存性
1

0

−1
21

ρ(1/ζ21) on S1

−1
21

2E17
量子論と半古典論による ICl 分子の光解離過程における量子干渉効果(慶大院理工) ○松岡 貴英，藪下 聡
Quantum Interference Effect in the Photodissociation Processes of ICl molecule with
Quantum and Semi-classical Approaches(Keio Univ.) ○ Takahide Matsuoka, Satoshi Yabushita【序】ICl の第 2 吸収帯(200~300 nm)に含まれる断熱電子状態のうちΩ=0+ 対称のものは，それぞれ
23Π0+(2341)と 3Σ-0+(2422)を 主 配 置 と す る 0+(III),
0+(IV)であり，その遷移双極子モーメント（TDM）
は二乗値 0.0179, 0.0099 bohr2 を持ち，図 1, 2 のように各解離チャネルに相関する[1]．ここで(pqrs)
は(σpπqπ∗rσ∗s)の電子配置を意味する．第 2 吸収帯
	
  	
  	
  	
  	
  	
  	
  	
  	
  	
  	
  	
  	
  	
  	
  	
  	
 
の 光 解 離 生 成 物 は ， 振 動 電 場 に 垂 直 方 向 で は 図１．0+(III)と 0+(IV)の相関図．Av-1 において+
I+Cl*が，また平行方向では I*+Cl が支配的に観 2 つの状態に分岐する．0 (III)に分岐した成分
は再び Av-2 において分岐する．
測されている（I=I(2P3/2), I*=I(2P1/2)など）[2,3]．
つまり前者はΩ=1 の対称性の励起状態で解離した生成物で，後者はΩ=0+の励起状態で解離した生成物である．I*+Cl に透熱的に相関するのは 0+(III)であるので，0+(III)の吸収強度は 0+(IV)のそれよ
りも非常に大きい，と実験から予想された．
	
  以前の理論計算[1]で，0+(III), 0+(IV)状態はそれぞれ断熱的に I+Cl*と I*+Cl に相関することを踏
まえ，吸収断面積とこれら 2 状態間に存在する回避交差(Av-1)における非断熱遷移確率から解離生
成物の分岐比を計算したところ，I*+Cl と I+Cl の!
比率は 2:1 となり実験値と矛盾した．ただしこの 図２．𝜴 = 𝟎 状態のポテンシャルエネルギー
理論計算においては断熱状態間の干渉効果を含め 曲線．
ていなかった．しかし本来は回避交差 Av-1 における 0+(III)と 0+(IV)状態の de Broglie 波の位相差Δϕに依存して、解離波束の間に量子干渉が生じるため(図 1)，これを含めた評価が必要である．Av-1 に至るまでの 2 つの断熱状態の deBroglie 位相は，Young の 2 重スリットにおける波の位相に対応する．
	  本研究では，この干渉効果を含めることによって実験と理論計算の間の矛盾が解消されることを示す．特に，Av-1 において 0+(III)状態の確率振幅が弱められ 0+(IV)状態のそれが強められることについて，その詳細を量子論と半古典論を用いて分析を行う．
【 理 論 と 計 算 】 RECP 法 に よ り ス ピ ン 軌 道 (SO) 相 互 作 用 を 含 ん だ 電 子 状 態 計 算 を ，COLUMBUS を用いて行った．ポテンシャルエネルギー曲線(PEC)と 1 次非断熱結合項(NACT)は SO 相互作用を含む多参照 2 電子励起配置間相互作用法により計算し，TDM は 1
電子励起配置までを考慮して計算した．これらの情報を用いて，Chebychev 展開法による波束計算，還元散乱行列を用いた半古典論による計算，および半古典的取り扱いに基づいたClassical Path Method[4]による計算を行った．
【結果と考察】0+(IV)状態の主配置は 3Σ-0+(2422)であり，1Σ+(2440) が主配置の基底状態 X(0+)
からは 2 電子励起である．しかしスピン軌道相互作用による 23Π0+(2341)からの intensity
borrowing のため 0+(III)と同程度の吸光係数が計算から得られた．第 2 吸収帯の同じ励起エネルギーにおいて，これら 2 状態の吸光係数は共に有限の大きさを持ち強く重なり合い，また Av-1 における非断熱遷移確率は p1 ≅ 0.80であった．第 2 吸収帯における 0+(III)と 0+(IV)の散乱部分断面積(それぞれ断熱的に I+Cl*と I*+Cl に相関)は，回避交差 Av-1, Av-2 における非断熱遷移確率 p1, p2 と基底状態から 0+(III), 0+(IV)状態への TDM𝜇III , 𝜇IV を用いて，次のよ
うに表現される．
𝜎I+Cl∗ = (1 − 𝑝! ) 𝑝! 𝜇IV + 1 − 𝑝! 𝜇III − 2𝜇III 𝜇IV 𝑝! 1 − 𝑝! cos  (Δ𝜙)

(1)

(2)
𝜎I* +Cl = 1 − 𝑝! 𝜇IV + 𝑝! µμ!!! + 2𝜇III 𝜇IV 𝑝! (1 − 𝑝! )cos  (Δ𝜙)+
解離断面積に寄与する干渉効果は2µμ!!! 𝜇!" 𝑝! (1 − 𝑝! )cos  (Δ𝜙)と表される( Δφ は 0 (III)と0+(IV)上の波が Av-1 通過時に持つ位相差（図１）)．特にこの振幅は，TDM の二乗値と遷移!
!
確率の積𝜇IV
𝑝! ，𝜇III(1 − 𝑝! )の幾何平均である．同じ平行成分の TDM を持ち，スペクトルが+強く重なる 0 (III), 0+(IV)は，非断熱遷移が完全に透熱でも断熱でもない，という干渉効果が強く生じる条件を全て満たす．Franck-Condon 領域で 0+(III)状態に遷移し回避交差を(1 − 𝑝! )!の確率で断熱的に通過した成分𝜇III(1 − 𝑝! )と，0+(IV)状態に遷移したのち𝑝! の確率で非断熱遷!
移した成分𝜇IV
𝑝! が，ほぼ同程度の大きさであるため，強い干渉効果が生じることになる．
位相差Δϕは，2 つの断熱状態の波数𝑘III,IV を
	
  	
  	
  	
  	
  	
  	
  	
  	
  	
  	
  	
  	
  	
 
座標𝑅について積分した作用積分の差として表される．2 つの断熱状態のポテンシャルエネ
ルギーの差が小さいほど，(3)式の位相差Δϕの
励 起 エ ネ ル ギ ー 依 存 性 ∂Δϕ/ ∂ΔEは 小 さ い ．
∂Δϕ/ ∂ΔEは，Av-1 を通過する前の断熱 PEC の
勾配差によって(4)式のように表現される．ここで，RIII , RIV は各励起エネルギーΔEにおける0+(III), 0+(IV)ポテンシャル曲線(PEC)の転向点
である．
!!
!!
𝜕Δ𝜙
𝜕
=
𝑘III d𝑅 −𝑘IV d𝑅(3)
𝜕Δ𝐸 𝜕𝛥𝐸 !III図３．I+Cl チャネル（赤）と I+Cl*チャネル
!IV
=

!!
!!

𝜇
2 Δ𝐸 − 𝑉

d𝑅III d𝑅IV
−
d𝑉
d𝑉
d𝑉

(4)

SOCI 計算結果によると 0+(III)と 0+(IV)の PEC
は FC 領域で近接し，回避交差 Av-1 までほぼ平行であった（図 2）．干渉効果が長波長側で顕著に寄与するため，I+Cl の断面積の平行成分が減少（異方性パラメータβ値も減少）し，実験値との不一致は解消された（図 3）．通常，
原子などの重粒子の干渉効果は励起エネルギ
ー依存性が強く，乱雑位相的に打ち消される場
合が多いが，0+(III), 0+(IV)の PEC の特異性により， (4)式の被積分関数は小さく，励起エネル
ギー依存性が弱いため顕在化したとみなせる．
	  Av-1 通過後における 0+(III)と 0+(IV)のノルムについて，いずれの計算方法も 0+(III)が弱められ，0+(IV)が強められる結果を示した（図 4）．
この干渉効果は 2 状態間の位相差Δϕを反映したものである．講演では，その位相差Δϕを決定する要因について述べる．

（青）の異方性パラメータβの計算結果（実線：干渉効果有り．点線：干渉効果無し）と実験値[2]（x）．

図４．断熱基底（実線）と透熱基底（破線）を用いて評価した Av-1 通過後における 0+(III)（黒・青）と 0+(IV)（赤・橙）のノルム．黒と赤の太線：波束計算．黒と赤の細線：半古典還元散乱行列．青と橙の細線：Classical Path Method.

[1] S. Yabushita, J. Mol. Struct.: THEOCHEM 461-462, 523(1999). [2] N. Diamantopoulou et al., J.
Chem. Phys. 134, 194314(2011). [3] L. Rogers et al., Chem. Phys. Lett. 258, 159(1996). [4] E. E.
Nikitin, Theory of Elementary Atomic and Molecular Processes in Gases (Clarendon, Oxford, 1974).

2E18
異種金属原子による窒素分子を挟んだ
逆サンドイッチ型錯体の特異なスピン状態(京大福井謙一研セ) ○中垣 雅之，榊 茂好

Unique Spin State of Inverted Sandwich Type Complexes with Hetero
Transition Metals Bridged by Dinitrogen Molecule(Kyoto univ. FIFC) ○Masayuki Nakagaki, Shigeyoshi Sakaki【緒言】β-ジケチミナト骨格を有する DDP (DDPH = 2-{(2,6diiso-propylphenyl)amino}pent-2-en; SCHEME 1)が配位した金属錯体は、有機分子を挟んだ逆サンドイッチ型二核錯体(ISTC)を形成する。これら ISTC は金属原子やサンドイッチされる分子の種類
によって多様なスピン状態をとり 1、その理解は理論化学のみならず材料科学の面からも興味深い。窒素分子を挟んだ錯体は特に興味深く、Cr 錯体では窒素分子は η2-side-on で配位し、１重項状態が最安定であるが 1、Fe 錯体では窒素分子は η1-end-on で配位し、７重項状態をとる 2。また Mn 錯体も Fe 錯体同様に高スピン(9 重項)状態を取る。これまでの研究で、金属の種類による構造、スピン状態の違いが、電子配置の違いによる金属原子‐窒素分子間の結合の違いに起因することを明らかにして来た。Cr 原子と Mn,Fe 原子では窒素分子との結合の性質が大きく異なることから、それらを組み合わせた異核錯体の電子状態は興味深い。本研究では Cr, Mn, Fe 原子を組み合わせた異核錯体の理論計算を行い、スピン状態及び窒素分子の配位構造の違いを検討した。
【計算方法】DDP の置換基を水素原子に置き換えた AIP (AIPH = 1-amino-3-imino-prop1-ene; SCHEME 1)を配位子とした錯体、(μ-N2)[M1(AIP)][M2(AIP)]2 {(M1,M2) = (Cr, Mn),(Cr,Fe), (Mn,Fe)}の各スピン状態について、CASSCF 法を用いて構造最適化を行い、CASPT2 法を用いてスピン状態の相対安定性を評価した。CASSCF 法の active space には、金属原子の3d 軌 道 と 窒 素 分 子 の ２ つ の * 軌 道 を 用 い た 。 遷 移 金 属 の 内 殻 電 子 は
Stuttgart-Dresden-Bonn の ECP で置き換え、原子価軌道には triple zeta 基底関数を用いた。
AIP の N 原子には aug-cc-pVDZ を、その他の原子には cc-pVDZ を用いた。
【結果及び考察】(Cr,Mn)錯体では B1 状態が基底状態であり η2-side-on 構造の 10 重項が最も安定である(表 1)。これは等核錯体で最も高いスピン多重度の Mn 錯体よりも高いスピン多重度であり興味深い。また A2 状態も存在するが、この場合、2 重項が基底状態である。A2
状態と B1 状態のエネルギー差は 4.4kcal/mol と比較的小さい。また、η1-end-on 構造では、
η2-side-on 構造の A2 状態に対応する B2 状態が、B1 状態よりも安定である。(Cr,Fe)錯体でも同様に 2 つの電子状態が得られたが、(Mn,Fe)錯体ではエネルギー差の小さい 2 つの電子状態は得られなかった。(Mn,Fe)錯体では 8 重項が基底状態であり、等核の(Mn,Mn)、(Fe,Fe)
錯体の中間のスピン多重度となる。

表 1 (μ-N2)[M1(AIP)][M2(AIP)]2 のスピン状態の相対エネルギー(kcal/mol)(M1,M2)(Cr,Mn)η1-end(Mn,Fe)η2-sideη1-endη2-sidemultiplicity B2 state B1 state A2 state B1 state A2 state B2 state
10tet
8tet
6tet
4tet
2let

13.5
7.6
8.2
9.5
11.4

14.1
15.5
16.0
16.2
16.1

18.3
6.2
5.7
5.5
4.4

0.0
1.4
2.3
2.8
3.0

26.5

12.1

0.0

1.7

5.7

6.4

9.6

9.7

20.2

16.8

(Cr,Mn)錯体の η2-side-on 構造の B1 状態では、11 電子が Scheme 2 のように MO を占める。(yz)結合は主に Cr の dyz 軌道と N2 の*z 軌道から形成され、最も安定なので 2 電子占有となるが、Mn の dyz 軌道はこの MO に関与せず非結合性軌道として１電子占有されている(図1)。結果として 1 つの 2 電子占有軌道と、Mn に 5 つ、Cr に 4 つの計 9 つの 1 電子占有軌道が現れ、全体として 10 重項状態をとる。一方 A2 状態では、Cr dyz-N2 *z 結合性 MO の代わりに、Mn の dyz 軌道と N2 の*z 軌道から形成される結合性 MO が 2 電子占有され、Mn のdxy 軌道と N2 の*x 軌道の δ 結合性 MO も 2 電子占有される。結果として 2 つの 2 電子占有
軌道と、Mn に 3 つ、Cr に 4 つの計 7 つの 1 電子占有軌道が現れるが、架橋 N2 のスピン分極の結果、金属間のスピンが Low spin カップリングした 2 重項状態が安定になる。(Mn,Fe)錯体の η1-end-on 構造の A2 状態では*z 軌道との結合性軌道、非結合性軌道共に 2 つの金属
原子の dyz 軌道を同等に含む。この結果、反結合性軌道が大きな占有数をもつなど、その電子状態は Mn,Fe の等核錯体と類似している。
金属原子‐窒素分子間の結合の性質が似ている Mn,Feの異核錯体では等核錯体に類似した 1 つの電子状態をとるのに対し、(Cr,Mn)、(Cr,Fe)錯体は、そのスピン多重度や結合が等核錯体と大きく異なる２つの電子状態が存在する。それぞれの電子状態は、N2 が片方の金属とより強く結合した状態に相当し Cr 原子との結合が強い B1 状
態がより安定である。これは Cr が Mn,Fe より d 軌道エ
ネルギーが高いためである。
(Cr,Mn) η2-side-on 10B1 状態(1.96)(1.00)(Cr,Mn) η2-side-on 2A2 状態(0.02)(1.50)(0.45)(0.51)(1.24)

0.76)(0.00)(Mn,Fe) η1-end-on 8A2 状態(1.49)(1.00)

図 1 (μ-N2)[M1(AIP)][M2(AIP)]の CASSCF 自然軌道 ( )内は占有数を表す[1] Y. I. Kurokawa, Y. Nakao, S. Sakaki, J. Phys. Chem. A, 2010, 114, 1191, 2012, 116, 2292.
[2] W. H. Monillas et al., J. Am. Chem. Soc. 2007, 129, 8090.
[3] J. M. Smith et al., J. Am. Chem. Soc. 2006, 128, 756.

2E19
ランタノイドキレート触媒を用いる向山アルドール反応の配位構造と立体選択性の理論的研究(京大福井謙一研究セ)○畑中美穂、諸熊奎治

Theoretical study of coordination structure and stereoselectivity on
Mukaiyama-Aldol reaction catalyzed by lanthanide chelate complex(FIFC Kyoto University) ○Miho Hatanaka, Keiji Morokuma【序】ランタノイドを触媒とする水中での C-C 結合生成反応は、環境負荷が低いことから、様々な種類への応用が期待されているが、これまで限られた基質でしか実現しておらず、な
かでも、立体選択的な反応の報告例は、非常に少ない。我々はこれまでに、Eu3+(H2O)8 を用いる向山アルドール反応の機構について調べ、(i)Eu3+がアルデヒドを活性化すること、(ii)反応が C-C 結合生成、プロトン移動、SiMe3 基の脱離の順に段階的に起こること、(iii)水の構造ゆ
らぎのために、構造が少しずつ異なる多数の遷移状態が存在し、立体選択性の定量的な議論のためには、これらの遷移状態をサンプリングする必要がある、ということを明らかにしていた。1,2 そこで本研究では、Eu3+周りに配位子を付けることで、水中で高い立体選択性を実現した Scheme 1 の反応 3 に着目する。

Scheme 1. Eu キレート触媒を用いる向山アルドール反応 3
次に、本反応に用いる配位子について説明する。ランタ
ノイドの水溶性配位子として広く用いられているものに、
DOTA 誘導体がある。この配位子は、カルボン酸エステルの「腕」の部分の回転と、12 員環の異性化により、様々な異Figure 1. Eu キレート触媒性体が存在することが知られている。また、実験から、アルデヒド 1 分子と水 2 分子が Eu3+に配位することが確認されているが、その配位構造は自明ではない。3
そこで、本研究では、まず始めに、Figure 1 の Eu キレート錯体の安定構造を非調和下方歪み追跡(ADDF)法によって調べる。得られた構造の中からボルツマン分布の大きい安定構造を選び出し、それらにアルデヒドと水が配位した構造を、人工力誘起反応(AFIR)法を用いて求める。最後に、立体選択性を決める C-C 結合生成の遷移状態を AFIR 法によって、サンプリ
ングし、立体選択性がどのような機構で発現するかを明らかにする。
【計算方法】ADDF 法による局所安定構造の探索では、Eu3+(イオン半径：1.07Å)を Lu3+(イオン半径：0.97Å)に置換し、PM6 を用いた。また、AFIR 法による近似反応経路探索には、B3LYP

法を用い、Eu3+には 4f 電子を内殻に含めた相対論的 ECP(RECP)及びその基底関数(7s6p5d)/[2s1p1d]を、他の原子には 6-31G を用いた。得られた全ての構造を B3LYP-D3 法で再度構造最適化した。ここでは、Eu3+には RECP 基底関数(7s6p5d)/[5s4p3d]を、O,N には 6-31+G*を、C,Hには 6-31G*を用い、溶媒和効果は PCM(水)で考慮した。最後に、Eu3+には RECP 基底関数
(8s7p56d3f2g)/[6s5p4d3f2g]を、他の原子には cc-pVTZ を基底関数に用いて B3LYP-D3 法で一点計算を行った。
【結果】Eu キレート錯体に ADDF 法を用いることで、26 種類の局所安定構造を得た。その中から、分布の大きい構造に着目した。これらに対してアルデヒド 1 分子と水 2 分子が配位した構造のうち、代表的なものを Figure 2 に示す。注目すべき点は、(a)のような 12 員環が正方形になる異性体だけでなく、(b)のように 12 員環が歪んだ構造を持ち、一辺の原子数が異なる四角形になるような異性体も安定に存在することである。(a)の場合、配位子の部分だけに着目すると C2 構造になっているため、±x 方向からエノールエーテルが近づくことができる
が、アルデヒドが配位子との間の水素結合によって固定されるため、エノールエーテルがア
ルデヒドの反応点に接近できるのは、
フェニル基に阻害されない＋x 方向のみに限定される。
(b)のような歪んだ 12 員環構造を持つ場合は、左右のカルボン酸エステルの「腕」の傾きが異なる。このため、エノールエーテルが接近できるのは、カルボン酸エステルの「腕」に阻害されない＋x 方向のみに限定される。このように、配位子の構
造毎に、立体選択性が発現する
原因が異なることが分かった。
発表では、アルデヒドとエノー
ルエーテルの間の C-C 結合生成
段階の遷移状態をサンプリング
し、配位子の異性体の内、どの構造が、低い反応障壁を持つのかを明らかにすると共に、ジアステレオ選択性とエナンチオ選
択性の両方が発現する理由について議論する。
Figure 2. 安定な配位子の異性体構造： (a,b-1)は z 方向から、(a,b-2)、(a,b-3)はそれぞれ＋x、－x 方向から見た図【参考文献】[1] M. Hatanaka, et al., J. Chem. Theory Comput., 2013, 9, 2882. [2] M. Hatanaka, et al., J. Am.
Chem. Soc., Submitted. [3] Y. Mei, J. Am. Chem. Soc., 2010, 132, 12871.

2E20
Computational Insight into the C-H functionalization Jungle
Djamaladdin G. Musaev
Cherry L. Emerson Center for Scientific Computation, Emory University,
Atlanta, Georgia, U.S.A.
I will present our collaborative approaches to the transition metal catalyzed C-H bondfunctionalization. At first, I will elaborate our efforts on understanding the di-Rh, di-Ru,(pybox)Ru catalyzed C-H bond alkylation and amination reactions. I will analyze thefactors controlling the reactivity of these complexes and make intriguing predictions. I
will apply of those principles to predict the (PDI)Fe-catalyst for the C-H bond alkylation.
In the second part of my presentation, I will demonstrate results from our recently
discovered mono-protected amino acid ligands (MPAA) promoted Pd(II)-catalyzed
enantioselective C–H activation reaction with both pyridine and carboxylic acid directing
groups. Presented computational investigations allowed us to gain insights into themechanisms and important elementary steps of the reaction, nature of active species, aligand coordination mode to the Pd(II), transition state structure of the C–H activation
step and base effect.

	
  

1