# 原子電荷と結合次数

「酸素にどれだけ負電荷があるか」「これは単結合か二重結合か」は化学で最もよくある問いの 2 つ — そして最も
*微妙な*2 つでもあります。なぜなら原子電荷も結合次数も量子力学的な観測量ではないからです。これらは**分割の
スキーム**:分子の密度を原子と原子対に分ける処方箋です。qc-rs は多数を実装しており、本章はそれらの計算方法と、
違いの読み方を示します。

## 理論:なぜ電荷スキームがこんなに多いのか

電子密度 $\rho(\mathbf r)$ は実在で一義的ですが、「*原子 A の*電荷」を言うには、**原子 A がどこで終わり原子 B が
どこで始まるか**を決めねばならず — 一意な答えはありません。スキームごとに選択が違います:

- **軌道ベース**（Mulliken, Löwdin）:密度をどの原子の基底関数が担うかで分ける。単純で速いが、**基底依存**
  （特に Mulliken）。
- **密度ベース**（Hirshfeld, MBIS, VDD）:原子参照密度を使って実空間を分割する。より頑健で基底安定。
- **静電ポテンシャルフィット**（MK/ChelpG/RESP）:分子の ESP を最もよく再現する電荷を選ぶ — 力場パラメータ化に
  最適。
- **自然母集団**（NPA）:自然原子軌道から;基底安定で広く使われる。

「正しい」スキームはなく、各々わずかに異なる問いに答えます。教訓は**目的に合ったものを 1 つ選び、一貫させる**
こと、そして最後の桁を過剰解釈しないことです。

## 原子電荷

各スキームは `qc.prop.chrg` の leaf です。2 つの軌道ベーススキームが最も単純な式を持ち — どちらも AO ごとの
**総母集団** $P_\mu$ を原子ごとの電荷 $q_A = Z_A - \sum_{\mu\in A} P_\mu$ に簡約します — 違いは*どの*母集団を
使うかだけです:

$$
\text{Mulliken: } P_\mu = (\mathbf{DS})_{\mu\mu} = \sum_\nu D_{\mu\nu}S_{\mu\nu},
\qquad
\text{Löwdin: } P_\mu = \big(\mathbf S^{1/2}\mathbf D\,\mathbf S^{1/2}\big)_{\mu\mu}.
$$

Mulliken の $\mathbf{DS}$ 母集団は $\mu,\nu$ について対称でなく、拡散した重なりの乏しい AO では*負*にさえ
なりえます — 下の表が示す基底感受性です。Löwdin の**対称直交化**（両側に $\mathbf S^{1/2}$）はこの非対称性を
直し、これが Löwdin が Mulliken より穏やかになりがちな理由の一部です。密度ベースのスキーム（Hirshfeld, MBIS）
は代わりに*実空間*密度 $\rho(\mathbf r)$ を、参照原子密度を使って直接分割し、AO 母集団を一切必要としません —
根本的に異なる（そしてより基底安定な）処方箋で、[電荷分割解析の理論](../../10-foundations/population-analysis-theory.md)
で詳述します。leaf によっては素の原子ごと配列を、あるいは `"charges"` フィールドを持つより豊かなレコードを
返します:

```python
import qc, numpy as np
water = "O 0 0 0.117; H 0 0.757 -0.469; H 0 -0.757 -0.469"
m = qc.chk.new(atom=water, ao="cc-pvdz", unit="angstrom").scf(ref="r").run()

np.asarray(qc.prop.chrg.hirshfeld(m))          # 配列: [-0.3219, 0.161, 0.161]
qc.prop.chrg.mulliken(m)["charges"]            # レコード: [-0.3056, 0.1528, 0.1528]
qc.prop.chrg.npa(m)["charges"]                 # [-0.9136, 0.4568, 0.4568]
```

水の酸素は**どのスキームでも**負、水素は正 — 定性的な描像は頑健です — が、*大きさ*は処方箋で大きく変わります:

| スキーム | `qc.prop.chrg.` | O 電荷 | H 電荷 | 特徴 |
|---|---|---|---|---|
| Löwdin | `lowdin` | −0.094 | +0.047 | 軌道、対称化 |
| Mulliken | `mulliken` | −0.306 | +0.153 | 軌道、基底に敏感 |
| Hirshfeld | `hirshfeld` | −0.322 | +0.161 | 密度、「ストックホルダー」 |
| CM5 | `cm5` | −0.658 | +0.329 | Hirshfeld + 経験的補正 |
| ADCH | `adch` | −0.731 | +0.365 | 原子双極子補正 Hirshfeld |
| MBIS | `mbis` | −0.862 | +0.431 | 最小基底反復ストックホルダー |
| NPA | `npa` | −0.914 | +0.457 | 自然母集団 |

*同じ*酸素で −0.09 から −0.91 という広がり — これこそ、どのスキームを使ったかを明記すべき理由です。一連の分子で
比較するなら、**密度ベース**（Hirshfeld/MBIS）または **NPA** が Mulliken より安全な選択です。

## 結合次数

**結合次数**は 2 原子間で共有される電子対の数を定量します。2 つの主力が `qc.prop.bond` にあり、各々、原子–原子の
`"matrix"` 全体と原子ごとの `"valence"` を返します:

```python
mayer = qc.prop.bond.mayer(m)
mayer["matrix"][0, 1]     # 1.0207   O–H の Mayer 結合次数
mayer["valence"][0]       # 2.041    酸素の総原子価（≈ 2 結合）

qc.prop.bond.wiberg(m)["matrix"][0, 1]   # 1.1886   O–H の Wiberg 結合次数
```

O–H の **Mayer** 結合次数は ~1.02 — 化学が期待するとおりの単結合で、酸素の原子価 ~2.04 は 2 本の O–H 結合を正しく
反映します。2 つの結合次数は、上で Mulliken 電荷を導いたのと*同じ*密度×重なり行列積を 2 乗し、和の取り方だけを
変えたものです:

$$
B_{AB}^{\text{Mayer}} = 2\!\!\sum_{\mu\in A,\,\nu\in B}\!\!\Big[(\mathbf P^\alpha\mathbf S)_{\mu\nu}(\mathbf P^\alpha\mathbf S)_{\nu\mu} + (\mathbf P^\beta\mathbf S)_{\mu\nu}(\mathbf P^\beta\mathbf S)_{\nu\mu}\Big],
\qquad
B_{AB}^{\text{Wiberg}} = 2\!\!\sum_{\mu\in A,\,\nu\in B}\!\!\Big[(\tilde{\mathbf P}^\alpha)_{\mu\nu}^2 + (\tilde{\mathbf P}^\beta)_{\mu\nu}^2\Big],
$$

$\mathbf P^\alpha,\mathbf P^\beta$ はスピン密度行列、$\tilde{\mathbf P}=\mathbf S^{1/2}\mathbf P\,\mathbf
S^{1/2}$ は Löwdin 直交化した密度です。**Wiberg はまさに直交化基底で評価した Mayer** — そこでは重なりが恒等
なので式は密度行列要素の 2 乗和に簡約されます — これがわずかに異なる数値（ここでは Mayer の 1.02 に対し
1.19）を与える理由です:2 つは同じ*考え*（$A$ と $B$ の間でどれだけ密度が共有されるか）を、2 つの異なる基底で
測ったものです。各原子の総**原子価** $\sum_{B\ne A}B_{AB}$ は対角を除いた行の和です。`qc.prop.bond` には
fuzzy 原子・multicenter（非局在/芳香族結合用）・IBSI・非局在化指数の結合尺度もあります。

## 総合例:電荷と結合を一緒に

```python
import qc, numpy as np
water = "O 0 0 0.117; H 0 0.757 -0.469; H 0 -0.757 -0.469"
m = qc.chk.new(atom=water, ao="cc-pvdz", unit="angstrom").scf(ref="r").run()

print("Hirshfeld charges:", np.round(np.asarray(m.prop.chrg.hirshfeld()), 3))   # [-0.322 0.161 0.161]
print("NPA charges      :", np.round(np.asarray(m.prop.chrg.npa()["charges"]), 3))
print("O–H Mayer BO     :", round(m.prop.bond.mayer()["matrix"][0, 1], 3))       # 1.021
```

:::{exercise}
:label: ex-charges

1. 2 本の論文が水の酸素電荷を −0.31 と −0.91 と報告しています。両方正しくありえますか。それらを整合させる 1 つの
   情報は何ですか。
2. 古典力場をフィットしていて点電荷が必要です。どの*ファミリ*の電荷スキームを使うべきで、なぜ Mulliken では
   ないのですか。
3. 水の酸素の Mayer 原子価は ≈ 2.04 になります。この数は何を教え、なぜちょうど 2 でないのですか。
:::

:::{solution} ex-charges
:class: dropdown

1. はい — 原子電荷は観測量ではありません。−0.31 は Mulliken 系、−0.91 は NPA 系です;**分割スキーム**を明記すれば
   整合します。比較は 1 つのスキーム内でのみ意味があります。
2. **静電ポテンシャルフィット**スキーム（MK / ChelpG / RESP、`qc.prop.chrg` 内） — 電荷は力場がモデル化すべき
   ESP を再現するよう選ばれます。Mulliken 電荷は ESP にフィットされず基底に敏感なので、転用性が低いです。
3. 酸素が ~2 結合分の共有電子対に関与している（2 本の O–H 結合）ことを示し — その二価性と整合します。ちょうど 2
   でないのは、実際の結合がわずかに分極/イオン性を持ち、それを Mayer の共有結合次数が完全な共有として数えない
   ためです。
:::

電荷と結合次数は密度を*数値的に*分割します。[次章](qtaim-elf.md) は密度を*トポロジー的に*分割します — QTAIM と
ELF が、密度そのものの形から原子・結合・孤立電子対を見つけます。