CI-NEB法: ミニマムエネルギーパスと鞍点探索における理論的進歩とアルゴリズムの詳説

[cite_start]理論化学および凝縮系物理学において、化学反応や拡散現象などの「遷移速度（Transition Rates）」を算出することは中心的な課題の一つです。原子の振動スケール（フェムト秒オーダー）に比べ、これらの遷移現象は通常、何桁も遅い時間スケールで発生します。例えば、活性化エネルギーが 0.5 eV 程度のプロセスであっても、古典的な分子動力学（MD）シミュレーションで遷移を1回観測するためには、現代の計算機を用いても極めて長い計算時間（論文発表当時の試算で $10^4$ 年オーダー）を要する場合があります [cite: 33, 34]。

[cite_start]このような「希少事象（Rare Event）」問題に対処するため、直接的な動力学シミュレーションの代わりに、統計的アプローチである 遷移状態理論（Transition State Theory: TST） が広く用いられています [cite: 35][cite_start]。特に、結晶中や表面上の原子が密に詰まった系においては、調和近似を用いた Harmonic TST (hTST) が有効です。hTSTにおいて、遷移速度 $k_{hTST}$ は以下の式で記述されます [cite: 38, 45]。

k_{hTST} = \frac{\prod_{i}^{3N} \nu_{i}^{init}}{\prod_{i}^{3N-1} \nu_{i}^{\dagger}} \exp\left( -\frac{E^{\dagger} - E^{init}}{k_B T} \right)

[cite_start]ここで、 $E^{\dagger}$ は鞍点（Saddle Point）のエネルギー、 $E^{init}$ は初期状態の極小点のエネルギー、 $\nu_i$ はそれぞれの基準振動数です [cite: 47][cite_start]。この式が示す通り、反応速度の推定において最もクリティカルな要素は、ポテンシャルエネルギー曲面（PES）上の 鞍点（遷移状態） を正確に特定することです [cite: 49]。

[cite_start]初期状態と最終状態を結ぶ経路の中で、統計的に最も支配的な経路は 最小エネルギー経路（Minimum Energy Path: MEP） と呼ばれます。MEP上の極大値が鞍点に対応するため、MEPを探索することは、反応機構の解明と反応速度の算出において不可欠な手続きとなります [cite: 50, 53]。

本稿では、MEP探索の標準的手法である Nudged Elastic Band (NEB) 法と、その精度を飛躍的に向上させる改良手法としてHenkelmanらが2000年に提案した Climbing Image NEB (CI-NEB) 法について、その数学的背景とアルゴリズムの詳細を解説します。

2. Nudged Elastic Band (NEB) 法の理論的枠組み#

[cite_start]NEB法は、初期状態 $R_0$ と最終状態 $R_N$ の間に、 $N-1$ 個の中間イメージ（Replicas） $[R_1, R_2, ..., R_{N-1}]$ を配置し、それらをバネで連結した「弾性バンド」として最適化を行う手法です [cite: 99]。

2.1 従来のChain-of-States法の問題点#

[cite_start]NEB以前の手法や、単純な弾性バンド法では、以下の2つの相反する問題に直面していました [cite: 79]。

Corner-cutting（ショートカット問題）: バネ定数が強すぎる場合、バンドがMEPの曲率に追従できず、エネルギー障壁の内側をショートカットしてしまい、鞍点エネルギーを過大評価する。
Sliding-down（滑り落ち問題）: バネ定数が弱すぎる場合、イメージがポテンシャルの勾配に従って初期状態や最終状態のベイスン（窪み）に滑り落ちてしまい、最も重要な鞍点付近の解像度が低下する。

2.2 NEBにおける力の射影（Nudging）#

[cite_start]NEB法は、これらの問題を「力の射影（Force Projection）」、通称 Nudging によって解決しました。各イメージ $i$ に働く全探索力 $F_i$ は、以下の2つの成分の和として定義されます [cite: 101, 102]。

F_i = F_i^s|_{||} - \nabla E(R_i)|_{\perp}

ここで、各項は以下のように定義されます。

垂直方向のポテンシャル力 ( $-\nabla E(R_i)|_{\perp}$ ): [cite_start]真のポテンシャル力 $-\nabla E(R_i)$ のうち、パスの接線ベクトル $\hat{\tau}_i$ に垂直な成分のみを抽出したものです [cite: 104]。
$-\nabla E(R_i)|_{\perp} = -\nabla E(R_i) + (\nabla E(R_i) \cdot \hat{\tau}_i) \hat{\tau}_i$
この項により、イメージはMEPに向かって収束します。
平行方向のバネ力 ( $F_i^s|_{||}$ ): [cite_start]隣接するイメージ間のバネ力のうち、接線ベクトル $\hat{\tau}_i$ に平行な成分のみを抽出したものです [cite: 109]。
$F_i^s|_{||} = k (|R_{i+1} - R_{i}| - |R_{i} - R_{i-1}|) \hat{\tau}_i$
この項は、イメージ間の間隔を制御する役割のみを果たし、パスの形状収束には干渉しません。

[cite_start]この構成により、NEB法は「バネの強さがパスの形状に影響を与えず」、かつ「ポテンシャル力がイメージの分布を偏らせない」という特性を獲得しました [cite: 80]。

2.3 従来のNEB法の限界#

[cite_start]NEB法はMEPの形状を特定するには極めて有効ですが、離散的なイメージでパスを表現するため、厳密な鞍点 にイメージが配置される保証はありません [cite: 86]。 [cite_start]通常、鞍点のエネルギーは、最高エネルギーを持つイメージとその隣接イメージを用いた補間（例えば3次多項式補間）によって推定されます [cite: 141][cite_start]。しかし、エネルギー障壁がパス全体の長さに比べて狭い場合、鞍点近傍に配置されるイメージ数が不足し、補間精度が著しく低下するという問題がありました [cite: 89, 116]。

3. Climbing Image NEB (CI-NEB) 法#

[cite_start]Henkelmanらが提案したCI-NEB法は、従来のNEB法にわずかな修正を加えるだけで、計算コストを増大させることなく、厳密に鞍点へ収束する ことを可能にした手法です [cite: 26, 90]。

3.1 アルゴリズムと数理#

CI-NEB法の手順は以下の通りです。

通常のNEB法で数回の反復を行い、大まかなパスを形成する。
[cite_start]最も高いエネルギーを持つイメージ $i_{max}$ を特定する [cite: 149]。
[cite_start]この $i_{max}$ （Climbing Image）に対してのみ、働く力を以下の式に書き換える [cite: 151]。

F_{i_{max}} = -\nabla E(R_{i_{max}}) + 2 \nabla E(R_{i_{max}})|_{||}

[cite_start]この式をベクトル成分で展開すると、その物理的意味が明確になります [cite: 152]。

F_{i_{max}} = \underbrace{-\nabla E(R_{i_{max}})|_{\perp}}_{\text{垂直方向の最小化}} + \underbrace{\nabla E(R_{i_{max}})|_{||}}_{\text{平行方向の最大化}}

通常の勾配降下法では、力はポテンシャルの逆勾配 $-\nabla E$ です。これを接線方向成分と垂直成分に分けると $-\nabla E = -\nabla E|_{||} - \nabla E|_{\perp}$ となります。 CI-NEBの式（Eq. [cite_start]5）は、ポテンシャルの勾配ベクトル $-\nabla E$ に対し、平行成分 $2 \nabla E|_{||}$ を加算することで、平行成分の富豪を反転させています [cite: 154]。

垂直成分: $-\nabla E(R_{i_{max}})|_{\perp}$ はそのまま維持されるため、イメージはMEPに沿ってエネルギーの谷底へ向かおうとします。
平行成分: $\nabla E(R_{i_{max}})|_{||}$ は、ポテンシャルエネルギーを 上昇させる 方向に働きます（通常の最適化とは逆方向）。

3.2 CI-NEBの収束特性#

[cite_start]この力の定義により、Climbing Imageは「バンドに垂直な方向にはエネルギー最小化」を行いながら、「バンドに沿った方向にはエネルギー最大化」を行います [cite: 156]。これはまさに、一次の鞍点の数学的定義（ある一方向に対して極大であり、他の全方向に対して極小である点）と一致します。

[cite_start]重要な点は、Climbing Imageはバネ力の影響を一切受けないということです [cite: 155][cite_start]。そのため、他のイメージがバネによって等間隔（あるいは指定された分布）を保とうとする中で、Climbing Imageだけは自由に動き回り、最終的に厳密な鞍点位置に収束します。計算コストに関しては、全てのイメージを同時に緩和させるため、通常のNEB法と比較して追加の計算量は無視できる範囲（同等あるいは10%以内の差）に収まります [cite: 160, 163]。

4. 可変バネ定数 (Variable Spring Constants)#

[cite_start]さらに本論文では、鞍点付近の解像度を向上させるための「可変バネ定数」スキームも提案されています [cite: 28, 175]。

4.1 均一バネ定数の課題#

[cite_start]通常のNEBでは全てのバネ定数 $k$ を均一に設定しますが、反応経路が長く、かつ遷移領域が狭い場合、限られたイメージ数の多くが反応に関与しない平坦な領域に浪費されてしまいます [cite: 179]。

4.2 エネルギー依存型のバネ定数#

[cite_start]そこで、イメージのエネルギーに応じてバネ定数を線形にスケーリングする手法が導入されました [cite: 184]。

k_i = \begin{cases} k_{max} - \Delta k \left( \frac{E_{max} - E_i}{E_{max} - E_{ref}} \right) & \text{if } E_i > E_{ref} \\ k_{max} - \Delta k & \text{if } E_i \le E_{ref} \end{cases}

[cite_start]ここで、 $E_{max}$ はバンド全体の最大エネルギー、 $E_{ref}$ は参照エネルギー（通常は初期・最終状態のうち高い方のエネルギー）です [cite: 208, 211]。 [cite_start]この設定により、エネルギーの高い（鞍点に近い）イメージ間のバネ定数 $k$ が強くなり、結果として鞍点付近にイメージが集中的に配置されるようになります。Nudgingのおかげで、バネ定数を不均一にしてもMEPへの収束性には悪影響を与えません [cite: 182]。

5. 適用事例と検証#

[cite_start]論文では、密度汎関数法（DFT/PW91, Ultrasoft Pseudopotentials）を用いた2つの表面反応系への適用を通じて、CI-NEB法の有効性が実証されています [cite: 29, 93]。

5.1 $CH_4$ の Ir(111) 表面での解離吸着#

[cite_start]この系は、メタン分子が表面から約 4 Å 離れた位置（反応座標 1.0）から、表面上で H と $CH_3$ フラグメントに解離吸着する（反応座標 0.0）プロセスです [cite: 114, 115]。

特徴: エネルギー障壁がMEP全体の長さに比べて非常に狭い。
[cite_start]通常NEBの結果: 鞍点付近にイメージが配置されず、補間による活性化エネルギーの推定値は過小評価され、大きな不確実性が残りました [cite: 141]。
[cite_start]CI-NEBの結果: Climbing Imageが自動的に鞍点の頂点へ移動しました。これにより得られた活性化エネルギー（約 0.4 eV）は、補間に頼ることなく厳密に決定されました [cite: 165]。
[cite_start]構造的知見: 遷移状態において、メタンに最も近いイリジウム原子が表面から 0.5 Å も引き出されるという、自明ではない構造緩和が確認されました [cite: 168]。

5.2 $H_2$ の Si(100) 表面での解離吸着#

[cite_start]水素分子がシリコン表面のダイマー上で解離する反応です。この系では、実験値と理論値の間に長年の不一致が存在していました [cite: 218]。

可変バネ定数の効果: $H_2$ が表面に接近する初期段階はポテンシャルが平坦であり、解像度を必要としません。可変バネ定数を導入した結果、エネルギーの低い領域のイメージ密度を下げ、代わりにバリア付近にイメージを集中させることに成功しました（Fig. [cite_start]2参照） [cite: 222]。
[cite_start]効率: 均一バネ定数を用いた場合（190回の力評価）と比較して、可変バネ定数を用いた場合（178回）は、同等以上の解像度を維持しつつ、計算コストが増加しない（むしろわずかに減少する）ことが確認されました [cite: 226]。

6. 実装と計算コストに関する考察#

6.1 並列化効率#

[cite_start]NEB法およびCI-NEB法は、各イメージの力計算が独立しているため、並列計算に極めて適しています [cite: 97]。本研究では、ワークステーションクラスタ上で通常のEthernet接続を用いて計算が行われましたが、通信ボトルネックは問題になりませんでした。

6.2 CI-NEBのオーバーヘッド#

CI-NEB法は、通常のNEB法の最適化ループの中に組み込まれます。Climbing Imageに対する特別な力計算（Eq. 5）は、通常の勾配計算（Eq. 2）と計算量は変わらないため、アルゴリズム的なオーバーヘッドは皆無です。 [cite_start]また、鞍点探索のために別のアルゴリズム（例えばHessian行列を必要とするNewton-Raphson法など）を併用する場合と比較して、CI-NEBはエネルギーと一次微分（力）のみを使用するため、大規模な系（第一原理計算など）において圧倒的に有利です [cite: 90]。

7. 結論#

Henkelmanらによって提案されたCI-NEB法は、既存のNEB法が抱えていた「鞍点エネルギーの補間精度」という問題を、エレガントかつ低コストな方法で解決しました。

厳密性: Climbing Imageの導入により、バンドはMEPを維持しつつ、その最高点（鞍点）に厳密に収束します。
効率性: 追加の計算コストは不要であり、可変バネ定数を併用することで、重要な遷移領域の記述能力を最大化できます。
汎用性: $CH_4/Ir(111)$ や $H_2/Si(100)$ といった複雑な表面反応系においても、その有効性が実証されました。

今日において、CI-NEB法は、VASPをはじめとする主要な第一原理計算コードや、ASE (Atomic Simulation Environment) などのライブラリに標準実装されており、計算化学における反応経路探索のデファクトスタンダードとしての地位を確立しています。

参考文献#

[1] Henkelman, G.; Uberuaga, B. P.; Jónsson, H. “A climbing image nudged elastic band method for finding saddle points and minimum energy paths.” J. Chem. Phys. 2000, 113, 9901–9904. DOI: 10.1063/1.1329672 [2] Eyring, H. J. Chem. Phys. 1935, 3, 107. [3] Wigner, E. Trans. Faraday Soc. 1938, 34, 29. [4] Keck, J. C. Adv. Chem. 1967, 13, 85. [5] Pechukas, P. In Dynamics of Molecular Collisions, Part B; Miller, W. H., Ed.; Plenum: New York, 1976. [6] Voter, A. F.; Doll, J. D. J. Chem. Phys. 1984, 80, 5832. [7] Vineyard, G. H. J. Phys. Chem. Solids 1957, 3, 121. [8] Mills, G.; Jónsson, H.; Schenter, G. K. Surf. Sci. 1995, 324, 305. [9] Jónsson, H.; Mills, G.; Jacobsen, K. W. In Classical and Quantum Dynamics in Condensed Phase Simulations; Berne, B. J., Ciccotti, G., Coker, D. F., Eds.; World Scientific: Singapore, 1998; p 385. [10] Henkelman, G.; Jónsson, H. J. Chem. Phys. 2000, 113, 9978. (Improved Tangent Method) [11] Elber, R.; Karplus, M. Chem. Phys. Lett. 1987, 139, 375. [12] Kresse, G.; Furthmüller, J. Phys. Rev. B 1996, 54, 11169.