韦尔莱积分法

韦尔莱算法是一种用于求解牛顿运动方程的数值方法，被广泛应用于分子动力学模拟以及视频游戏中。韦尔莱算法的优点在于：数值稳定性比简单的欧拉方法高很多，并保持了物理系统中的时间可逆性与相空间体积元体积守恒的性质。

Carl Størmer首次应用韦尔莱算法求解磁场中运动粒子的轨迹，因此韦尔莱算法又被称为Størmer算法。1967年法国物理学家Loup Verlet将其应用于分子动力学计算，从此韦尔莱算法流行起来。

基本韦尔莱算法

牛顿运动方程为

a (t) = \frac{f (t)}{m}

代入到粒子的坐标关于时间步的Taylor展式中

r (t + Δ t) = r (t) + v (t) Δ t + \frac{a (t)}{2} Δ t^{2} + \frac{1}{3!} \frac{d^{3} r}{d t^{3}} Δ t^{3} + 𝒪 (Δ t^{4})

得

r (t + Δ t) = r (t) + v (t) Δ t + \frac{f (t)}{2 m} Δ t^{2} + \frac{1}{3!} \frac{d^{3} r}{d t^{3}} Δ t^{3} + 𝒪 (Δ t^{4})

同理

r (t - Δ t) = r (t) - v (t) Δ t + \frac{f (t)}{2 m} Δ t^{2} - \frac{1}{3!} \frac{d^{3} r}{d t^{3}} Δ t^{3} + 𝒪 (Δ t^{4}) :

两式相加，得

r (t + Δ t) + r (t - Δ t) = 2 r (t) + \frac{f (t)}{m} Δ t^{2} + 𝒪 (Δ t^{4})

则

r (t + Δ t) ≃ 2 r (t) - r (t - Δ t) + \frac{f (t)}{m} Δ t^{2}

新位置的计算误差为四阶， $Δ t$ 为时间步。因为韦尔莱算法中不涉及速度，如果希望得到速度，须从轨线中推导速度表达式：

v (t) = \frac{r (t + Δ t) - r (t - Δ t) + 𝒪 (Δ t^{3})}{2 Δ t} = \frac{r (t + Δ t) - r (t - Δ t)}{2 Δ t} + 𝒪 (Δ t^{2})

一般地，速度表示下的韦尔莱算法更为常用，它可以给出同一时间变量下的速度和位置。它实际上与基本的韦尔莱算法等价，精度相同。

首先对位置进行 Taylor 展开

r (t + Δ t) = r (t) + v (t) Δ t + \frac{f (t)}{2 m} Δ t^{2}

r (t + 2 Δ t) = r (t + Δ t) + v (t + Δ t) Δ t + \frac{f (t + Δ t)}{2 m} Δ t^{2}

对两式相减可得

r (t + 2 Δ t) + r (t) = 2 r (t + Δ t) + [v (t + Δ t) - v (t)] Δ t + \frac{f (t + Δ t) - f (t)}{2 m} Δ t^{2}

将最初的 Verlet 公式中的 $t$ 换为 $t + Δ t$ ，

r (t + 2 Δ t) ≃ 2 r (t + Δ t) - r (t) + \frac{f (t + Δ t)}{m} Δ t^{2}

代入前式，可得

v (t + Δ t) = v (t) + \frac{f (t + Δ t) + f (t)}{2 m} Δ t

此式即为速度表示的韦尔莱算法。实际常用的计算步骤为，

1.首先通过 Taylor 展开 $r (t + Δ t) = r (t) + v (t) Δ t + \frac{f (t)}{2 m} Δ t^{2}$ 计算得到位置 $r (t + Δ t)$

2.由 $r (t + Δ t)$ 和系统的相互作用势条件（如果相互作用仅依赖位置 $r$ ）可以求的力场 $f (t + Δ t)$

3.由速度表示的韦尔莱公式求出新的速度 $v (t + Δ t)$ 。