查看“︁韦尔莱积分法”︁的源代码

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

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

== 基本韦尔莱算法 ==

牛顿运动方程为
:<math>
a(t)=\frac{f(t)}{m}
</math>

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

:<math>
r(t+\Delta t)=r(t)+v(t)\Delta t+\frac{a(t)}{2}\Delta t^{2}+\frac{1}{3!}\frac{d^3r}{dt^3}\Delta t^{3}+\mathcal O(\Delta t^{4})
</math>

得

:<math>
r(t+\Delta t)=r(t)+v(t)\Delta t+\frac{f(t)}{2m}\Delta t^{2}+\frac{1}{3!}\frac{d^3r}{dt^3}\Delta t^{3}+\mathcal O(\Delta t^{4})
</math>

同理

:<math>
r(t-\Delta t)=r(t)-v(t)\Delta t+\frac{f(t)}{2m}\Delta t^{2}-\frac{1}{3!}\frac{d^3r}{dt^3}\Delta t^{3}+\mathcal O(\Delta t^{4})
:</math>

两式相加，得

:<math>
r(t+\Delta t)+r(t-\Delta t)=2r(t)+\frac{f(t)}{m}\Delta t^{2}+\mathcal O(\Delta t^{4})
</math>

则

:<math>
r(t+\Delta t)\simeq2r(t)-r(t-\Delta t)+\frac{f(t)}{m}\Delta t^{2}
</math>

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

:<math>
v(t)=\frac{r(t+\Delta t)-r(t-\Delta t)+\mathcal O(\Delta t^{3})}{2\Delta t}=\frac{r(t+\Delta t)-r(t-\Delta t)}{2\Delta t} + \mathcal O(\Delta t^{2})
</math>

== 速度表示的韦尔莱算法 ==

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

首先对位置进行 Taylor 展开

:<math>
r(t+\Delta t)=r(t)+v(t)\Delta t+\frac{f(t)}{2m}\Delta t^{2}
</math>
:<math>
r(t+2\Delta t)=r(t+\Delta t)+v(t+\Delta t)\Delta t+\frac{f(t+\Delta t)}{2m}\Delta t^{2}
</math>
对两式相减可得

:<math>
r(t+2\Delta t)+r(t)=2r(t+\Delta t)+\left[v(t+\Delta t)-v(t)\right]\Delta t+\frac{f(t+\Delta t)-f(t)}{2m}\Delta t^{2}
</math>
将最初的 Verlet 公式中的 <math>t</math> 换为 <math>t+\Delta t</math> ，

:<math>
r(t+2\Delta t)\simeq2r(t+\Delta t)-r(t)+\frac{f(t+\Delta t)}{m}\Delta t^{2}
</math>

代入前式，可得

:<math>
v(t+\Delta t)=v(t)+\frac{f(t+\Delta t)+f(t)}{2m}\Delta t
</math>

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

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

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

3.由速度表示的韦尔莱公式求出新的速度 <math>v(t+\Delta t)</math> 。

== 参考文献 ==
*{{cite book|author=Daan Frenkel|chapter=第四章|title="Understanding Molecular Simulation - From Algorithms to Applications" Academic Press|year=2002}}
{{常微分方程数值方法}}
[[Category:数值微分方程]]