分子动力学-粒子相互作用

Question

我试图模拟基于Lennar势的惰性气体相互作用。我遇到了这样的问题:粒子之间的相互作用显然是以离散的步骤计算出来的，所以加速度往往会跳到不合理的高值，或者当两个粒子相互接近时不会有明显的变化。我使用的检测算法如下所示：

def accerlation(self):
   index = 0
   for particle, particle2 in itertools.combinations(self.particles, 2):
       index += 1
       
       x_diff = (particle.r[0] - particle2.r[0])*S
       y_diff = (particle.r[1] - particle2.r[1])*S
       distance = np.sqrt(x_diff**2 + y_diff**2)
       
       if distance < criticaldistance1 and distance> criticaldistance2:
           print("Interaction") 
           particle.a[0] -= 1/(Particle().m) *(LennardJones_Force(x_diff))
           particle.a[1] -= 1/Particle().m *(LennardJones_Force(y_diff))
           print(particle.a)

(打印命令用于调试) Lennard只是一个返回(24*epsilon*sigma**6*(Distance**6 - 2*sigma**6))/Distance**13的函数，其中epsilon和sigma是常量值。

如果您在那里找不到错误，那么问题也可能出现在以下部分：

def Changeovertime(self): 
       self.Collision_Box() #elaxtic collsions with box
       for particle in self.particles:
           particle.r += self.dt * particle.v + self.dt**2*particle.a 
           particle.a_prior = particle.a 
           self.accerlation() #
           particle.a = np.zeros(2) #
           particle.v = (particle.v + self.dt/2 * (particle.a + particle.a_prior)) *self.scaling()

当我在动画中打印一对粒子的加速度时，它看起来如下：

[-3.21599405e-01 -1.05489024e-18]\
Interaction\
[-3.35299415e-14  2.61407475e-19]\
Interaction\
[-2.52825200e+31 -1.05489024e-07]\
Interaction\
[-6.70598831e-14  5.22814950e-19]\
Interaction\
[ 1.57188229e-01 -5.51566856e-19]\
Interaction\
[-6.70598831e-03  5.22814950e-08]\
Interaction\
[ 3.14376458e-01 -1.10313371e-18]\
Interaction\
[-3.35299416e-14  2.72195193e-19]\
Interaction\
[-6.70598831e-14  5.44390386e-19]\
(Note the jump from entry 2 to entry 3)

Answer 1

在不对代码的其余部分进行注释的情况下，您所看到的问题是由以下代码行造成的：

particle.a[0] -= 1/(Particle().m) *(LennardJones_Force(x_diff))
particle.a[1] -= 1/Particle().m *(LennardJones_Force(y_diff))

力和加速度是向量:你不能把它们分解成x和y组件。你必须计算力并将它投射到x和y中(如果你想的话，你可以使用加速)。

force = LennardJones_Force(distance)
particle.a[0] -= 1/(particle.m) * force * x_diff/distance
particle.a[1] -= 1/(particle..m) * force * y_diff/distance

请检查这些标志( particle和particle2应该在对面)。

你看到的巨大跳跃是因为两个粒子有一个非常相似的x或y坐标(即使它们相距很远)，而计算错误的力是巨大的。

Answer 2

物理似乎不被遵守。

验证牛顿第三定律是否已经实现。

当这些粒子相互作用时，形成一个作用-反应对，因此为了确保能量的守恒，作用在粒子 j上的力必须是作用在粒子i (方向和绝对值)上的相同的力：Fx[i] = -Fx[j] (想象i红色和j蓝色！)。

，这里有一个注意事项，：如果这些粒子的质量相等，它们的值​​是1，你可以假设加速度等于力Fx[i] = ax[i]/1。然后，剩下的就是将这些力(Fri和Frj)分解成它们的分量(Fxi、Fyi和Fzi表示Fri和Fxj，Fyj和Fzj表示Frj)，始终记住粒子j对粒子i的作用力等于粒子i对粒子j的作用力，只是它们的符号不同(不同的感官)。

在计算机实践中，你是这样做的

def acceleration():
    for i in range(0,N_particles):
        for j in range(i+1,N_particles):
            rij = np.sqrt((x[i] - x[j])**2 + (y[i] - y[j])**2 + (z[i] - z[j])**2)
            if (rij < rmaxlj):
                flj = (48*epsilon*(sigma**12/rij**13) + 24*self.epsilon*(sigma**6/rij**7))/rij**2
                Fx[i] += flj*(x[i] - x[j])
                Fy[i] += flj*(y[i] - y[j])
                Fz[i] += flj*(z[i] - z[j])
                Fx[j] -= flj*(x[i] - x[j])
                Fy[j] -= flj*(y[i] - y[j])
                Fz[j] -= flj*(z[i] - z[j])
    return Fx, Fy, Fz

注意:所有添加到粒子j i**.的每个组件的力都从粒子的每个组件中“移除”。这实际上是牛顿的第三定律！**

这些rmaxlj criticaldistance 工作起来就像是在使用 0 ，也可以使用 i 作为0和E 248 j E 152作为E 253 1。

测试这个，稍后告诉我！