Practical 07: Verlet Integrator in 3D

Harmonic bonds

Your task is to simulate a carbon-carbon bond.

We model it as a simple harmonic oscillator with masses 6u and force constant 1860e2 kJ/mol/nm^2 with a reference bond length of 0.15 nm. We integrate the equations of motion with the Verlet algorithm.

• We will treat the two particles as points in 3-space.
• Outline of code is provided as hoscmd_v0.py.

Task

1. Discuss the problem with your neighbors in your row. Derive the equations for the potential energy and the force (see Assignment 04). Note that this should be treated as a 3D problem: input and output are points with three cartesian coordinates.

2. Download hoscmd_v0.py and mdIO.py from http://becksteinlab.physics.asu.edu/pages/courses/2013/SimBioNano/07/

   Rename hoscmd_v0.py to hoscmd.py and work on this code. The file mdIO.py contains updated functions for trajectory output (namely mdIO.write_xyz_frame()).

3. Make the code work. See the comments in the file.

   You may work with your neighbors.

4. Produce a graph of the potential energy over 1000 steps at a time step of 0.001 ps (1 fs), initial atom separation of 0.2 nm, and initial velocities 0.

5. Visualize the resulting trajectory hosc.xyz in VMD :

   1. Open /Applications/VMD
   2. Menu: File -> New Molecule: browse to hosc.xyz and Load
   3. Menu: Graphics -> Representations:
      • Drawing Method: CPK
   4. Play trajectory: Main window lower right corner ‘>’
   5. Rotate with mouse

   You can try to make a movie with Extensions -> Visualization -> Movie Maker:

   • Movie Settings: Trajectory
   • Format: try different ones, maybe one works

6. As a bonus: also calculate the kinetic and total energy of the system and plot.