如何在python中实现适用于谐振子的velocity Verlet积分器？

Question

我是python新手，我正在尝试实现一个适用于谐振子的velocity Verlet积分器。 从我下面的笔记本（摘自：http://hplgit.github.io/prog4comp/doc/pub/._p4c-solarized-Python022.html）中可以看到，Euler的方法有效，但verlet积分器无效。我错过了什么

## Euler's method

from numpy import zeros, linspace, pi, cos, array
import matplotlib.pyplot as plt
import numpy as np

omega = 1
m=1
N=500
dt=0.08
t = linspace(0, N*dt, N+1)

u = zeros(N+1)
vi = zeros(N+1)

X_0 = 2
u[0] = X_0
vi[0] = 0

for n in range(N):
    u[n+1] = u[n] + dt*vi[n]
    vi[n+1] = vi[n] - dt*omega**2*u[n]

fig = plt.figure()
l1, l2 = plt.plot(t, u, 'b-', t, X_0*cos(omega*t), 'r--')
fig.legend((l1, l2), ('numerical', 'exact'), 'upper left')
plt.xlabel('t')
plt.show()

## Velocity Verlet

for n in range(N):
    v_next = vi[n] - 0.5*omega**2*u[n]*dt
    u[n+1]= u[n]+ dt*v_next
    vi[n+1] = v_next - 0.5*omega**2*u[n]*dt
    np.append(u, x_next)

plt.plot(t, X_0*cos(omega*t), label = 'Analytical')
plt.plot(t, u, 'r--', label = 'Verlet')
plt.title('Verlet', fontweight = 'bold', fontsize = 1)
plt.xlabel('t', fontweight = 'bold', fontsize = 14)
plt.ylabel('X', fontweight = 'bold', fontsize = 14)
plt.legend()
plt.show()

Answer 1

您没有正确执行velocity Verlet步骤，第二次速度更新使用新计算的位置，而不是旧位置

vi[n+1] = v_next - 0.5*omega**2*u[n+1]*dt

这一微小变化应能恢复二阶和能量/振幅保持

---

另外，如果删除以前更改的剩余内容，则没有x_next可追加