埃瓦尔德求和

埃瓦尔德求和（），是一种计算中长程力（如静电力）的方法，以德国物理学家保罗·彼得·埃瓦尔德命名。埃瓦尔德求和最初用于计算离子晶体的电势能，现在用于计算化学中计算长程力。埃瓦尔德求和是泊松求和公式的特殊形式，用倒空间中的等效求和代替实空间中的总和。埃瓦尔德求和将分为短程力和无奇点的长程力两部分，短程力在实空间中计算，长程力用傅里叶变换计算。与直接求和相比，此方法的优势为能量能够快速收敛，这意味着此方法在计算长程力时具有较高的精度和合理的速度，是计算中长程力的标准方法。此方法需要分子系统的电中性，以准确计算总库仑力。

推导
埃瓦尔德求和将表示为两部分之和：

:\varphi(\mathbf{r}) \ \stackrel{\mathrm{def}}{=}\ \varphi_{sr}(\mathbf{r}) + \varphi_{\ell r}(\mathbf{r}),

其中，\varphi_{sr}(\mathbf{r})表示实空间中和值快速收敛的短程势，\varphi_{\ell r}(\mathbf{r})表示倒空间中和值快速收敛的长程势。所有量（如r）的长程部分是有限的，但可能有简易的数学形式，如高斯分布。该方法假设短程势容易求和，因此需要重点考虑的是长程势。由于使用了傅里叶级数，该方法将周期性边界条件作为假设，此周期性系统的重复单元称为原胞，选择一个原胞作为中央原胞作为参考，其余单元称为镜像。

长程力的能量是中央原胞的电荷与晶格所有电荷间之和，因此可以表示为原胞和晶格的电荷密度的双重积分：

:
E_{\ell r} = \iint d\mathbf{r}\, d\mathbf{r}^\prime\, \rho_\text{TOT}(\mathbf{r}) \rho_{uc}(\mathbf{r}^\prime) \ \varphi_{\ell r}(\mathbf{r} - \mathbf{r}^\prime)

其中原胞的电荷密度\rho_{uc}(\mathbf{r})是中央原胞中位置\mathbf{r}_k上的电量q_k之和：

:
\rho_{uc}(\mathbf{r}) \ \stackrel{\mathrm{def}}{=}\ \sum_{\mathrm{charges}\ k} q_k \delta(\mathbf{r} - \mathbf{r}_k)

总电荷密度\rho_\text{TOT}(\mathbf{r})是原胞及其镜像电量q_{k}之和：

:
\rho_\text{TOT}(\mathbf{r}) \ \stackrel{\mathrm{def}}{=}\ \sum_{n_1, n_2, n_3} \sum_{\mathrm{charges}\ k}
q_k \delta(\mathbf{r} - \mathbf{r}_k - n_1 \mathbf{a}_1 - n_2 \mathbf{a}_2 - n_3 \mathbf{a}_3)

这里，\delta(\mathbf{x})表示狄拉克δ函数，\mathbf{a}_1、\mathbf{a}_2、\mathbf{a}_3表示晶格矢量，n_1、n_2、n_3的范围为所有整数。总电荷密度\rho_\text{TOT}(\mathbf{r})可以表示为\rho_{uc}(\mathbf{r})与晶格函数L(\mathbf{r})的卷积：

:
L(\mathbf{r}) \ \stackrel{\mathrm{def}}{=}\ \sum_{n_1, n_2, n_3}
\delta(\mathbf{r} - n_1 \mathbf{a}_{1} - n_{2} \mathbf{a}_2 - n_3 \mathbf{a}_3)

由于\rho_\text{TOT}(\mathbf{r})为卷积，其傅里叶变换为一个积：

:
\tilde{\rho}_\text{TOT}(\mathbf{k}) = \tilde{L}(\mathbf{k}) \tilde{\rho}_{uc}(\mathbf{k})

其中晶格函数的傅里叶变换是狄拉克δ函数的另一个和：

:
\tilde{L}(\mathbf{k}) =
\frac{\left(2\pi \right)^{3}}{\Omega} \sum_{m_1, m_2, m_3}
\delta(\mathbf{k} - m_1 \mathbf{b}_1 - m_2 \mathbf{b}_2 - m_3 \mathbf{b}_3)

其中定义倒空间向量为\mathbf{b}_{1} \ \stackrel{\mathrm{def}}{=}\ \frac{\mathbf{a}_{2} \times \mathbf{a}_{3}}{\Omega}（周期性排列），其中\Omega \ \stackrel{\mathrm{def}}{=}\ \mathbf{a}_{1} \cdot \left( \mathbf{a}_{2} \times \mathbf{a}_{3} \right)为中心原胞的体积（几何形状通常为平行六面体），L(\mathbf{r})和\tilde{L}(\mathbf{k})为实函数和偶函数。

为了简洁起见，定义有效单粒子势能：

:
v(\mathbf{r}) \ \stackrel{\mathrm{def}}{=}\ \int d\mathbf{r}^{\prime}\, \rho_{uc}(\mathbf{r}^\prime) \ \varphi_{\ell r}(\mathbf{r} - \mathbf{r}^\prime)

因为其亦为卷积，其傅里叶变换是一个积：

:
\tilde{V}(\mathbf{k}) \ \stackrel{\mathrm{def}}{=}\ \tilde{\rho}_{uc}(\mathbf{k}) \tilde{\Phi}(\mathbf{k})

其中定义了傅里叶变换：

:
\tilde{V}(\mathbf{k}) = \int d\mathbf{r} \ v(\mathbf{r}) \ e^{-i\mathbf{k} \cdot \mathbf{r}}

现在，长程力的能量可以表示为单个电荷密度的积分：

:
E_{\ell r} = \int d\mathbf{r} \ \rho_\text{TOT}(\mathbf{r}) \ v(\mathbf{r})

使用帕塞瓦尔定理，能量亦可于倒空间中求和：

:
E_{\ell r} =
\int \frac{d\mathbf{k}}{\left(2\pi\right)^3} \ \tilde{\rho}_\text{TOT}^*(\mathbf{k}) \tilde{V}(\mathbf{k}) =
\int \frac{d\mathbf{k}}{\left(2\pi\right)^3} \tilde{L}^*(\mathbf{k}) \left| \tilde{\rho}_{uc}(\mathbf{k})\right|^2 \tilde{\Phi}(\mathbf{k}) =
\frac{1}{\Omega} \sum_{m_1, m_2, m_3} \left| \tilde{\rho}_{uc}(\mathbf{k})\right|^2 \tilde{\Phi}(\mathbf{k})

其中\mathbf{k} = m_1 \mathbf{b}_1 + m_2 \mathbf{b}_2 + m_3 \mathbf{b}_3是最终的和值。

计算出\tilde{\rho}_{uc}(\mathbf{k})后，\mathbf{k}的和值或积分是显然的，可以很快地收敛。不能收敛的最常见原因是原胞不太明确，其必须为电中性，以避免无穷大的和。

粒子网格埃瓦尔德（PME）方法
在计算机普及前，埃瓦尔德求和是理论物理的理论。然而，自20世纪70年代以来，埃瓦尔德求和在粒子系统的计算机模拟中被广泛使用，尤其是遵守平方反比定律的粒子相互作用，如重力和静电力。最近，粒子网格埃瓦尔德方法也用于计算兰纳-琼斯势的r^{-6}部分，以消除产生的。其应用包括等离子体、星系及分子的模拟。

在粒子网格埃瓦尔德方法中，和标准埃瓦尔德求和相同，被分为两部分\varphi(\mathbf{r}) \ \stackrel{\mathrm{def}}{=}\ \varphi_{sr}(\mathbf{r}) + \varphi_{\ell r}(\mathbf{r})，其基本思想是用实空间中短程力的直接求和E_{sr}（粒子部分），及倒空间中长程力的求和（埃瓦尔德部分），代替点粒子间相互作用的能量的直接求和：

:
E_\text{TOT} = \sum_{i,j} \varphi(\mathbf{r}_{j} - \mathbf{r}_i) = E_{sr} + E_{\ell r}

:
E_{sr} = \sum_{i,j} \varphi_{sr}(\mathbf{r}_j - \mathbf{r}_i)

:
E_{\ell r} = \sum_{\mathbf{k}} \tilde{\Phi}_{\ell r}(\mathbf{k}) \left| \tilde{\rho}(\mathbf{k}) \right|^2

其中\tilde{\Phi}_{\ell r}和\tilde{\rho}(\mathbf{k})表示力和电荷密度的傅里叶变换。由于两个求和分别在实空间和倒空间中迅速收敛，它们可能被精确，且所需计算时间大幅减少。计算电荷密度的傅里叶变换\tilde{\rho}(\mathbf{k})可使用快速傅里叶变换，需在空间中的上（即网格部分）估计电荷密度。

由于埃瓦尔德方法隐含的周期性假设，粒子网格埃瓦尔德方法于物理系统中的应用需施加周期性。因此，该方法最适合用于空间范围内可以模拟为无限的系统。在分子动力学模拟中，常构造可以无限平铺形成镜像的电中性原胞；然而，为了正确解释这种近似效应，这些镜像被重新并入原始模拟原胞中，这种整体效应被称为周期性边界条件。 想象一个单位立方体，上表面与下表面有效接触，右侧面与左侧面有效接触，前表面与后表面有效接触。因此，原胞的尺寸必须足够大，以避免两个接触面间不正确的运动相关性，但仍需足够小以便计算。短程力与长程力间的定义也可以引入。

电荷密度对网格的限制，使得粒子网格埃瓦尔德方法对电荷密度或势函数平滑变化的系统更有效。利用可以更有效地处理局部系统或电荷密度波动较大的系统。

偶极子
极性晶体（即原胞中具有净偶极子\mathbf{p}_{uc}的晶体）的静电能為条件收敛，即取决于求和顺序。例如，若中央原胞的偶极与不断增加的立方体上的原胞偶极相互作用，则其能量收敛值並不會与考慮不斷增大的球面時相等。大致来说，这种条件收敛是因为在半径为R的壳上的偶极子数約為R^{2}；偶极-偶极相互作用的强度約為\frac{1}{R^{3}}；而兩者相乘的結果是發散的调和级数\sum_{n=1}^{\infty} \frac{1}{n}。

這看似令人驚訝的結果並不與現實晶體能量有限的事實相違背，因為現實晶體並非無限，具有特定邊界。具体而言，极性晶体的边界的有效表面电荷密度为\sigma = \mathbf{P} \cdot \mathbf{n}，其中\mathbf{n}为表面法向量，\mathbf{P}为单位体积的净偶极矩。則中央原胞之偶极子與表面电荷密度\sigma的相互作用能U可寫為：

:
U = \frac{1}{2V_{uc}} \int
\frac{\left( \mathbf{p}_{uc}\cdot \mathbf{r} \right)
\left( \mathbf{p}_{uc} \cdot \mathbf{n} \right)dS}{r^3}

其中，\mathbf{p}_{uc}和V_{uc}分别为原胞的净偶极矩和体积，dS为晶面上的无穷小区域，\mathbf{r}为中央原胞到无穷小区域的向量。此公式来自于对能量 dU = -\mathbf{p}_{uc} \cdot \mathbf{dE}积分，其中d\mathbf{E}表示无穷小电场，由无穷小的表面电荷dq \stackrel{\mathrm{def}}{=}\sigma dS产生（库仑定律）：

:
d\mathbf{E} \ \stackrel{\mathrm{def}}{=}\
\left( \frac{-1}{4\pi\epsilon} \right) \frac{dq \ \mathbf{r}}{r^3} =
\left( \frac{-1}{4\pi\epsilon} \right)
\frac{\sigma\, dS \ \mathbf{r} }{r^3}

负号来自于\mathbf{r}的定义，其指向电荷方向为正方向。

历史
埃瓦尔德求和由德国物理学家保罗·彼得·埃瓦尔德于1921年发表，用于确定离子晶体的静电能及马德隆常数。

复杂度
不同的埃瓦尔德求和具有不同的时间复杂度。直接求和的时间复杂度为O(N^2)，其中N为系统中原子数。粒子网格埃瓦尔德方法的时间复杂度为O(N\,\log N)。

参见

• 保罗·彼得·埃瓦尔德
• 馬德隆常數
• 泊松求和公式
• 分子建模

*

参考文献