Home
4584 words
23 minutes
【計算化学】NDDO近似(Neglect of Diatomic Differential Overlap)の数理的詳細と実装論:現代半経験的分子軌道法の基盤

最終更新:2025-12-31

注意: この記事はAIによって自動生成されたものです。正確な情報については、引用元の原著論文(Pople et al., 1965; Dewar et al., 1977等)をご確認ください。

序論:量子化学における計算コストの壁と近似の必要性#

1951年にRoothaanによって提唱されたLCAO-MO-SCF法(Roothaan-Hall方程式)は、分子の電子状態を非経験的(ab initio)に求めるための基礎方程式を確立した。しかし、この方程式を厳密に解くためには、基底関数数 NN に対して N4N^4 のオーダーで増大する二電子積分 (μνλσ)(\mu\nu|\lambda\sigma) をすべて計算する必要がある。1960年代から70年代の計算機資源においては、これは極めて小規模な分子を除いて実行不可能なタスクであった。

この計算コストの壁を突破するために考案されたのが、「微分重なりの無視(Neglect of Differential Overlap: NDO)」に基づく一連の近似手法である。これらは、値が小さいと予想される積分を計算せずにゼロとみなすことで、計算量を劇的に削減する。

その中でも NDDO (Neglect of Diatomic Differential Overlap) 法は、最も洗練された近似レベルに位置する。CNDO (Complete Neglect of Differential Overlap) や INDO (Intermediate Neglect of Differential Overlap) が、電子反発積分における方向依存性を大きく損なう近似を行っていたのに対し、NDDOは原子内の軌道の異方性(s, p軌道の区別や方向性)を保持しつつ、異なる原子間の微分重なりだけを無視する。この特性により、NDDOは孤立電子対の相互作用や立体的な反発を記述する能力を持ち、後のMNDO、AM1、PM3といった実用的な半経験的手法の直接的な基盤となった。

本稿では、NDDO法の数理的構造を、基礎方程式からプログラム実装レベルまで詳細に解剖する。


1. 理論的枠組み:Roothaan-Hall方程式とNDO近似#

1.1 基礎方程式#

閉殻系のHartree-Fock法において、分子軌道係数行列 C\mathbf{C} を決定するRoothaan-Hall方程式は以下で与えられる。

FC=SCϵ\mathbf{FC} = \mathbf{SC\epsilon}

ここで、F\mathbf{F} はFock行列、S\mathbf{S} は重なり行列、ϵ\mathbf{\epsilon} は軌道エネルギーの対角行列である。Fock行列の要素 FμνF_{\mu\nu} は次のように定義される。

Fμν=Hμν+λσPλσ[(μνλσ)12(μλνσ)]F_{\mu\nu} = H_{\mu\nu} + \sum_{\lambda}\sum_{\sigma} P_{\lambda\sigma} \left[ (\mu\nu|\lambda\sigma) - \frac{1}{2}(\mu\lambda|\nu\sigma) \right]
  • HμνH_{\mu\nu}: 一電子積分(運動エネルギー+核引力ポテンシャル)
  • PλσP_{\lambda\sigma}: 密度行列 (Pλσ=2ioccCλiCσiP_{\lambda\sigma} = 2 \sum_{i}^{occ} C_{\lambda i} C_{\sigma i})
  • (μνλσ)(\mu\nu|\lambda\sigma): 二電子積分(電子間反発)
(μνλσ)=ϕμ(1)ϕν(1)1r12ϕλ(2)ϕσ(2)dτ1dτ2(\mu\nu|\lambda\sigma) = \iint \phi_\mu^*(1) \phi_\nu(1) \frac{1}{r_{12}} \phi_\lambda^*(2) \phi_\sigma(2) d\tau_1 d\tau_2

1.2 NDO近似のクラス#

NDO近似とは、異なる原子にある基底関数 ϕμ\phi_\mu (原子A) と ϕν\phi_\nu (原子B, ABA \neq B) の積(微分重なり)を体積要素 dτd\tau 全域でゼロとみなすことである。

ϕμ(1)ϕν(1)dτ10(if AB)\phi_\mu(1)\phi_\nu(1) d\tau_1 \approx 0 \quad (\text{if } A \neq B)

これを二電子積分に適用するレベルによって手法が分類される。

  1. CNDO: すべての微分重なりを無視。(μνλσ)(\mu\nu|\lambda\sigma)μ=ν\mu=\nu かつ λ=σ\lambda=\sigma の場合のみ残り、かつ方向性を無視して (μμλλ)=γAB(\mu\mu|\lambda\lambda) = \gamma_{AB} とする。
  2. INDO: 一中心積分(原子A上の μ,ν,λ,σ\mu, \nu, \lambda, \sigma)については交換積分などを考慮するが、二中心積分はCNDOと同様。
  3. NDDO: 異なる原子間の微分重なりのみを無視する。同一原子上の微分重なり(例:原子A上の ss 軌道と pxp_x 軌道の積)は保持する。

2. NDDO近似の数学的定式化#

NDDO法において、二電子積分 (μνλσ)(\mu\nu|\lambda\sigma) が非ゼロとなる条件は以下の通りである。

μA,νA,λB,σB\mu \in A, \nu \in A, \lambda \in B, \sigma \in B

すなわち、電子1は原子A上の分布 ϕμϕν\phi_\mu \phi_\nu に属し、電子2は原子B上の分布 ϕλϕσ\phi_\lambda \phi_\sigma に属する場合である。ここで A=BA=B の場合(一中心積分)と ABA \neq B の場合(二中心積分)が含まれる。

一方、ハイブリッドな積分、例えば (μAνBλCσD)(\mu_A \nu_B | \lambda_C \sigma_D) のような項はすべてゼロとなる。これにより、計算すべき二電子積分の数は N4N^4 から N2M2N^2 M^2MMは原子あたりの基底関数数)程度まで劇的に削減される。

2.1 回転不変性(Rotational Invariance)の要請#

Popleらが1965年の論文で強調した重要な点は、**「座標軸の回転に対して、エネルギーや電子密度などの物理量が不変でなければならない」**という条件である。 一般に、原子 pp 軌道などは座標軸に依存して定義される(px,py,pzp_x, p_y, p_z)。局所的な座標軸を回転させたとき、基底関数は線形変換される。NDDO近似を行った結果の積分値が、この変換に対して整合性が取れていない場合、計算結果は分子の向きに依存してしまう(空間対称性の破れ)。

これを満たすために、NDDO法では以下の制約が課される。

  1. 一中心二電子積分: 原子内の積分は、純粋な原子スペクトルデータから決定されるか、あるいはSlater-Condonパラメータ(F0,F2,G1F^0, F^2, G^1 等)を用いて解析的に表現されなければならない。
  2. 二中心二電子積分: 原子Aと原子Bの間の相互作用は、それぞれの原子上の局所座標系における多重極子間相互作用として記述され、かつその相互作用テンソルは回転に対して適切に変換されなければならない。

3. Fock行列の詳細な導出#

NDDO近似下でのFock行列要素 FμνF_{\mu\nu} を、対角ブロック(同一原子内)と非対角ブロック(原子間)に分けて導出する。

3.1 同一原子内要素 (μA,νA\mu \in A, \nu \in A)#

Fμν=Hμν+λ,σAPλσ[(μνλσ)12(μλνσ)]+BAλ,σBPλσ(μνλσ)F_{\mu\nu} = H_{\mu\nu} + \sum_{\lambda, \sigma \in A} P_{\lambda\sigma} \left[ (\mu\nu|\lambda\sigma) - \frac{1}{2}(\mu\lambda|\nu\sigma) \right] + \sum_{B \neq A} \sum_{\lambda, \sigma \in B} P_{\lambda\sigma} (\mu\nu|\lambda\sigma)

この式は3つの部分から成る。

  1. 一電子項 HμνH_{\mu\nu}:

    Hμν=UμνBAVμν,BH_{\mu\nu} = U_{\mu\nu} - \sum_{B \neq A} V_{\mu\nu, B}
    • UμνU_{\mu\nu}: 原子AのCore Hamiltonian(運動エネルギー + 原子Aの核引力)。μ=ν\mu=\nu の場合は原子価状態イオン化エネルギーなどに相当。
    • Vμν,BV_{\mu\nu, B}: 原子Bの核(コア)による引力ポテンシャル積分 (μνZB/r1B)(\mu\nu | Z_B/r_{1B})。NDDOでは、これを核電荷 ZBZ_B と電子分布 ϕμϕν\phi_\mu\phi_\nu の静電相互作用として扱う。しばしば、核をs軌道電子雲とみなして二電子積分で近似する技法が使われる(Vμν,BZB(μνsBsB)V_{\mu\nu, B} \approx Z_B (\mu\nu | s_B s_B))。
  2. 一中心電子反発: λ,σAPλσ[(μνλσ)12(μλνσ)]\sum_{\lambda, \sigma \in A} P_{\lambda\sigma} [ (\mu\nu|\lambda\sigma) - \frac{1}{2}(\mu\lambda|\nu\sigma) ]

    • これは原子A内部の電子間相互作用である。これらは実験値(原子スペクトル)からパラメタライズされるか、解析的に計算される。
  3. 二中心クーロン反発: BAλ,σBPλσ(μνλσ)\sum_{B \neq A} \sum_{\lambda, \sigma \in B} P_{\lambda\sigma} (\mu\nu|\lambda\sigma)

    • 原子A上の電子分布 ϕμϕν\phi_\mu\phi_\nu と、原子B上の電子分布 ϕλϕσ\phi_\lambda\phi_\sigma の間のクーロン相互作用。NDDOの最大の特徴であり、分布が球対称(s-s)だけでなく、双極子(s-p)や四重極子(p-p)の分布を持つ場合も考慮される。

3.2 異なる原子間要素 (μA,λB,AB\mu \in A, \lambda \in B, A \neq B)#

NDDO近似を適用すると、二電子積分の交換項のみが残る(クーロン項は μA,λB\mu \in A, \lambda \in B なので (μλ)(\mu\lambda | \dots) となり、微分重なり無視により消える… というのは早計で、正しくは以下の通りである)。

元のFock行列の式:

Fμλ=Hμλ+ρσPρσ[(μλρσ)12(μρλσ)]F_{\mu\lambda} = H_{\mu\lambda} + \sum_{\rho\sigma} P_{\rho\sigma} \left[ (\mu\lambda|\rho\sigma) - \frac{1}{2}(\mu\rho|\lambda\sigma) \right]

NDDO条件 ϕaϕb0(aA,bB)\phi_a \phi_b \approx 0 \, (a \in A, b \in B) を適用する。

  • クーロン積分項 (μλρσ)(\mu\lambda|\rho\sigma): μA,λB\mu \in A, \lambda \in B なので、被積分関数内の ϕμϕλ\phi_\mu \phi_\lambda がゼロとなる。よってこの項は消滅する。
  • 交換積分項 (μρλσ)(\mu\rho|\lambda\sigma): これが非ゼロになるには、μ\muρ\rho が同一原子、λ\lambdaσ\sigma が同一原子でなければならない。すなわち ρA,σB\rho \in A, \sigma \in B である必要がある。

よって、

Fμλ=Hμλ12ρAσBPρσ(μρλσ)F_{\mu\lambda} = H_{\mu\lambda} - \frac{1}{2} \sum_{\rho \in A} \sum_{\sigma \in B} P_{\rho\sigma} (\mu\rho|\lambda\sigma)
  1. 共鳴積分 Hμλ=βμλH_{\mu\lambda} = \beta_{\mu\lambda}:

    • 電子が原子Aから原子Bへ移動するエネルギーを表す。通常、重なり積分 SμλS_{\mu\lambda} に比例する形で経験的に定められる。
    • DewarのMNDO/AM1では βμλ=12(βA+βB)Sμλ\beta_{\mu\lambda} = \frac{1}{2}(\beta_A + \beta_B) S_{\mu\lambda} とされる。
  2. 交換相互作用:

    • 12ρAσBPρσ(μρλσ)-\frac{1}{2} \sum_{\rho \in A} \sum_{\sigma \in B} P_{\rho\sigma} (\mu\rho|\lambda\sigma)
    • 密度行列の非対角要素 PρσP_{\rho\sigma}(結合次数)と、二中心二電子積分の積で表される。

4. 実装のためのアルゴリズム詳細#

NDDO法を実装する際の最大の山場は、二中心二電子積分 (μνλσ)(\mu\nu|\lambda\sigma) の計算と、その回転操作である。

4.1 局所座標系と多重極展開#

原子Aと原子Bの間の積分を計算する場合、一般のグローバル座標系 (X,Y,Z)(X, Y, Z) で直接計算するのは複雑である。そこで、以下の手順を採る。

  1. 局所座標系の定義: 原子Aから原子Bに向かうベクトルを zz' 軸とし、適切な x,yx', y' 軸を定義する局所座標系を設定する。
  2. 基底関数の回転: グローバル座標系の基底関数(px,py,pzp_x, p_y, p_z など)を、局所座標系の基底関数(pσ,pπp_\sigma, p_\pi など)の線形結合で表す。これには回転行列 R\mathbf{R} を用いる。
  3. 局所系での積分評価: 局所系では、原子AとBが zz' 軸上に並ぶため、対称性により多くの積分がゼロになる。非ゼロの積分のみを計算する。
  4. 積分の逆変換: 計算された局所積分を、回転行列を用いてグローバル座標系に戻すか、あるいは密度行列の方を局所系に変換してFock行列を構築し、最後にFock行列を逆変換する。一般的には後者が計算効率が良い。

4.2 二中心二電子積分の近似式#

NDDOでは、二中心積分を点電荷同士の相互作用として近似することが多い。ただし、単なる点電荷ではなく、軌道分布 ΩA=ϕμϕν\Omega_A = \phi_\mu \phi_\nu を表現する「多重極配置」として扱う。

DewarのMNDO/AM1では、原子A上の分布 μν\mu\nu (例:sspzp_zの積は双極子的な分布を持つ)を、微小距離離れた点電荷のセットで表現する手法(Klopman-Ohno式の拡張)を用いる。

基本的な積分近似式(Klopman-Ohno):

(μνλσ)1(RAB+δμν+δλσ)2+(multipole corrections)(\mu\nu|\lambda\sigma) \approx \frac{1}{\sqrt{ (R_{AB} + \delta_{\mu\nu} + \delta_{\lambda\sigma})^2 + (\text{multipole corrections}) }}

ここで RABR_{AB} は原子間距離、δ\delta は原子の大きさに関連するパラメータである。より厳密には、多重極展開された各構成要素ごとの距離を用いて総和を取る。

4.3 コア間反発エネルギー#

全エネルギー EtotE_{tot} は電子エネルギー EelecE_{elec} とコア間反発 EcoreE_{core} の和である。

Ecore=A<BZAZBγAB(RAB)×(Screening Factor)E_{core} = \sum_{A < B} Z_A Z_B \gamma_{AB}(R_{AB}) \times (\text{Screening Factor})

NDDOの枠組み自体は電子部分の近似であるが、実用的な化学精度を出すためには、このコア反発項のパラメータ化が極めて重要である。AM1ではここにガウス関数を加えている。


5. Pythonによるプログラム実装例#

以下に、NDDO法の核心部分であるFock行列構築の擬似コード(Python/NumPy形式)を示す。ここでは可読性を重視し、回転操作や対称性は簡略化して概念を示す。

import numpy as np

class NDDO_Solver:
    def __init__(self, molecule, basis_set, params):
        self.mol = molecule
        self.basis = basis_set
        self.params = params
        self.num_orbs = basis_set.num_orbitals
        
    def get_rotation_matrix(self, atom_A, atom_B):
        """
        グローバル座標から、A->Bをz軸とする局所座標への回転行列を計算
        """
        # (実装省略: ベクトル演算により3x3または軌道数に応じた行列を生成)
        return R_mat

    def calc_two_center_integrals_local(self, atom_A, atom_B, R_dist):
        """
        局所座標系における二中心二電子積分 (mu nu | lambda sigma) を計算
        NDDO近似により、A上の分布(mu*nu)とB上の分布(lambda*sigma)の相互作用のみ
        """
        integrals = {}
        # A上の軌道ペア (mu, nu)
        for mu in atom_A.orbitals:
            for nu in atom_A.orbitals:
                # B上の軌道ペア (lam, sig)
                for lam in atom_B.orbitals:
                    for sig in atom_B.orbitals:
                        # 多重極近似などを用いて積分値を計算
                        # 例: Klopman-Ohno式
                        val = self.multipole_interaction(mu, nu, lam, sig, R_dist)
                        integrals[(mu, nu, lam, sig)] = val
        return integrals

    def build_fock_matrix(self, P):
        """
        密度行列 P から Fock行列 F を構築
        """
        F = np.zeros((self.num_orbs, self.num_orbs))
        
        # 原子ペアごとのループ
        for i, atom_A in enumerate(self.mol.atoms):
            for j, atom_B in enumerate(self.mol.atoms):
                
                # インデックス範囲の取得
                idx_A = self.basis.get_indices(i) # 例: [0, 1, 2, 3] (s, px, py, pz)
                idx_B = self.basis.get_indices(j)
                
                if i == j: # 同一原子内 (One-center block)
                    # 1. 一電子項 (Core Hamiltonian)
                    H_core = self.get_one_center_H(atom_A, P)
                    
                    # 2. 電子反発項 (Coulomb & Exchange)
                    G_intra = np.zeros((len(idx_A), len(idx_A)))
                    P_sub = P[np.ix_(idx_A, idx_A)]
                    
                    for mu in range(len(idx_A)):
                        for nu in range(len(idx_A)):
                            sum_val = 0.0
                            for lam in range(len(idx_A)):
                                for sig in range(len(idx_A)):
                                    # 一中心二電子積分 (mu nu | lam sig)
                                    g_val = self.get_one_center_2e(atom_A, mu, nu, lam, sig)
                                    # P_lam,sig * [ (mu nu | lam sig) - 0.5 * (mu lam | nu sig) ]
                                    exc_val = self.get_one_center_2e(atom_A, mu, lam, nu, sig)
                                    sum_val += P_sub[lam, sig] * (g_val - 0.5 * exc_val)
                            G_intra[mu, nu] = sum_val
                    
                    # 3. 他核からのポテンシャルと他核上の電子との反発
                    # NDDOではこれらをまとめて扱うことが多い
                    V_eff = np.zeros((len(idx_A), len(idx_A)))
                    for k, atom_C in enumerate(self.mol.atoms):
                        if k == i: continue
                        idx_C = self.basis.get_indices(k)
                        P_C = P[np.ix_(idx_C, idx_C)]
                        
                        # A-C間の二中心積分を取得 (回転考慮)
                        R_AC = np.linalg.norm(atom_A.coord - atom_C.coord)
                        # (実装簡略化: グローバル座標での積分が得られたとする)
                        ERIs = self.get_two_center_ERIs(atom_A, atom_C) 
                        
                        for mu in range(len(idx_A)):
                            for nu in range(len(idx_A)):
                                # 核引力 V_mu,nu,C
                                v_nuc = -atom_C.nuclear_charge * ERIs[(mu, nu, 's_C', 's_C')] # 核をs軌道と近似
                                # 電子反発 sum P_lam,sig (mu nu | lam sig)
                                v_elec = 0.0
                                for lam in range(len(idx_C)):
                                    for sig in range(len(idx_C)):
                                        v_elec += P_C[lam, sig] * ERIs[(mu, nu, lam, sig)]
                                V_eff[mu, nu] += v_nuc + v_elec

                    F[np.ix_(idx_A, idx_A)] += H_core + G_intra + V_eff

                else: # 異なる原子間 (Two-center block)
                    # 1. 共鳴積分 (Beta * Overlap)
                    S_sub = self.basis.overlap[np.ix_(idx_A, idx_B)]
                    beta_A = self.params[atom_A.type]['beta']
                    beta_B = self.params[atom_B.type]['beta']
                    H_res = 0.5 * (beta_A + beta_B) * S_sub
                    
                    # 2. 交換積分項 -0.5 * sum P (mu lam | nu sig)
                    # NDDO近似: (mu nu | lam sig) の形のみ非ゼロ
                    # ここで mu in A, lam in B. 
                    # 式: F_mu,lam += -0.5 * sum_{nu in A} sum_{sig in B} P_nu,sig * (mu nu | lam sig)
                    
                    G_exch = np.zeros((len(idx_A), len(idx_B)))
                    P_cross = P[np.ix_(idx_A, idx_B)] # P_nu,sig
                    
                    ERIs = self.get_two_center_ERIs(atom_A, atom_B)
                    
                    for mu in range(len(idx_A)):
                        for lam in range(len(idx_B)):
                            sum_val = 0.0
                            for nu in range(len(idx_A)):
                                for sig in range(len(idx_B)):
                                    # (mu nu | lam sig)
                                    # mu, nu は A上, lam, sig は B上
                                    val = ERIs[(mu, nu, lam, sig)]
                                    sum_val += P_cross[nu, sig] * val
                            G_exch[mu, lam] = -0.5 * sum_val
                    
                    F[np.ix_(idx_A, idx_B)] += H_res + G_exch

        return F

このコードにおける get_two_center_ERIs は、内部で回転行列を用いて局所座標系での積分をグローバル系に変換する処理を含む。NDDO法の効率性は、この積分計算が解析的な近似式で高速に行える点にある。

6. 歴史的意義と実利的な成果#

6.1 Popleによる定式化からDewarによる実用化へ#

1965年のPopleらの論文は、NDDOの理論的整合性(不変性)を確立したが、実用的なパラメータセットを提供するまでには至らなかった。1970年代、DewarらはこのNDDO近似を採用し、実験値を再現するようにパラメータを徹底的にチューニングする方針(半経験的アプローチ)を採った。これが MNDO (Modified Neglect of Diatomic Overlap) である。MNDOは有機化学反応の遷移状態探索や生成熱の予測において、当時のab initio法を遥かに凌ぐコストパフォーマンスを示し、計算化学を実験化学者のツールへと変貌させた。

6.2 CNDO/INDOに対する優位性#

CNDOやINDOは、孤立電子対の方向性や、原子間の反発における立体的な異方性を記述できないという致命的な欠点があった。NDDOは、同一原子上の s,ps, p 軌道間の分布(双極子)を保持するため、例えば水の孤立電子対同士の反発や、アミド基の平面性などを正しく記述できるポテンシャルを持っていた。6.3 AM1への接続しかし、MNDO(純粋なNDDOのパラメータ化)にも限界があった。特に原子間の反発が強すぎ、水素結合を再現できないという問題である。これを解決するために、DewarらはNDDOの枠組みはそのままに、コア反発項(Core Hamiltonianの一部ではないが、全エネルギー計算に使われる項)にガウス型関数を加える改良を行った。これが AM1 である。つまり、NDDOはAM1の「エンジン(電子状態計算部分)」であり、AM1はそのエンジンに「サスペンション(コア反発補正)」と「チューニング(パラメータ)」を施した完成車であると言える。

7. 結論#

NDDO法は、量子化学の第一原理方程式であるRoothaan-Hall方程式に対し、物理的直観に基づいた大胆かつ合理的な近似(二原子微分重なりの無視)を導入することで、計算量を N4N^4 から N2N^2 オーダーへと削減した。その数理的本質は、分子内の相互作用を「原子内構造」と「原子間多重極相互作用」へと分離することにあり、これにより回転不変性を保ちつつ、化学的に重要な方向性を維持することに成功した。今日、MNDO、AM1、PM6、さらには密度汎関数タイトバインディング(DFTB)法に至るまで、多くの高速計算手法がこのNDDOの系譜に連なっていることは、この近似の物理的妥当性と工学的有用性を如実に物語っている。

参考文献#

  • J. A. Pople, D. P. Santry, and G. A. Segal, “Approximate Self-Consistent Molecular Orbital Theory. I. Invariant Procedures,” J. Chem. Phys. 43, S129 (1965).
  • M. J. S. Dewar and W. Thiel, “Ground states of molecules. 38. The MNDO method. Approximations and parameters,” J. Am. Chem. Soc. 99, 4899 (1977).
【計算化学】NDDO近似(Neglect of Diatomic Differential Overlap)の数理的詳細と実装論:現代半経験的分子軌道法の基盤
https://ss0832.github.io/posts/20251231_compchem_nddo_overview/
Author
ss0832
Published at
2025-12-31