Difference between revisions of "CompSciHW9"
From Predictive Chemistry
(Created page with "Introduction to Scientific Computing, HW 9. Due Friday, Apr. 15, 2016. <pre> 1) Write a complete code to simulate the Lennard-Jones gas in a 2D periodic box with box leng...") |
|||
Line 32: | Line 32: | ||
r = closest_distance(z, x) # matrix of closest distances (n x 3) |
r = closest_distance(z, x) # matrix of closest distances (n x 3) |
||
r2 = sum(r*r, 1) # vector of squared distances (n) |
r2 = sum(r*r, 1) # vector of squared distances (n) |
||
− | u = r2**3 |
+ | u = r2**-3 |
print 6 * sum(r * ((2*u*u - u)/sqrt(r2))[:,newaxis], 0) # sum over other atoms |
print 6 * sum(r * ((2*u*u - u)/sqrt(r2))[:,newaxis], 0) # sum over other atoms |
||
</source> |
</source> |
Revision as of 10:06, 4 April 2016
Introduction to Scientific Computing, HW 9. Due Friday, Apr. 15, 2016.
1) Write a complete code to simulate the Lennard-Jones gas in a 2D periodic box with box length L = 14 and n=100 particles. Start them off on a 10x10 grid with Gaussian distributed velocities. Ignore units and assume beta = m = 1. The Hamiltonian is given by H = sum_j m v_j^2/2 + 1/2 sum_{i != j} u_ij^2 - u_ij where u_ij = |x_i - x_j|**-6 The force on each particle, i, is therefore F_i = sum_{j != i} 6 (x_i - x_j) / |x_i - x_j| ( 2 u_ij^2 - u_ij ) 2) Run the simulation for 100 steps, and create a plot showing the locations of the atoms every 10 steps. 3) For every timestep, calculate the kinetic and potential energies. What do you observe about the behavior of the potential energy? 4) Make a plot of the total energy vs. time for your 100 step simulation. Overlay these plots for several different values of the numerical timestep, dt.
Hints: <source lang="python">
- Wrap all coordinates in an array to the range [0,L)
print x - L*floor(x/L)
- Find the closest distance between two points, r_ij
y = x[i] - x[j] print y - L*floor(y/L + 0.5)
- Calculate the LJ force on an atom at point z
r = closest_distance(z, x) # matrix of closest distances (n x 3) r2 = sum(r*r, 1) # vector of squared distances (n) u = r2**-3 print 6 * sum(r * ((2*u*u - u)/sqrt(r2))[:,newaxis], 0) # sum over other atoms </source>