几何积分

在常微分方程的数值计算中，几何积分是一种保留微分方程的流的精确几何特性的数值方法。

以摆为例

可考虑单摆运动以引出几何积分的研究。

设摆锤质量为$m=1$，摆杆长度为$\ell =1$。设重力加速度为$g=1$。用$q(t)$表示杆偏移垂直方向的角位移，并用$p(t)$表示摆的动量，则系统的哈密顿量（动能与势能之和）为

$H(q,p)=T(p)+U(q)=\frac{1}{2}{p}^2-\cos q,$

其给出哈密顿方程

$(\dot{q},\dot{p})=(\partial H/\partial p,-\partial H/\partial q)=(p,-\sin q).$

很自然，可将所有$q$的位形空间$Q$看做单位圆${𝕊}^1$，这样$(q,p)$就位于圆柱体${𝕊}^1\times ℝ$上。取$(q,p)\in {ℝ}^2$只是因为$(q,p)$空间会更方便绘制。定义$z(t)=(q(t),p(t){)}^{\mathrm{T}}$、$f(z)=(p,-\sin q{)}^{\mathrm{T}}$。让我们用一些简单的数值方法对这个系统进行积分。像往常一样，选择常数步长$h$，对任意非负整数$k$记$z_k:=z(kh)$。 我们用以下方法：

$z_{k+1}=z_k+hf(z_k)$（显式欧拉）；

$z_{k+1}=z_k+hf(z_{k+1})$（隐式欧拉）；

$z_{k+1}=z_k+hf(q_k,p_{k+1})$（辛欧拉）；

$z_{k+1}=z_k+hf((z_{k+1}+z_k)/2)$（隐式中点法则）。

（注意，辛欧拉法用显式欧拉法处理q，用隐式欧拉法处理$p$。）

观察到$H$在哈密顿方程的解曲线上是常数，于是可以描述系统的精确轨迹，是$p^2/2-\cos q$的水平曲线。在${ℝ}^2$中绘制了系统的精确轨迹和数值解。对显式、隐式欧拉法，分别取$h=0.2$；z₀ = (0.5, 0)及(1.5, 0)；对其他两种方法，分别取$h=0.3$、z₀ = (0, 0.7)；(0, 1.4)及(0, 2.1)。

单摆轨迹

显式（或隐式）欧拉法是从原点向外（或向内）的螺旋运动。另两种方法显示了正确的定性行为，隐式中点法则与精确解的吻合程度高于辛欧拉法。

回顾一下，具有1自由度的哈密顿系统的精确流${\varphi }_t$是保面积的，即

$\det \frac{\partial {\varphi }_t}{\partial (q_0,p_0)}=1$ for all $t$.

此式很容易手动验证。对我们的单摆例子，可以发现，显式欧拉法的数值流${\mathrm{\Phi }}_{\mathrm{e}\mathrm{E},h}:z_k\mapsto z_{k+1}$不保面积；即

$\det \frac{\partial }{\partial (q_0,p_0)}{\mathrm{\Phi }}_{\mathrm{e}\mathrm{E},h}(z_0)=\left|\begin{array}{cc}1& h\\ -h\cos q_0& 1\end{array}\right|=1+h^2\cos q_0.$

隐式欧拉法也可进行类似计算，行列式为

$\det \frac{\partial }{\partial (q_0,p_0)}{\mathrm{\Phi }}_{\mathrm{i}\mathrm{E},h}(z_0)=(1+h^2\cos q_1{)}^{-1}.$

辛欧拉法是保面积的：

$\left(\begin{array}{cc}1& -h\\ 0& 1\end{array}\right)\frac{\partial }{\partial (q_0,p_0)}{\mathrm{\Phi }}_{\mathrm{s}\mathrm{E},h}(z_0)=\left(\begin{array}{cc}1& 0\\ -h\cos q_0& 1\end{array}\right),$

于是$\det (\partial {\mathrm{\Phi }}_{\mathrm{s}\mathrm{E},h}/\partial (q_0,p_0))=1$。隐式中点法则具有类似的几何特性。

总结：单摆例表明，除显式、隐式欧拉法不是解决问题的好方法外，辛欧拉法和隐式中点法则与系统的精确流非常吻合，后者更精确。而且后两种方案与精确流都保面积，是几何积分（实际上是辛积分）的两个例子。

活动标架法

活动标架法可用于构建保持ODE李对称性的数值方法。龙格-库塔法等现有方法可用活动标架法进行修改，以产生不变版本。

另见

• 能量漂移

阅读更多

• Template:Cite book
• Template:Cite book
• Template:Cite news
• Template:Cite news