哈密顿正则方程的数值积分

zhailiangjun

各位大侠请帮忙看看。小弟正在做一个哈密顿系统的正则运动方程积分的问题，采用了matlab中的ode函数对微分方程组进行积分运算（正则运动方程是微分方程组），但是得到的解并不能保持哈密顿系统的能量守恒，多数是表现为能量逐渐耗散的形式，如图能量的随时间变化的曲线，不知道大家有没有遇到过类似的情形，还请大家多多帮忙啊。
嗯。可能我说的不太明白吧。也可以还原成这样的一个数学问题，就是说有一个满足一个约束条件的微分方程组，如何求解才能使求得的这样一组数值解尽可能的满足这么一个约束条件。举个例子，有这个一个微分方程组；
dx1/dt=x2；
dx2/dt=-x1;
满足的约束条件 E=1/2*x1^2+1/2*x2^2,E为一个常数，这是一个谐振子的哈密顿量和正则方程，根据哈密顿原理，由微分方程组解出来的随时间变化的数值解，应满足约束条件， E=1/2*x1^2+1/2*x2^2，但事实上得到的结果并非如此。我想是我的计算方法出现了问题，不知道大家有没有遇到过这样的关于约束条件的问题，有什么好的解决办法没有，请大家不吝赐教

[ 本帖最后由 ChaChing 于 2010-4-25 16:37 编辑 ]

yangzhanwen

关于这类问题不能用一般的RK方法，要用基于哈密顿系统建立的辛方法（冯康提出的）或者是李群上的数值方法。没做过具体的问题，详细的内容可以在网上查到。

zhailiangjun

多谢提醒，

Rainyboy

回复 3 # zhailiangjun 的帖子

问一下啊，你开始使用辛方法了么？如果是的话，用的是什么编程语言？

zhailiangjun

回复 4 # Rainyboy 的帖子

没有，因为我能够用龙哥库塔的方法来做一部分，暂时的回避了这些问题，但是现在又不得不面对了。我不会这个辛算法，你有辛算法的程序吗？

Rainyboy

回复 5 # zhailiangjun 的帖子

我还没到程序那一步呢，现在还在学习，呵呵。
我看的书是《力功能量与幸数学》，还有《弹性力学的幸数学方法》。

mistuning

matlab里的rk方法是不能的，你可以查查文献，有保辛的rk方法，但好像是隐式的，也有改进的保辛性质比较好的（但不完全保辛）的伪辛算法。