東京大学 生産技術研究所
北條 博彦

化学生命工学専攻

有機機能材料学特論II

分子をモデリングするモデリングとは？
modeling --- 模型を作ること・模型の観察によって実物を理解する
ex. 飛行機と気流、船と波、etc.

分子モデリング≒計算（機）化学・分子模型から分子のはたらきが理解できるのか？
→生体分子システムの機構解明や、有機材料のa priori（予見的）設計が可能に！？

化学の世界のスケール
物理量X

日常世界

分子の世界

電子の世界

E  RT ln p

F  mamv  h

長さ（m）
質量（kg）

時間（s）
速度（m s-1）
運動量（kg m s-1）

エネルギー（J）
log10X50-5-10-15-20-25-30-35記号からモデルへ
紀元前のギリシャ：Democritosの原子論→物質の性質の違いはatomosの形の違いが生む

17-18世紀：Daltonの原子説→化学の要素を記号で表した

水素

窒素

リン硫黄

炭素

酸素マグネシウム 石灰

OH

19世紀以降：分子の概念→共有結合を実線で表した
元素周期表の確立→元素記号

20世紀：分子の立体構造→CPK (Corey-PaulingKultun)モデル
20世紀後半～：巨大な分子の構造モデル→コンピュータグラフィックス

NHNN

HNN

N

NO

モデリングの方法論
実験的方法

理論的方法

反応生成物や活性分子力学法（静的構造）

X線回折法（結合距離、分子配列）

分子動力学法（動的構造）

核磁気共鳴法（電子状態、原子間距離）

モンテカルロ法（動的構造）

紫外･可視吸光分光法（電子遷移）

分子軌道法（電子状態）

赤外吸光分光法（分子の振動）

密度汎関数法（電子状態）

マイクロ波吸収分光法（分子の回転）

QM/MM（巨大系）

実験値をもっともよく説明できるようにモデルを組み立てる（帰納的）一定の計算条件下での正確な解を与える（演繹的）

結果は実験値を解釈するための理論にもとづいている（規約主義）

現実の単純化によって逆説的に主要因子を探る

理論計算は正しい解を与えるかよくある誤解：実験結果は事実で、計算結果は架空である
よくある疑問：計算結果は実験値を再現できるのか？

計算化学の捉え方（ゴールの設定）
第1段階：現象の理解
第2段階：実験値の再現、説明
第3段階：実験値の補間、予測

理論は現象のモデル化である問題を単純化しよう
計算条件を設定しよう
モデルの限界を知っておこう応用例：視物質の機能解明
信号伝達タンパクを刺激

視神経へ

光刺激

N

H

Lys296

11-cis-retinal（Lys残基とのSchiff 塩基）

光受容タンパク“ロドプシン”

桿体錐体(rhodopsin)(cone opsin)
感度 高い
低い
応答 遅い
速い
λ / nm500400ロドプシンの結晶構造（二量体）

レチナールの周辺構造
バクテリオロドプシンの光反応中心●単結合のねじれ（コンホメーションの変化）
2~3000 cm-1

●対イオンとの相互作用（静電力，水素結合）4~5000 cm-1

●芳香族残基との相互作用（分散力） ~1000 cm-1

●極性基との相互作用（配向力，誘起力）~1000 cm-1

Figure cited from: Gao, et al., J. Comput. Chem. 23, 96-105 (2002)

分子間にはたらく力
相互作用の種類

力の起源

特徴

静電力

電荷―電荷，永久双極子

長距離まではたらく

永久双極子―永久双極子

距離の－6乗で減衰し、近距離でのみ重要な寄与をする
無極性分子では分散力が支配的

配向力
van derWaals力 誘起力

電荷―誘起双極子

分散力

誘起双極子―誘起双極子

電荷移動相互作用

電位供与体―電子受容体

電荷移動吸収による
着色

水素結合

OH, NHなどの官能基

特異的，指向性強い疎水性相互作用

極性溶媒中の無極性基

エントロピー支配オプシンシフト問題
発色団の吸収波長を決める様々な要因
分子内：（１）コンホメーション変化（～101 kJ mol-1）
分子間：（２）周辺残基との相互作用(a) 静電力（～102 kJ mol-1）(b) vdW力（～100 kJ mol-1）(c) 電荷移動力（～100 kJ mol-1）(d) 水素結合（～101 kJ mol-1）（３）非特異的相互作用(a) 溶媒効果（～101 kJ mol-1）(b) 疎溶媒効果（～101 kJ mol-1）

最近の論文から

Figure 1 Graphical representation of thechromophore and its surrounding residues in bR.
The residues and water molecules labeled withthe blue letters are included as the QM region in
the ONIOM geometry optimization for model 3.

Figure 3↑ Contribution of each residue or water molecule tothe opsin shift. Each ΔE value was obtained by subtracting theexcitation energy for the wild-type bR from that for each “Switchedoff” mutant. The figure shows only the cases for which the absolutevalues of ΔE were greater than 0.01 eV.
Figure 4→ The charge distribution along the retinal polyene chain.
Each value of charge corresponds to the sum of the charge on thecarbon atom of interest and those on the attached H atoms. (a) The
data for the gas-phase PRSB, whose geometry was optimized atthe B3LYP/6-31G(d,p) level. (b) The data for the supermolecularcomplex composed of PRSB and the counterion residues (Asp85,Asp212 and W402). Its geometry was taken from the optimized bR(model 3). (c) The data for PRSB in bR (model 3).

Figure 7 ↑ Schematic representation of the molecular orbital
diagram for the excited state responsible for the main
absorption band (λ max ). (a) the wild-type bR, (b) W138G andW189G, and (c) Y185G. In all the cases, LUMO was localized
only on the chromophore PRSB (see Figure 4(b)).
←Figure 6 (a) Graphical representation of the electron
density difference between the ground state and the excitedstate responsible for the main absorption band (λ max ): an
increase in the electron density upon excitation is shownby the yellow lobes and a decrease by the green lobes.
(b) The occupied and unoccupied MOs to form the singly
excited configurations whose absolute values of the CI
coefficients are larger than 0.3. The data were taken fromthe calculated results for model 3.

どんな方法で何を計算するか
計算手法
分子力学法

基礎となる方程式
Newton方程式

分子動力学法

半経験的分子軌道法

最安定構造
溶媒効果

Hartree-FockRoothaan方程式

ab initio分子軌道法
密度汎関数法

計算結果

基準振動
反応経路

Kohn-Sham方程式

励起状態摂動法

磁気遮蔽定数

配置間相互作用

静電ポテンシャル構造最適化●非効率的な方法
二原子分子の場合

RE（R)

R0R

●効率的な方法R

量子化学の基礎
光の伝搬
Maxwell方程式より
  ckck

2π 2πPlanckの量子仮説よりE  h  
Einsteinの相対性理論より
Ec m c  p
2 22if m  0
 E  cp
E h hp   kcc


波長



波数k

振動数

2π


角振動数   2π

粒子と波動の二重性～量子力学の公理～

速度（群速度）

運動量

エネルギー粒子

波動pm


kp

kp22m波動を表す式
波動方程式

2k 2 2
 ( x, t )  2 2 2  ( x, t )2 2 x t
一般解 ( x, t )   0 {cos( kx  t )  i sin( kx  t )}  0 exp( ikx) exp( it )変数分離した方程式

22

(x,t)

k
 ( x, t )2 2
 x
運動量の演算子表現i 

 pˆ
 x

 p  k 

運動量の固有状態と固有値pˆ  ( x, t ) i  0 exp( ikx) exp( it )x k 0 exp( ikx) exp( it ) k ( x, t )波の重ね合わせと期待値
波数kの正弦波 k ( x, t )   k exp( ikx) exp( it )二種類のkの波の重ね合わせ ( x)   1 exp( ik1 x)  2 exp( ik 2 x)

これは、状態Ψ(x)において運動量がki となる確率が

i22

1   22

Ψ(x)は運動量の固有状態にはなっていない であることを意味している。ここpˆ  ( x)  k1 1 exp( ik1 x)  k 2 2 exp( ik 2 x) で、

 (k1  k 2 ) ( x)運動量の期待値（平均値）の定義

pˆ * ( x) pˆ ( x) dx*
 ( x)( x) dx2k1  1  k 2  22

1   22222

1   2  1となるようにΨ1, Ψ2を選んでおく（規格化する）と式が簡単になる。
また、k1 ≠ k2ならば、

 exp( ik x) exp(ik x) dx  δijである（規格直交化）。

ij

波束と不確定性原理
波数kの正弦波 k ( x, t )   k exp( ikx) exp( it )さまざまなkの波の重ね合わせ ( x, t )   k exp( ikx) dk exp( it )kの重みが正規分布だったら・・・

k2

 (k  k0 ) 2 1exp  
3/ 22(2π) s2sxの重みも正規分布になる

p, x の不確定さ（～標準偏差） ( x, t )   k exp( ikx) dk exp( it )22s22exp2s(xx)0(2 π)1/ 2

2

（Fourier変換）

p ~ s , x ~
不確定性原理 xp ~1212s円環状に拘束された波～ボーアモデル～
中心力のポテンシャル

Ze 2U (r )  4π 0 r

中心力と遠心力のつりあい

p2Ze 2
mr 4π 0 r 2mZe2r
4 π 0 p 2
角運動量の量子化（ボーアの量子条件）

rn p  n

軌道半径は量子化される

h 2 0rn  nπmZe2np22πrn  nh p  

力を受けて運動する波
力学的全エネルギーp2E U (r )2m
全エネルギーの演算子表現
2pˆ1ˆH U (rˆ)  
 2  U (rˆ)
2m2m
 2 22
2    2  2  2 xyz 
ポテンシャルを感じながら運動する波動が満たす式（Schrödinger方程式）

Hˆ  ( x)  E ( x)水素原子1s軌道の解

1s (r ) 1201exp( r )
π

100
0.05

800.2

600.1

40
2000

20

406080100

シュレーディンガー方程式
水素原子の場合

ヘリウム原子の場合11Hˆ    2 2r

1 22Hˆ   1   22
2 2
1     r1 r2  r1  r21
 hˆ(r1 )  hˆ(r2 ) r121exp  r ,π12  r exp   r ,
4 2π 2

二電子演算子
解析的に解ける例はわずか一電子演算子近似を用いなくてはならない

電子の存在を確率で表す
2.0y / a00.0-1.0-2.0-3.0-3.0

N = 100-2.0-1.0

0.0

1.0

3.0
Bohr半径

2.0

2.0

1.0

1.0

0.0

0.0

-1.0-1.0-2.0-2.0-3.0-3.0-2.0
2.0
3.0

N = 500-1.002N = 5000

-1.0

0.0

1.0

2.0

3.0

x / a01.00.80.8f (r )   4πr 2  (r ) dr -2.0

0.60.6

0.4

0.4

0.2

0.2

0.0

0.00123r / a0456f (r)r-3.0
3.0 -3.01.02

動径分布の累積2.024 πr 2  (r )1.0x / a0x / a0

動径分布関数0.0

vdW半径4r  / a0-1

y / a0

1.0

3.0y / a03.0

電子雲の広がり
H原子の電子軌道1s

Na原子の電子軌道

3s
1s

2s

2p5Å波動関数と物理量物理量を表す演算子Wに対し，期待値は量子化学における波動関数→主に定在波を扱う+

W －*ˆ
Ψ WΨ dv*Ψ Ψ dvex. 双極子モーメント m
軌道角運動量 L2Lのz成分 Lzスピン角運動量 S2Sのz成分 Sz：
Schrödinger方程式の解はエネルギーの固有状態E* ˆΨ H Ψ dv*Ψ Ψ dv