陰的溶媒和の理論:PCM

陰的溶媒和の理論:PCM#

溶質を明示的な溶媒分子の殻で取り囲んでモデル化することは高価であり、既に高価な計算に相当な配座 サンプリングを加えます。**分極連続体モデル(PCM)**は代わりに、溶媒を、分子形の空洞の外側にある 構造を持たない誘電体連続体で置き換え、その連続体が溶質自身の静電ポテンシャルに応答して発達させる 表面電荷を解きます — その反応場が今度は溶質自身を分極させ、SCF 自己無撞着ループに結合します。 本章は IEF-PCM と C-PCM の背後にある境界積分静電学、分子表面を離散化する GePol 空洞、そして非静電的 な溶媒和物理を加える SMD 拡張を導出します — qc-rs の pcm-core クレート(純 Rust の、PCMSolver に 忠実な移植)と AGENTS.md の pcm-core 規則に基づきます。

境界積分の描像#

PCM の静電学は完全な 3 次元誘電体連続体問題を 2 次元問題に還元します:空洞の外側のいたるところで Poisson 方程式を解く代わりに、連続体が生成するはずの正しい反応ポテンシャルを再現する見かけの 表面電荷 \(q\) — 離散化された空洞表面上の点電荷の集合 — を見つければ十分です。空洞表面はテッセラ と呼ばれる小さな多角形パッチに分割され、それぞれが面積 \(a_k\)、中心、外向き法線を持ちます;未知数は テッセラごとに 1 つの見かけの電荷です。真空 Green 関数 \(G(\mathbf r,\mathbf r')=1/|\mathbf r-\mathbf r'|\) から構築される 2 つの古典的境界積分演算子が、テッセラ \(i,j\) の間で作用します:

\[ S_{ij} = G(\mathbf r_i,\mathbf r_j), \qquad D_{ij} = \frac{\partial G}{\partial\mathbf n_i}\bigg|_{\mathbf r_j} = \hat{\mathbf n}_i\cdot\frac{\mathbf r_j-\mathbf r_i}{|\mathbf r_j-\mathbf r_i|^3} \ \ (i\ne j), \]

単一層(\(S\)、テッセラ間の通常の Coulomb 様カーネル)と 2 重層(\(D\)、その法線微分 — テッセラ \(j\) にある点電荷がテッセラ \(i\) の外向き法線に沿って投影する場)演算子です — PCMSolver の符号 規約であり、qc-rs はこれに厳密に一致します:\(\hat{\mathbf n}_{\text{probe}}\cdot(\text{source}- \text{probe})/|\text{source}-\text{probe}|^3\)。対角(\(i=j\))要素は、カーネルがそこで特異なので特別な 扱いが必要です;qc-rs が実装する重心法(collocation)は、素朴な(発散する)点評価の代わりに閉形式の 対角公式を使います:

\[ S_{ii} = c\sqrt{\frac{4\pi}{a_i}}, \qquad D_{ii} = -c\sqrt{\frac{\pi}{a_i}}\cdot\frac{1}{R_i}, \]

\(c\approx1.07\)(PCMSolver の既定 collocation 係数)、\(R_i\) はテッセラ \(i\) が属する球の半径です — これらの対角公式こそが、あらゆるテッセラで真の特異積分求積を必要とするのではなく、collocation を 頑健で安価な既定にしているものです。

IEF-PCM:一般境界積分解法#

**積分方程式形式(IEF)**は一般的な PCM 定式化であり、等方性誘電体だけでなくより特殊な外部環境 (イオン液体、異方性媒質)にも有効です。空洞外側の誘電率 \(\varepsilon\) の等方性誘電体と、内側の 真空(\(\varepsilon=1\))について、\(A\) をテッセラ面積の対角行列とすると、系の行列は

\[ \mathbf T(\varepsilon) = \Bigl[2\pi f(\varepsilon)\,\mathbf I - \mathbf D_i\mathbf A\Bigr]\,\mathbf S_i, \qquad f(\varepsilon) = \frac{\varepsilon+1}{\varepsilon-1}, \qquad \mathbf R_\infty = 2\pi\,\mathbf I - \mathbf D_i\mathbf A, \]

内側(真空)Green 関数演算子 \(S_i,D_i\) から構築され、見かけの表面電荷は線形系を解きます:

\[ \mathbf T(\varepsilon)\,\mathbf q = -\mathbf R_\infty\,\mathbf v, \]

\(\mathbf v\) は各テッセラでサンプリングされた溶質の静電ポテンシャルです(核と電子密度の両方からの 分子静電ポテンシャル、まさに本マニュアルの他所で導出される ESP そのものであり、一般のグリッドではなく空洞表面で評価されているだけです)。因子 \(f(\varepsilon)\) は 周囲の媒質が溶質の場をどれほど強く遮蔽するかを符号化します — \(\varepsilon\to\infty\)(完全導体) では \(f(\varepsilon)\to1\) となり、IEF 定式化は次節 C-PCM で導出される導体様の極限になめらかに帰着 します。真に非一様な外部環境(イオン液体、異方性誘電体)については、異方性一般化が単一の Green 関数を別々の内側/外側演算子 \(S_o,D_o\) で置き換えます(等方性の場合の厳密な上位集合であり、異なる 定式化ではありません):

\[ \mathbf T = (2\pi\mathbf I - \mathbf D_o\mathbf A)\mathbf S_i + \mathbf S_o(2\pi\mathbf I + \mathbf A\mathbf D_i^{\mathsf T}), \qquad \mathbf R = (2\pi\mathbf I - \mathbf D_o\mathbf A) - \mathbf S_o\mathbf S_i^{-1}(2\pi\mathbf I - \mathbf D_i\mathbf A), \]

同じ電荷方程式 \(\mathbf T\mathbf q=-\mathbf R\mathbf v\) を伴います。qc-rs は非一様な外側 Green 関数 を IEF-PCM に限定します — C-PCM のより単純な導出(下記)は最初から一様な外側誘電体を仮定しており、 類似の異方性一般化を持ちません。

C-PCM:導体様の極限#

C-PCM(COSMO とも呼ばれる)は物理的により単純な描像から出発します:溶媒を有限 \(\varepsilon\) 誘電体ではなく完全導体だと見なし、はるかに単純な導体境界条件を解き、その後経験的なスケーリング 補正を適用して有限 \(\varepsilon\) の振る舞いを回復します。理想導体表面では反応ポテンシャルが溶質 自身のポテンシャルを厳密に相殺し、非スケール導体電荷方程式 \(\mathbf S\mathbf q=-\mathbf v\) を 与えます;経験的な修正は単一層行列を \(\varepsilon\) に依存する因子でスケールします:

\[ \mathbf S_{\text{scaled}} = \frac{\mathbf S}{f_\varepsilon}, \qquad f_\varepsilon = \frac{\varepsilon-1}{\varepsilon+k}, \qquad \mathbf S_{\text{scaled}}\,\mathbf q = -\mathbf v, \]

\(k\) は小さな経験的補正です(qc-rs は PySCF の規約に従います:素の C-PCM では \(k=0\)、密接に関連する COSMO 変種では \(k=0.5\) — 同じ背後のソルバ、異なるスケーリング定数だけ)。C-PCM の導出は 2 重層 演算子 \(D\) をまったく必要としません — 導体境界条件は純粋に単一層(\(S\))の関係だからです — これが IEF-PCM より組み立てが著しく安価である理由であり、非一様な外部環境への自然な拡張を持たない理由でも あります(そもそも一般化すべき \(D_o/S_o\) の分割がありません)。通常の誘電率のほとんどの溶媒では、 C-PCM と IEF-PCM は近く一致します(どちらも同じ物理的反応場を、異なる閉形式経路で近似している だけです)が、これらは同一の定式化ではなく、非常に小さな \(\varepsilon\) や特殊な幾何ではより大きく 分岐しえます — これは同じ答えへの単なる 2 つの数値経路ではなく、解かれる方程式における真の違いです。

GePol 空洞:原子球からテッセレーションされた表面へ#

どちらのソルバも、入力として位置・法線・面積を持つテッセラの集合、すなわち離散化された空洞表面を 必要とします。qc-rs のGePol空洞(PCMSolver プロジェクトの PEDRA アルゴリズムの Rust 移植)は これを、原子中心の球の和集合(標準テーブルからの半径 — Bondi あるいは類似の van der Waals 半径、 通常 \(\sim1.2\) の因子でスケール)から、特定の意図的な構成を使って構築します:アルゴリズムは常に 単一球の \(D_{2h}\) テッセレーションから始まります — 球の 8 分の 1(1 つのオクタント)を直接 テッセレーションし、その後 3 つの鏡映面(\(Oyz\)、\(Oxz\)、\(Oxy\))を使ってそのパッチを完全な球に 複製します。これは、分子全体としてたまたま何らかの対称性を持つかどうかとは独立に、Abel 点群構造を 空洞構成アルゴリズム自体の中で直接利用します — これは 1 つの球のテッセラを効率的に生成する ための計算上の便宜であり、対称な分子だけが空洞を得るという主張ではありません。互いに十分近い原子 から来る重なり合う球(van der Waals 球が交差する)は、最終的な空洞が分子全体を、内部の隙間や偽の 再入表面なしに囲むように、その交差領域を除外する必要があります — 上記のテッセレーションと鏡映の 手順と組み合わされて物理的に妥当な分子表面を生成しなければならない、GePol 系空洞構成における標準的 で確立された幾何ステップです。

SCF への結合:反応場 Fock 項#

見かけの表面電荷 \(\mathbf q\) は静的な一度限りの補正ではありません — それは(\(\mathbf v\) を通じて) 溶質の電子密度に依存し、それが生む反応場は今度はその同じ密度を分極させるので、PCM は真に密度 依存の Fock 項として SCF ループに参加しなければなりません、qc-scf の PcmTerm です。各 SCF サイクル:現在の密度と原子核から \(\mathbf v\)(各テッセラでの静電ポテンシャル)を構築し、 どちらかのソルバ(IEF-PCM では \(\mathbf T\mathbf q=-\mathbf R_\infty\mathbf v\)、C-PCM では \(\mathbf S_{\text{scaled}}\mathbf q=-\mathbf v\))で \(\mathbf q\) を解き、その結果の反応場ポテンシャル を、ESP の章 が任意の点でのポテンシャル評価に使うのと 同じ int1e_rinv 型テッセラポテンシャル積分を通じて Fock 行列に加え戻します(ここでは一般のグリッド ではなく各テッセラで評価されます)。これは構造的に、線形応答の章 が 真に密度依存の環境項について記述するのと同じ継ぎ目パターンです:固定 ECP 演算子とは異なり、PCM の 反応場はどんな密度にも応答します — 安定性検査や CPHF 右辺で使われる遷移密度を含めて — なので、 溶媒和系で密度依存の摂動を探るあらゆる安定性解析・CPHF 解法・(原理的には)TDDFT 計算において、 単に基底状態 Fock 構築だけでなく、応答カーネル \(\sigma\) 構築にも自身の項を寄与させる必要が真に あります。

SMD:非静電的溶媒和物理を加える#

IEF-PCM と C-PCM は溶媒和の静電的部分だけ — 電荷分極からの反応場 — を捉えます。実際の溶媒和 自由エネルギーには相当な非静電的寄与もあります:溶媒中に空洞を形成するのに必要な仕事(空洞化)、 溶質-溶媒間の分散相互作用、そして溶質-溶媒界面での構造/交換反発効果です。SMD(密度に基づく 溶媒和モデル)は、幾何だけから構築される単一の追加項でこれを捉えます:

\[ E_{\text{solv}} = \underbrace{E_{\text{electrostatic}}(\text{IEF-PCM}, R_{\text{SMD}})}_{\text{同じ IEF-PCM 機構、SMD 自身の固有 Coulomb 半径を使う}} + \underbrace{G_{\text{CDS}}}_{\text{空洞-分散- 溶媒構造}}, \]

\(G_{\text{CDS}}\) は幾何のみのスカラーです — 原子および分子の表面張力様パラメータと溶媒アクセス 可能表面積(SASA)の積の和であり、追加の Fock 行列項も、それが包む IEF-PCM 部分を超える追加の SCF 自己無撞着性もありません。静電部分は SWIG IEF-PCM ソルバ経路を真に変更なく再利用します — SMD の 特徴的な寄与は、その共有機構に供給する半径(既定で IEF-PCM/C-PCM が使う素の van der Waals 半径とは 異なる、SMD 自身の固有 Coulomb 半径)と、その後加えられる純粋に幾何学的な \(G_{\text{CDS}}\) 補正 だけです。

検証済み例 — 水、真空中と 3 つの溶媒和処理:

import qc
water = "O 0 0 0.117; H 0 0.757 -0.469; H 0 -0.757 -0.469"

m_gas = qc.chk.new(atom=water, ao="sto-3g", unit="angstrom").scf(ref="r").run()
m_gas.scf.energy   # -74.96294665653868

m_ief = qc.chk.new(atom=water, ao="sto-3g", unit="angstrom").scf(ref="r",
                    pcm={"solvent": "water", "model": "IEF-PCM"}).run()
m_ief.scf.energy, m_ief.scf.energy - m_gas.scf.energy
# (-74.9690349700267, -0.006088313488021413)

m_cpcm = qc.chk.new(atom=water, ao="sto-3g", unit="angstrom").scf(ref="r",
                     pcm={"solvent": "water", "model": "C-PCM"}).run()
m_cpcm.scf.energy, m_cpcm.scf.energy - m_gas.scf.energy
# (-74.96907206979104, -0.006125413252362932)

m_smd = qc.chk.new(atom=water, ao="sto-3g", unit="angstrom").scf(ref="r",
                    pcm={"solvent": "water", "model": "smd"}).run()
m_smd.scf.energy, m_smd.scf.energy - m_gas.scf.energy
# (-74.97017942182167, -0.007232765282992659)

3 つの処理すべてが真空に対してエネルギーを安定化させます、極性溶媒中の極性溶質として物理的に 期待される通りです。IEF-PCM と C-PCM は近く一致しますが同一ではありません(ここでは \(\sim4\times10^{-5}\ E_h\) の差)— 上で導出したまさに「密接に関連するが真に異なる方程式」の 関係であり、数値精度の副産物ではありません。SMD のより大きな安定化は、自身の異なる半径で構築 された静電部分に加えて、追加の非静電 \(G_{\text{CDS}}\) 寄与を反映しています — 単に再スケール された IEF-PCM の数値ではありません。

練習 14

  1. C-PCM の導出は 2 重層演算子 \(D\) を一切導入せず、\(S\) だけです。導体境界条件(C-PCM の導出が 出発する物理的出発点)が \(D\) を不要にする理由、そして IEF-PCM の有限 \(\varepsilon\) 境界条件が それを避けられない理由を一文で説明しなさい。

  2. GePol 空洞構成は、実際の分子が何らかの対称性をまったく持たない場合でも、常に単一球の \(D_{2h}\)(1 オクタント、3 回鏡映)をテッセレーションします。これがなぜ矛盾ではないのか — 対称 な単一球テッセレーションを構築することが、点群対称性をまったく持たない分子の最終的な空洞を 生成することとどう両立するのか説明しなさい。

  3. PCM は(基底状態 Fock 構築だけでなく)線形応答 \(\sigma\) 構築にも項を寄与させなければなりません が、固定 ECP 演算子はそうしません。線形応答の章 自身の、外部 ポテンシャルが応答寄与を必要とするときの基準を使って、この違いを一文で説明しなさい。

溶媒和と ECP は Part II における qc-rs の 2 つの Hamiltonian 修正の章です;分散理論 は、それらのどちらか(あるいはどちらでもない)とも組み合わさる、幾何のみの 3 つ目の補正を扱います。 Hartree-Fock もほとんどの密度汎関数も、それ自体では長距離 van der Waals 引力を捉えないため必要と なるものです。