Lennard Jones interaction between particles. Particles moving to one point

Question

import numpy as np
import random
import pygame
background_colour = (255,255,255)
width, height = 300, 325
eps = 1
sigma = 1
dt = 0.05

class Particle():
    def __init__(self):
        self.x = random.uniform(0,400)
        self.y = random.uniform(0,500)
        self.vx = random.uniform(-.1,.1)
        self.vy = random.uniform(-.1,.1)
        self.fx = 0
        self.fy = 0
        self.m = 1
        self.size = 10
        self.colour = (0, 0, 255)
        self.thickness = 0

    def bounce(self):
        if self.x > width - self.size:
            self.x = 2*(width - self.size) - self.x

        elif self.x < self.size:
            self.x = 2*self.size - self.x

        if self.y > height - self.size:
            self.y = 2*(height - self.size) - self.y

        elif self.y < self.size:
            self.y = 2*self.size - self.y
    def getForce(self, p2):
        dx = self.x - p2.x
        dy = self.y - p2.y
        self.fx = 500*(-8*eps*((3*sigma**6*dx/(dx**2+dy**2)**4 - 6*sigma**12*dx/(dx**2+dy**2)**7)))
        self.fy = 500*(-8*eps*((3*sigma**6*dy/(dx**2+dy**2)**4 - 6*sigma**12*dy/(dx**2+dy**2)**7)))
        return self.fx, self.fy

    def verletUpdate(self,dt):
        self.x = self.x + dt*self.vx+0.5*dt**2*self.fx/self.m
        self.y = self.y + dt*self.vy+0.5*dt**2*self.fy/self.m
    def display(self):
       pygame.draw.circle(screen, self.colour, (int(self.x), int(self.y)), self.size, self.thickness)
screen = pygame.display.set_mode((width, height))
screen.fill(background_colour) 

partList = []
for k in range(10):
    partList.append(Particle())

running = True
while running:
    for event in pygame.event.get():
        if event.type == pygame.QUIT:
            running = False
    screen.fill(background_colour)

    for k, particle in enumerate(partList):
        for p2 in partList[k+1:]:
            particle.getForce(p2)
        particle.verletUpdate(dt)
        particle.bounce()
        particle.display()

    pygame.display.flip()
pygame.quit()

Is my code correct? I tried to simulate particles in 2D move with Lennard Jones forces. I think calculating forces works okay but why my particles are moving to one point? Ocasionally I also get error OverflowError: Python int too large to convert to C long Any advice would be useful.

No, your code is not correct. Your force is not the force resulting from the Lennart-Jones *potential*. You do not accumulate the interaction forces. Your verlet function does not implement Verlet since it does not update the velocities. You can not mix interaction computation and position/velocity updates, you will need 3 or 4 separate loops over the particles. — Lutz Lehmann, Apr 22 '15 at 14:14
see https://stackoverflow.com/a/29379415/3088138 and https://stackoverflow.com/a/29796633/3088138 — Lutz Lehmann, Apr 22 '15 at 14:41

shaunakde · Answer 1 · 2015-04-22T14:40:27.913

I can not comment on the physics of the simulation, but as far as the display is concerned following are my observations:

Your particles move to one point because the update condition for the x and y parameter in your code in verletUpdate are slowly moving to values beyond the display area. Also to values out of the range of the int() function which is causing your error. You can see this with the statement:

def verletUpdate(self,dt):
    self.x = self.x + dt*self.vx+0.5*dt**2*self.fx/self.m
    self.y = self.y + dt*self.vy+0.5*dt**2*self.fy/self.m
    print self.x
    print self.y

Sample Output:

290.034892392
9.98686293664
290.028208837
9.99352484332
-2.55451579742e+19
1.12437640586e+19

Also they saturate and with iterations, the update gets smaller and smaller:

def display(self):
          print ' %s + %s '%(self.x,self.y)
          pygame.draw.circle(screen, self.colour, (int(self.x), int(self.y)), self.size, self.thickness)

Output:

 10.0009120033 + 10.0042647307 
 10.0009163718 + 10.0000322065 
 10.0009120033 + 10.0042647307 
 10.0009163718 + 10.0000322065 
 ...
 10.0009163718 + 10.0000322065 
 10.0009120033 + 10.0042647307 
 10.0009163718 + 10.0000322065

This is also why your bounce functions and your limit checking is not working. And after a lot of iterations on occasion your self.x and self.y are far exceeding the limits of int().

The code seems fine, but you can get rid of the overflow error by adding some checks above the draw line. For instance I initialized them randomly again to simulate a particle going off screen and us tracking a new one. Feel free to change it.

   def display(self):
          if(self.x<0 or self.x>height):
              self.__init__()
              print "reset"
          if(self.y<0 or self.y>width):
              self.__init__()
              print "reset"
          print ' %s + %s '%(self.x,self.y)
          pygame.draw.circle(screen, self.colour, (int(self.x), int(self.y)), self.size, self.thickness)

Also at one point you adress the array as [k+1:], and addressing the zero element caused a divide by zero error. You might want to look at that.

Lennard Jones interaction between particles. Particles moving to one point

1 Answers1