Examples 127the face directly opposite. The correct handling of periodic boundary conditions within
the force calculation is described in Appendix B. Numerical calculations in periodic
boxes rarely make use of the Lennard–Jones potential in the formgiven in eqn. (3.14.3)
but rather exploit the short-range nature of the functionu(r) = 4ǫ[(σ/r)^12 −(σ/r)^6 ]
by introducing a truncated interaction
̃u(r) =u(r)S(r), (3.14.4)whereS(r) is aswitching functionthat smoothly truncates the Lennard–Jones poten-
tial to 0 at a valuer=rc, wherercis typically chosen to be between 2.5σand 3.0σ.
A useful choice forS(r) is
S(r) =
1 r < rc−λ
1 +R^2 (2R−3) rc−λ < r≤rc
0 r > rc(3.14.5)
(Watanabe and Reinhardt, 1990), whereR= [r−(rc−λ)]/λ. The parameterλis
called thehealing lengthof the switching function. This switching function has two
continuous derivatives, thus ensuring that the forces, which require ̃u′(r) =u′(r)S(r)+
u(r)S′(r), are continuous.
It is important to note several crucial differences between a simplesystem such
as the harmonic oscillator and a highly complex system such as the Lennard–Jones
(LJ) fluid. First, the LJ fluid is an example of a system that is highlychaotic. A
key characteristic of a chaotic system is known assensitive dependence on initial
conditions. That is, two trajectories in phase space with only a minute difference
in their initial conditions will diverge exponentially in time. In order to illustrate this
fact, consider two trajectories for the Lennard–Jones potential whose initial conditions
differ in just a single particle. In one of the trajectories, the initial position of a
randomly chosen particle is different from that in the other by only 10−^10 %. In this
simulation, the Lennard–Jones parameters corresponding to fluidargon (ǫ= 119. 8
Kelvin,σ= 3. 405 ̊A,m= 39.948 a.u.). Each system containsN = 864 particles
in a cubic box of volumeV = 42811. 0867 ̊A^3 , corresponding to a density of 1.34
g/cm^3. The equations of motion are integrated with a time step of 5.0 fs using a cutoff
ofrc= 2. 5 σ. The value of the total Hamiltonian is approximately 0.65 Hartrees or
7.5× 10 −^4 Hartrees/atom. The average temperature over each run is approximately
227 K. The thermodynamic parameters such as temperature and density, as well as
the time step, can also be expressed in terms of the so-called Lennard–Jonesreduced
units, in which combinations ofm,σ, andǫare multiplied by quantities such as number
density (ρ=N/V), temperature and time step to yield dimensionless versions of these.
Thus, the reduced density, denotedρ∗, is given in terms ofρbyρ∗=ρσ^3. The reduced
temperature,T∗, isT∗=T/ǫ, and the reduced time step ∆t∗= ∆t
√
ǫ/mσ^2. For the
fluid argon parameters above, we findρ= 0. 02 ̊A−^3 andρ∗= 0.8,T∗ = 1.9, and
∆t∗= 2. 3 × 10 −^3.
Fig. 3.7 shows theyposition of this particle (particle 1 in this case) in both
trajectories as functions of time when integrated numerically usingthe velocity Verlet
algorithm. Note that the trajectories follow each other closely foran initial period
but then begin to diverge. Soon, the trajectories do not resembleeach other at all.