Tutorial 6.4: MD Simulation on Two Atoms of Lennard-Jonesium

INTRODUCTION:

"Lennard-Jonesium" is an idealized model for the interaction of two uncharged atoms. The interaction potential is described by the Lennard-Jones 6-12 potential. This potential combines the weakly attractive London dispersion forces and the strongly repulsive exchange (atomic overlap) forces.

This tutorial carries out an MD simulation on two atoms of Lennard-Jonesium, restricted to move in only one dimension. This will introduce you to the dynamic behavior of non-covalently-linked atoms, and to the Lennard-Jones potential.

For two atoms separated by a distance r,

In this exercise, the ideal distance is r0 = 3Å, and the well depth is E0 = 0.1 kJ/mol.

As seen in the graph, the negative r-6 potential is weakly attractive, and it dominates when r > r0, the ideal contact distance. The positive r-12 term is strongly repulsive and dominates when r < r0.

GETTING READY:

Create a new folder for your work on this tutorial.

Download the following files and put them in that folder:

PDB File: LJ2.pdb

FF Mod File: LJ_E0_0.1.frcmod

Recipe...: LJ2_recipe.json

PREDICTION:

An MD simulation mimics the behavior of real masses, in this case two atoms acting under the influence of the Lennard-Jones force between them. At the start of the simulation (t = 0), the two atoms are at rest on the x-axis, separated by a distance of 8Å. Their motion is confined to that single dimension: throughout the simulation, y = 0 and z = 0. What do you think they will do as time passes? To predict the behavior, just ask yourself, "What would two real masses do if they were sitting, at rest, on the Lennard-Jones energy surface shown above, and then released?" You are to express your prediction by describing the trajectory of the system, drawing a graph showing r(t).

What's important here is for you to think about the system for a while and make a guess about its behavior. It's not important that your prediction be right. You will learn more from this simulation if you make a prediction – even if it's wrong – than you will if you simply run it without thinking about what might happen.

The simplest way to think about this system is to imagine one atom as a ball sitting on the van der Waals energy surface produced by the other atom, as in the figure. The position of the single ball (r) represents the distance between the two atoms. (This is an example of a simple coordinate transformation. In this case the motions of the two atoms, x1(t) and x2(t) are reduced to a single coordinate, r(t). You may find it easier to think about the motion of this "ball" in a fixed potential than it is to think about the motions of the two atoms.)

Gravity pulls the ball downward. If you answer the following questions, you should be able to predict how the ball will move, that is, how the distance between the two atoms will change as time passes.

At the start of the simulation, do you expect the atoms to move, or not? In which direction? (That is, do you expect the two atoms to move toward one another, or away from one another, or to stay fixed?)

Do you expect the motion to be monotonic (in one direction, forever)? Or do you expect the motion to reverse direction? If so, how far will the atoms go before changing direction?

And what will happen then? Will the atoms return to their original positions or not?

Will the motion be periodic, or not? If it is, do you expect it to be damped by friction? Will the atoms eventually come to rest? If so, where?

Once you've made your final guess, make your prediction semi-quantitative by sketching a rough graph of r(t), that is, with time on the horizontal axis and r on the vertical axis. Again, it doesn't matter whether your guess is right or not – what's important is for you to think about the likely behavior of the system. This will help you understand the results of the simulation, once it's been run.

Throughout this tutorial, do NOT rotate the structure in the Chimera window.

This is a one-dimensional system. y=0 and z = 0 for all atoms throughout the trajectory.

If you rotate the structure, you will not be able to see motion in the x direction correctly..

INSTRUCTIONS:

1. Run an MD simulation in Sophia using the input files listed above. (If you need help setting up and running the simulation, follow instructions in this link.) The recipe should produce a simulation with

• a box length of 10 nm

• a cutoff of 1 nm

• a stepsize of 1 fs

• 10,000,000 steps (10 ns)

• an output interval of 100 steps (100 fs)

2. Create a plot of r(t), the distance between the two atoms. (You learned how to do this in the first Sophia tutorial, using Chimera to measure the distance, and the Sophia Trajectory Viewer window to plot the distance. If you need help with this, here's a link to those instructions.)

ANALYSIS:

1. Compare the plot of r(t) with the sketch you made of your predicted behavior of r(t). Write a brief description of any differences, explaining the reasoning behind your prediction, and the reason the system behaved differently than you expected.

2. Examine the time-dependence of Epot, Ekin and Etot. Explain the behavior.

3. Look at Etot in detail. How well is it conserved?

4. What is the period of the motion?

CHALLENGE:

Can you come up with an analytic expression for the period as a function of the Lennard-Jones parameters E0 and r0, and the initial position r(t=0)? This is a difficult problem, but it has a solution.

Show that the period predicted by your analytic expression is equal to the period observed in the simulation.