最好的方法是使用保证保持初始条件范数的 ODE 求解器，即$\|y_n\| = \|y_0\|$对于所有。这样的求解器存在，并且被称为几何积分器，因为它们保留了精确解的几何特性（在这种情况下，能量是守恒的，即） . （一个特殊的类别是哈密​​顿力学中 ODE 系统的辛积分器。）$n\in\mathbb{N}$$\frac{d}{dt}\|y(t)\| = 0$

最简单的这种积分器（并非特定于哈密顿系统）实际上是一种非常标准的方法，分别称为梯形法、隐式中点法或Crank-Nicolson方法。在你的情况下，它相当于在每一步求解线性系统 其中是单位矩阵，是当前时间步长。这种方法是隐式的，但二阶准确且 A 稳定。（我认为没有一阶几何积分器。）我不知道它是否包含在 ODEPACK 中（scipy 的基础

$$\left(I+\frac{i\,h_n}2 H(t_{n+1})\right) y_{n+1} = \left(I-\frac{i\,h_n}2 H(t_{n})\right) y_n,$$ $I$$h_n = t_{n+1}-t_n$odeint例程），但自己实现应该不会太难。

如果您想了解有关几何积分的更多信息：

• 关于几何积分最容易理解的讨论可能是A. Iserles 的第 5 章，A First Course in the Numerical Analysis of Differential Equations，第 2 版，剑桥大学出版社，2008 年。

• “圣经”当然是E. Hairer、G. Wanner、C. Lubich，Geometric Numerical Integration，第 2 版，Springer，2006 年（有关量子力学的应用，请参见第 VII.6 章）。Acta Numerica 12, 2003, p.中还有一篇调查文章。399-450（这里是预印本），虽然这似乎集中在二阶 ODE。

• 具体而言，薛定谔方程的应用在E. Faou，几何数值积分和薛定谔方程，苏黎世高等数学讲座，EMS，2012中进行了讨论。