【DFT】B2PLYP汎関数の数理と歴史：摂動論的第2次相関を取り入れた二重混成汎関数の先駆的展開

最終更新：2025-12-30

注意: この記事はAIによって自動生成されたものです。正確な情報については、必ず引用元の原著論文をご確認ください。

序論：DFTの限界突破と「第5の階層」への挑戦#

密度汎関数法（Density Functional Theory: DFT）は、現代の量子化学計算において最も成功した電子状態計算手法の一つである。特に、Hartree-Fock（HF）交換項を一定割合で混合したハイブリッド汎関数（B3LYP、PBE0など）の登場により、DFTは平均的な化学的精度において従来の波動関数理論を凌駕するコストパフォーマンスを実現した。 John Perdewは、DFTの近似レベルの向上を旧約聖書になぞらえて「ヤコブの梯子（Jacob’s Ladder）」と表現した。

局所密度近似 (LDA): 電子密度 $\rho$ のみに依存。
一般化勾配近似 (GGA): 密度勾配 $\nabla \rho$ を導入。
Meta-GGA: 運動エネルギー密度 $\tau$ またはラプラシアン $\nabla^2 \rho$ を導入。
ハイブリッド (Hyper-GGA): 占有軌道を用いたHF交換を導入。

この階層構造において、第4段のハイブリッド汎関数までは、電子相関（electron correlation）に関しては依然としてDFTの相関汎関数（半局所的な近似）に依存していた。DFTの相関汎関数は、動的相関（dynamic correlation）を効率的に取り込む一方で、長距離分散力やファンデルワールス相互作用のような非局所的な相関効果の記述には原理的な限界があった。また、自己相互作用誤差（Self-Interaction Error: SIE）の完全な除去も困難であり、反応障壁や電荷移動励起の記述に課題を残していた。

2006年、ミュンスター大学のStefan Grimmeは、この壁を破るために「ヤコブの梯子の第5段」に相当する新しいクラスの汎関数、すなわち二重混成汎関数（Double Hybrid Density Functional: DHDF）を提案した。その最初にして最も象徴的なモデルがB2PLYPである。 B2PLYPは、Kohn-Sham軌道を用いた**第2次Møller-Plesset摂動論（MP2）**による相関エネルギーを、DFTの相関エネルギーの一部と置き換えるという大胆な発想に基づいている。これにより、波動関数理論（WFT）の系統的な精度向上能力と、DFTの効率性を融合させることに成功した。

本稿では、Grimmeの原著論文 “Semiempirical hybrid density functional with perturbative second-order correlation” [1] に基づき、B2PLYPの数理的導出、断熱接続公式に基づく理論的背景、および化学計算における実利的な成果について、一万字程度の規模で詳細に解説する。

1. 数理的背景：断熱接続と摂動論的相関の融合#

B2PLYPの設計思想は、単なる経験的な「重ね合わせ」ではなく、DFTの厳密な理論的枠組みである**断熱接続（Adiabatic Connection: AC）**の拡張に基づいている。

1.1 断熱接続公式（Adiabatic Connection Formula）#

DFTにおける交換相関エネルギー $E_{xc}$ は、断熱接続公式を用いて以下のように厳密に表現される。

E_{xc} = \int_0^1 U_{xc,\lambda} \, d\lambda

ここで、 $\lambda$ は結合定数（coupling constant）であり、 $\lambda=0$ は相互作用のないKohn-Sham参照系（非相互作用系）、 $\lambda=1$ は現実の相互作用系（完全な物理系）に対応する。 $U_{xc,\lambda}$ は、結合定数が $\lambda$ であるときのポテンシャルエネルギー期待値から、古典的なクーロンエネルギーを引いたものである（交換相関ポテンシャルエネルギー）。

$\lambda=0$ の極限: このとき、系は非相互作用であり、波動関数は単一のスレーター行列式で表される。したがって、 $U_{xc, \lambda=0}$ は純粋なHartree-Fock交換エネルギー ( $E_x^{HF}$ ) に等しい。
$\lambda=1$ の極限: 完全な多体相互作用を含む状態であり、正確な交換相関エネルギーに対応する密度汎関数近似（DFT）で記述されるべき領域である。

1.2 ハイブリッド汎関数の再解釈#

従来のハイブリッド汎関数（例えばPBE0やB3LYP）は、この積分の被積分関数 $U_{xc,\lambda}$ を、両端（ $\lambda=0$ と $\lambda=1$ ）の線形補間、あるいは低次の多項式で近似したものと解釈できる。例えば、最も単純な “Half-and-Half” 汎関数は、以下のような台形公式による近似に対応する。

E_{xc} \approx \frac{1}{2} (U_{xc,0} + U_{xc,1}) = \frac{1}{2} E_x^{HF} + \frac{1}{2} E_{xc}^{DFT}

より一般的なハイブリッド汎関数は、係数 $a_x$ を用いて以下のように書ける。

E_{xc}^{hybrid} = a_x E_x^{HF} + (1 - a_x) E_x^{DFT} + E_c^{DFT}

（注：通常、相関エネルギーは全領域でDFT近似を用いることが多いが、厳密なACの観点からは相関も $\lambda$ 依存性を持つ）

1.3 二重混成（Double Hybrid）への拡張#

Grimmeは、このAC積分において、 $\lambda=0$ 付近での挙動をより正確に記述するために、**Görling-Levyの摂動論（GLPT）**の知見を取り入れた。 GLPTによれば、 $\lambda \to 0$ の極限において、交換相関エネルギーの勾配（初期勾配）は、第2次摂動エネルギーによって記述されるべきである。すなわち、相関エネルギーの展開は $\lambda$ の2次オーダーから始まる。

E_c = \lim_{\lambda \to 1} E_c^\lambda \approx E_c^{GL2}

Kohn-Sham軌道および軌道エネルギーを用いた場合、この $E_c^{GL2}$ は、通常のMøller-Plesset第2次摂動エネルギー（MP2）と形式的に同等となる（ただし、参照軌道がHF軌道ではなくKS軌道である点が異なる）。

Grimmeは、AC被積分関数 $U_{xc,\lambda}$ を $\lambda$ の関数としてより柔軟にモデル化するために、交換項だけでなく相関項についてもHF（摂動論的）成分とDFT成分を混合するというアイデアに至った。これが「二重混成（Double Hybrid）」と呼ばれる所以である。数式としては以下の形式をとる。

E_{xc}^{B2PLYP} = (1 - a_x) E_x^{DFT} + a_x E_x^{HF} + (1 - a_c) E_c^{DFT} + a_c E_c^{PT2}

ここで、各項の意味は以下の通りである。

$E_x^{DFT}$ : DFT交換汎関数（B2PLYPではB88交換）。
$E_x^{HF}$ : Exact Exchange（HF交換）。KS軌道を用いて計算される。
$E_c^{DFT}$ : DFT相関汎関数（B2PLYPではLYP相関）。
$E_c^{PT2}$ : 摂動論的第2次相関エネルギー（MP2様相関）。

1.4 第2次摂動エネルギー ( $E_c^{PT2}$ ) の定義#

B2PLYPにおける $E_c^{PT2}$ 項は、収束したKohn-Sham軌道 $\{\psi_p\}$ とその軌道エネルギー $\{\epsilon_p\}$ を用いて、通常のMP2エネルギー表式により計算される。

E_c^{PT2} = \frac{1}{4} \sum_{i,j}^{occ} \sum_{a,b}^{virt} \frac{|\langle ij || ab \rangle|^2}{\epsilon_i + \epsilon_j - \epsilon_a - \epsilon_b}

ここで、 $\langle ij || ab \rangle$ は二電子積分（反対称化積分）である。 $i, j$ は占有軌道、 $a, b$ は仮想軌道を走る。分母は軌道エネルギー差（HOMO-LUMOギャップに関連）である。重要な点は、この計算が自己無撞着（SCF）ループの外で行われることである。すなわち、まずハイブリッドDFT（ $a_c=0$ の状態）でKS方程式を解き、得られた軌道を使って一度だけPT2エネルギーを計算し、加算する。これを「post-Kohn-Sham」アプローチと呼ぶ。

2. 汎関数の構築とパラメータ決定#

B2PLYPの具体的な関数形とパラメータ決定プロセスについて詳述する。名称の “B2” は “Becke 2-parameter” に由来し、“PLYP” は摂動論（Perturbative）とLYP相関の組み合わせを示唆している。

2.1 関数形：B88交換とLYP相関の採用#

Grimmeは、ベースとなるDFT汎関数として、実績のあるB88（Becke 1988）交換汎関数とLYP（Lee-Yang-Parr）相関汎関数の組み合わせを採用した。これはB3LYPと同じ構成要素であるが、パラメータ数と構造が異なる。

B2PLYPのエネルギー式は、より具体的に書くと以下のようになる。

E_{xc}^{B2PLYP} = (1 - a_x) E_x^{B88} + a_x E_x^{HF} + (1 - a_c) E_c^{LYP} + a_c E_c^{PT2}

2.2 パラメータ $a_x, a_c$ の最適化#

Grimmeは、この2つの混合係数 $a_x$ （HF交換率）と $a_c$ （PT2相関率）を決定するために、広範な化学データベース（G2/97サブセットなど）を用いた半経験的なフィッティングを行った。最適化のターゲットとなったのは、以下の物性値である。

原子化エネルギー（Atomization Energies）: 分子の結合エネルギーの総和。熱化学的精度の基本指標。
イオン化ポテンシャル（IP）および電子親和力（EA）。
プロトン親和力（PA）。

広範な探索の結果、以下の最適値が決定された。

$a_x = 0.53$ (53%)
$a_c = 0.27$ (27%)

このパラメータセットには、従来のハイブリッド汎関数とは異なる際立った特徴がある。まず、HF交換率（53%）が非常に高い。B3LYP（20%）やPBE0（25%）と比較して倍以上の割合である。通常、DFTにおいてHF交換を増やしすぎると、静的相関の記述が不足し、結合エネルギーが過小評価される傾向がある。しかし、B2PLYPでは $E_c^{PT2}$ 項（27%）が導入されているため、これが長距離相関や静的相関の一部を効果的に補い、高いHF交換率による弊害を打ち消しているのである。この「高HF交換＋PT2相関」というバランスこそが、B2PLYPの精度の鍵である。

3. 歴史的背景：MP2とDFTの架け橋#

B2PLYPの登場は、計算化学の歴史において「DFTとWFT（波動関数理論）の本格的な融合」の始まりを告げるものであった。

3.1 従来のDFTの限界#

2000年代前半まで、DFTは基底状態の構造やエネルギーにおいて無類の強さを誇っていたが、いくつかの「弱点」が克服できずにいた。

分散力（Dispersion Force）の欠如: 局所的な相関汎関数は、長距離での電子相関（ $1/R^6$ の引力）を記述できない。
反応障壁の過小評価: 自己相互作用誤差により、遷移状態が過度に安定化される。
Rydberg励起・電荷移動励起の失敗: 交換ポテンシャルの漸近挙動の誤り。

3.2 MP2の復権と限界#

一方、MP2法は、分散力を（定性的に、ある程度定量的に）記述できる最も安価なWFT手法であったが、計算コストが $O(N^5)$ と高く（DFTは $O(N^3)$ 〜 $O(N^4)$ ）、また開殻系や静的相関が強い系（スピン汚染）に弱いという欠点があった。

3.3 Grimmeの洞察#

Stefan Grimmeは、DFTの「静的相関への強さ・低コスト」と、MP2の「動的相関・分散力の記述力」を組み合わせれば、互いの欠点を補完できると考えた。 B2PLYP以前にも、DFTとMP2を組み合わせる試みはあったが、多くは複雑な範囲分離を用いたり、パラメータ数が多かったりした。Grimmeの功績は、「全領域で一定割合を混ぜる」という極めてシンプルなグローバルハイブリッド形式で、実用的な精度が出ることを示した点にある。特に、PT2項を相関エネルギー全体の約1/4に留めることで、MP2の計算コストを正当化しつつ、MP2特有の過剰評価（Basis Set Superposition Errorなど）の影響を緩和することに成功した。

4. 実利的な成果と検証#

原著論文およびその後の研究において示されたB2PLYPの性能は、当時の標準（B3LYP）を大きく上回るものであった。

4.1 熱化学的精度 (G2/97セット)#

原子化エネルギーの平均絶対偏差（MAD）において、B2PLYPは卓越した性能を示す。

B3LYP: ~2-3 kcal/mol
B2PLYP: ~1.0 kcal/mol 以下 (G2セットの多くのサブセットで)

特に、B3LYPが苦手としていたフッ素化合物や、電子豊富な系、ラジカル反応などにおいて、B2PLYPはMP2単体やB3LYPよりも安定して高い精度を与えた。これは、PT2項による動的相関の補正と、高いHF交換による自己相互作用の低減が相乗的に働いた結果である。

4.2 反応障壁 (Barrier Heights)#

反応速度論において、遷移状態のエネルギー精度は極めて重要である。B2PLYPは53%という高いHF交換率を持つため、自己相互作用誤差が大幅に低減されている。

HTBH（水素移動反応障壁）: B3LYPが数kcal/molの誤差を出すのに対し、B2PLYPは 1 kcal/mol 前後の精度を達成する。これは、反応速度定数の予測において桁違いの信頼性をもたらす。

4.3 弱相関相互作用 (Weak Interactions)#

B2PLYPの最も画期的な点は、ファンデルワールス力（分散力）がある程度記述できるようになったことである。

希ガスダイマーやメタンダイマーなどの分散力支配系において、B3LYPは結合を全く再現しない（斥力のみ）。
B2PLYPは、PT2項（非局所相関成分）のおかげで、これらの系においてもポテンシャル井戸（結合状態）を再現する。ただし、MP2自体も分散力の記述には限界があり（特に $\pi$ -スタッキングなどで過大評価する傾向がある）、B2PLYPも完全ではない。後の研究で、これにさらに経験的な分散力補正を加えた B2PLYP-D や B2PLYP-D3 が開発され、これらは現在における非共有結合系の標準的な高精度手法となっている。

4.4 分子構造と振動数#

分子構造（結合長）に関しては、MP2は結合を短く見積もり、DFT（GGA）は長く見積もる傾向がある。B2PLYPはその中間（混合）であるため、実験値との一致が極めて良い。特に、強い相関が関与する多重結合などの記述において優れている。

5. 議論：メリットとデメリットの深掘り#

B2PLYPは万能ではなく、明確なトレードオフが存在する。

5.1 計算コストとスケーリング#

最大のデメリットは計算コストである。

DFT: 計算時間は系サイズ $N$ に対して $O(N^3)$ または $O(N^4)$ でスケールする。
B2PLYP: PT2項の計算に $O(N^5)$ のコストがかかる。また、仮想軌道の情報が必要なため、大きな基底関数系を必要とする。

したがって、数百原子を超えるような巨大なタンパク質やナノ粒子に対して、B2PLYPを日常的に適用するのは困難である。しかし、数十原子程度の中規模分子においては、CCSD(T)などの高精度WFTよりも遥かに高速であり、現実的な時間で「CCSD(T)に近い精度」を得られる選択肢として重宝される。

5.2 基底関数依存性#

MP2項を含むため、B2PLYPは基底関数の大きさに対して敏感である。DFTは比較的小さな基底関数（6-31G*など）でもエネルギーが収束しやすいが、MP2（およびB2PLYP）は相関エネルギーの収束が遅く、Triple-zeta (TZ) や Quadruple-zeta (QZ) クラスの基底関数（def2-TZVPPなど）が必要となる。これが実質的な計算コストをさらに押し上げる要因となる。

5.3 スピン汚染と開殻系#

B2PLYPは非制限Kohn-Sham (UKS) 軌道を参照してPT2計算を行うことができるが、MP2部分はスピン汚染（Spin Contamination）に脆弱である。遷移金属錯体やラジカル対などでスピン汚染が大きい場合、PT2項のエネルギーが発散したり、不物理的な値になったりするリスクがある。このような系では、通常のDFTの方がロバストな場合がある。

6. プログラム出力：B2PLYPの実装概念#

以下に、B2PLYPのエネルギー計算の流れを概念的なPythonコード（Psi4やPySCFなどの量子化学ライブラリを用いた際のロジック）として表現する。これはB2PLYPが「post-Kohn-Sham」法であることを理解するために有用である。

"""
Conceptual implementation of B2PLYP Energy Calculation
This script demonstrates the logic flow, not a standalone executable.
Libraries like Psi4 or PySCF are typically used for actual computation.
"""

class B2PLYP_Functional:
    def __init__(self):
        # Grimme's original parameters (2006)
        self.ax = 0.53  # Fraction of HF Exchange
        self.ac = 0.27  # Fraction of PT2 Correlation
        
        # DFT components weights implies:
        # Weight of DFT Exchange = 1 - ax = 0.47
        # Weight of DFT Correlation = 1 - ac = 0.73
        self.w_dft_x = 1.0 - self.ax
        self.w_dft_c = 1.0 - self.ac

    def compute_energy(self, molecule, basis_set):
        """
        Step-by-step calculation of B2PLYP energy
        """
        print(f"Calculating B2PLYP for {molecule} with {basis_set}")

        # --- Step 1: Self-Consistent Field (SCF) Calculation ---
        # The SCF is performed using a specific hybrid DFT functional.
        # This "parent" functional includes the DFT part and HF exchange,
        # but NOT the PT2 part. It defines the Kohn-Sham orbitals.
        #
        # Functional for SCF: 
        # E_SCF = 0.53 * E_HF + 0.47 * E_B88 + 0.73 * E_LYP
        # Note: Standard B2PLYP SCF does not include PT2.
        
        scf_result = self.run_scf(molecule, basis_set, 
                                  hf_fraction=self.ax, 
                                  dft_x_weight=self.w_dft_x, 
                                  dft_c_weight=self.w_dft_c)
        
        orbitals = scf_result.orbitals
        orbital_energies = scf_result.orbital_energies
        e_scf_converged = scf_result.total_energy

        print(f"SCF Converged Energy: {e_scf_converged:.8f} Ha")

        # --- Step 2: Perturbative Correction (Post-SCF) ---
        # Calculate MP2-like correlation energy using the converged KS orbitals.
        # Standard MP2 formula:
        # E(2) = (1/4) * sum_{ijab} |<ij||ab>|^2 / (e_i + e_j - e_a - e_b)
        
        e_pt2 = self.calculate_mp2_energy(orbitals, orbital_energies)
        
        print(f"Raw PT2 Correlation Energy: {e_pt2:.8f} Ha")

        # --- Step 3: Final Energy Assembly ---
        # Combine the SCF energy with the scaled PT2 energy.
        # B2PLYP Energy = E_SCF + ac * E_PT2
        
        e_b2plyp = e_scf_converged + self.ac * e_pt2

        return e_b2plyp

    def run_scf(self, mol, basis, hf_fraction, dft_x_weight, dft_c_weight):
        # Placeholder for SCF procedure calling external QC codes
        # In Psi4, this would be defining a custom functional or using 'B2PLYP' keyword
        pass

    def calculate_mp2_energy(self, orbitals, eps):
        # Placeholder for integral transformation and MP2 sum
        # Computational cost is O(N^5)
        # eps: orbital energies (eigenvalues of the KS matrix)
        pass

# Example Usage
# calc = B2PLYP_Functional()
# energy = calc.compute_energy("H2O", "def2-TZVPP")
# print(f"Final B2PLYP Energy: {energy:.8f} Ha")

コード解説#

SCFステップ: まず、ハイブリッド汎関数としてSCF計算を行う。この時の汎関数は、B2PLYPのDFT成分（B88交換 47%、LYP相関 73%）とHF交換（53%）のみを含む。PT2項は軌道や密度に影響を与えない（非自己無撞着）。

PT2ステップ: 収束したKohn-Sham軌道を用いて、MP2エネルギー $E^{(2)}$ を計算する。ここが計算のボトルネックとなる。

エネルギー合成: 最後に、SCFエネルギーに重み付けされたPT2エネルギー（ $0.27 \times E^{(2)}$ ）を加算して、最終的なB2PLYPエネルギーを得る。

結論#

Stefan Grimmeによって提案されたB2PLYPは、密度汎関数法と波動関数理論の境界を取り払い、両者の利点を融合させたマイルストーン的な汎関数である。断熱接続公式に基づく摂動論的展開という堅固な理論的背景と、 $a_x=0.53, a_c=0.27$ という絶妙なパラメータ設定により、B2PLYPは以下の成果を達成した。

熱化学精度の飛躍的向上: G2セットにおいて化学的精度（~1 kcal/mol）を達成。

反応速度論への適用: 高いHF交換率による反応障壁の正確な記述。

非局所相関の取り込み: PT2項による分散力の定性的な記述（後にB2PLYP-Dで定量化）。

計算コスト（ $O(N^5)$ ）と基底関数依存性という課題はあるものの、B2PLYPは「コストをかけてでもB3LYPより良い結果が欲しいが、CCSD(T)には手が届かない」という多くの研究者にとって、実用的かつ高精度な解を提供し続けている。また、その後に続く多数の二重混成汎関数（mPW2PLYP, DSD-PBEP86など）の基礎となり、現代計算化学における「第5の階層」を切り拓いた功績は計り知れない。

##参考文献

S. Grimme, “Semiempirical hybrid density functional with perturbative second-order correlation”, J. Chem. Phys. 124, 034108 (2006).