蛙跳算法相空间面积守恒疑问：库仑力下面积未保持恒定

根据相关讲义，蛙跳积分算法应在单体系统的相空间中保持面积。我针对库仑力实现了该算法，能量符合预期守恒，但由四个不同初始条件定义的平行四边形面积并未如预期保持恒定。

补充说明：若针对简谐振荡器运行该算法，面积是守恒的。下图展示了面积随时间的变化情况，有趣的是，使用欧拉法得到的曲线完全相同。

实现代码如下：

DT = 0.1
T = 20

def inverse_square_law(x, k=1):
    return k / (x ** 2)

def KE(v, m=1):
    return 0.5 * m * v**2

def PE(x, k=1):
    return k / x

def area(a, b):
    return np.abs(a[0] * b[1] - a[1] * b[0])

x0 = np.array([2, 2, 1, 1])
v0 = np.array([1, 0, 1, 0])

x = x0
v = v0

times = []
ke_list = []
pe_list = []
total_energy_list = []
all_positions = [[] for _ in range(len(x))]
all_velocities = [[] for _ in range(len(v))]
areas = []

v = v + 0.5 * DT * inverse_square_law(x)

for t in np.arange(0, T, DT):
    x = x + DT * v
    a = inverse_square_law(x)
    v = v + DT * a

    v_full = v - 0.5 * DT * a
    ke = KE(v_full[0])
    pe = PE(x[0])
    total_energy = ke + pe

    times.append(t)
    ke_list.append(ke)
    pe_list.append(pe)
    total_energy_list.append(total_energy)

    for i in range(len(x)):
        all_positions[i].append(x[i])
        all_velocities[i].append(v_full[i])

    p0 = np.array([x[0], v_full[0]])
    p1 = np.array([x[1], v_full[1]])
    p2 = np.array([x[2], v_full[2]])
    p3 = np.array([x[3], v_full[3]])

    area_val = area(p2 - p0, p1 - p0)
    areas.append(area_val)

问题根源 #

1. 对辛积分特性的理解偏差
   蛙跳积分是辛积分算法，它保证的是无穷小相空间体积元的辛结构守恒，而非你定义的由四个有限初始点构成的平行四边形面积守恒。

2. 线性与非线性系统的差异
   简谐振荡器的哈密顿量是二次型，其相流是线性变换，有限体积的平行四边形会被线性映射为同面积的平行四边形，因此蛙跳积分能保持这个有限面积。但库仑力的哈密顿量是非线性的，相流是非线性变换，有限区域的形状会被扭曲，面积自然不会保持恒定——这是系统本身的特性，和蛙跳积分的实现无关。

3. 有限区域面积的本质
   你用四个离散初始点的演化计算平行四边形面积，本质是观察有限区域的体积变化。辛积分只保证无穷小体积元的辛测度守恒，对有限区域，非线性相流会拉伸或压缩体积，只是整体的辛测度总和不变，但单个有限区域的面积可以变化。

4. 欧拉法的巧合
   你提到欧拉法得到的曲线和蛙跳一致，这是因为你测试的是有限区域的面积变化，在这个有限尺度下，两种算法的表现碰巧相似。但欧拉法不是辛积分，它连无穷小体积元的守恒都无法保证，只是你的测试场景下，有限面积的变化趋势刚好重合。

另外，你的代码实现是正确的：v_full的计算符合蛙跳积分的半步速度要求，能量守恒也验证了这一点，问题核心是对“相空间面积守恒”的概念误解。

内容的提问来源于stack exchange，提问作者AJB