MD -моделирование с использованием Velocity Verlet в Python ⇐ Python

[Anonymous]

Я пытаюсь реализовать простую моделирование MD в Python (я новичок в этом), я использую уравнения LJ потенциал и силу вместе с методом Verlet: < /p>
def LJ_VF(r):
#r = distance in Å
#Returns V in (eV) and F in (eV/Å)
V = 4 * epsilon * ( (sigma/r)**(12) - (sigma/r)**6 )
F = 24 * epsilon * ( 2 * ((sigma**12)/(r**(13))) - ( (sigma**6)/(r**7) ))
return V , F

def velocity_verlet(x, v, f_old, f_new): #setting m=1 so that a=f
x_new = x + v * dt + 0.5 * f_old * dt**2
v_new = v + 0.5 * (f_old + f_new) * dt
return x_new, v_new
< /code>
Теперь, чтобы убедиться, что он работает, я пытаюсь использовать это в самой простой системе, то есть 2 атомах, разделенных расстоянием 4 Å: < /p>
#Constants

epsilon = 0.0103
sigma = 3.4
m = 1.0
t0 = 0.0
v0 = 0.0
dt = 0.1
N = 1000

def simulate_two_atoms(p1_x0, p1_v0, p2_x0, p2_v0):
p1_x, p2_x = [p1_x0], [p2_x0]
p1_v, p2_v = [p1_v0], [p2_v0]
p1_F, p1_V, p2_F, p2_V = [], [], [], []

r = abs(p2_x0 - p1_x0)
V, F = LJ_VF(r)
p1_F.append(F)
p1_V.append(V)
p2_F.append(-F)
p2_V.append(V)

for i in range(N - 1):
r_new = abs(p2_x[-1] - p1_x[-1])
V_new, F_new = LJ_VF(r_new)
x1_new, v1_new = velocity_verlet(p1_x[-1], p1_v[-1], p1_F[-1], F_new)
x2_new, v2_new = velocity_verlet(p2_x[-1], p2_v[-1], p2_F[-1], -F_new)

p1_x.append(x1_new)
p1_v.append(v1_new)
p2_x.append(x2_new)
p2_v.append(v2_new)

p1_F.append(F_new)
p2_F.append(-F_new)

p1_V.append(V_new)
p2_V.append(V_new)

return np.array(p1_x), np.array(p1_v), np.array(p2_x), np.array(p2_v)

#Initial conditions

p1_x0 = 0.0
p1_v0 = 0.0
p2_x0 = 4.0
p2_v0 = 0.0

#Plot

p1_x, p1_v, p2_x, p2_v = simulate_two_atoms(p1_x0, p1_v0, p2_x0, p2_v0)
time = np.arange(N)

plt.plot(time, p1_x, label="Particle 1", color="blue")
plt.plot(time, p2_x, label="Particle 2", color="green")
plt.xlabel("Time (t)")
plt.ylabel("x (Å)")
plt.title("Particle Positions Over Time (Bouncing Test)")
plt.legend()
plt.grid(True)
plt.show()
< /code>
Но я получаю неправильные значения, и график показывает, что атомы вообще не отскакивают, в отличие от друг друга! < /p>

Я пытался найти, где я Прошел не так уже очень долго, но не продвигаясь ни на каком уровне! Думал, что второй глаз может помочь! Есть советы/советы?

Подробнее здесь: https://stackoverflow.com/questions/794 ... -in-python