查看“︁几何积分”︁的源代码

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

==以摆为例==
可考虑[[单摆]]运动以引出几何积分的研究。

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

:<math>H(q,p) = T(p)+U(q) = \frac{1}{2}p^2 - \cos q, </math>

其给出[[哈密顿方程]]

:<math>(\dot q,\dot p) = (\partial H / \partial p, -\partial H / \partial q ) = (p,-\sin q). \, </math>

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

: <math> z_{k+1} = z_k + hf(z_k) \, </math>（[[欧拉法|显式欧拉]]）；

: <math> z_{k+1} = z_k + hf(z_{k+1}) \, </math>（隐式欧拉）；

: <math> z_{k+1} = z_k + hf(q_k,p_{k+1}) \, </math>（辛欧拉）；

: <math> z_{k+1} = z_k + hf((z_{k+1}+z_k)/2) \, </math>（隐式中点法则）。

（注意，辛欧拉法用显式欧拉法处理''q''，用隐式欧拉法处理<math>p</math>。）

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

[[Image:pendulumtrajectories.png|thumb|300px|单摆轨迹]]
显式（或隐式）欧拉法是从原点向外（或向内）的螺旋运动。另两种方法显示了正确的定性行为，隐式中点法则与精确解的吻合程度高于辛欧拉法。

回顾一下，具有1自由度的哈密顿系统的精确流<math>\phi_t</math>是保面积的，即
:<math>\det\frac{\partial\phi_t}{\partial (q_0,p_0)} = 1</math> for all <math>t</math>.
此式很容易手动验证。对我们的单摆例子，可以发现，显式欧拉法的数值流<math>\Phi_{{\mathrm{eE}},h}:z_k\mapsto z_{k+1}</math>'''不'''保面积；即

:<math>\det\frac{\partial}{\partial (q_0,p_0)}\Phi_{{\mathrm{eE}},h}(z_0)
   = \begin{vmatrix}1&h\\-h\cos q_0&1\end{vmatrix}
   = 1+h^2\cos q_0.</math>

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

:<math>\det\frac{\partial}{\partial (q_0,p_0)}\Phi_{{\mathrm{iE}},h}(z_0)
    = (1+h^2\cos q_1)^{-1}.</math>

辛欧拉法'''是'''保面积的：

:<math>
    \begin{pmatrix}1&-h\\0&1\end{pmatrix}\frac{\partial}{\partial (q_0,p_0)}\Phi_{{\mathrm{sE}},h}(z_0)
    = \begin{pmatrix}1&0\\-h\cos q_0&1\end{pmatrix},</math>

于是<math>\det(\partial\Phi_{{\mathrm{sE}},h}/\partial (q_0,p_0)) = 1</math>。隐式中点法则具有类似的几何特性。

总结：单摆例表明，除显式、隐式欧拉法不是解决问题的好方法外，辛欧拉法和隐式中点法则与系统的精确流非常吻合，后者更精确。而且后两种方案与精确流都保面积，是几何积分（实际上是[[辛积分]]）的两个例子。
==活动标架法==
[[活动标架]]法可用于构建保持ODE[[李群|李]][[对称群|对称性]]的数值方法。[[龙格-库塔法]]等现有方法可用活动标架法进行修改，以产生不变版本。<ref>Pilwon Kim (2006),  " [http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.527.93&rep=rep1&type=pdf Invariantization of Numerical Schemes Using Moving Frames] {{Wayback|url=http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.527.93&rep=rep1&type=pdf |date=20170202162211 }}"</ref>

==另见==
* [[能量漂移]]

==参考文献==
{{Reflist}}

==阅读更多==
* {{cite book
  | first1 = Ernst | last1 = Hairer
  | first2 = Christian | last2 = Lubich
  | first3 = Gerhard | last3 = Wanner
  | title = Geometric Numerical Integration: Structure-Preserving Algorithms for Ordinary Differential Equations
  | year = 2002
  | publisher = Springer-Verlag
  | isbn = 3-540-43003-2
  }} 
* {{cite book
  | last1 = Leimkuhler | first1 = Ben
  | last2 = Reich | first2 = Sebastian
  | title = Simulating Hamiltonian Dynamics
  | year = 2005
  | publisher = Cambridge University Press
  | isbn = 0-521-77290-7
  }}
* {{cite news
  | first1 = C.J. | last1 = Budd
  | first2 = M.D. | last2 = Piggott
  | title = Geometric Integration and its Applications
  | series = Handbook of Numerical Analysis
  | volume = 11
  | pages = 35–139
  | year = 2003
  | publisher = Elsevier
  | doi = 10.1016/S1570-8659(02)11002-7
  | isbn = 9780444512475
 }} 
* {{cite news
  | first1 = Pilwon | last1 = Kim
  | title = Invariantization of Numerical Schemes Using Moving Frames
  | series = BIT Numerical Mathematics
  | volume = 47
  | pages = 525–546
  | year = 2007
  | issue = 3
 | publisher = Springer
  | doi = 10.1007/s10543-007-0138-8
  }}

[[Category:数值微分方程]]