Molecular dynamics thesis

2.3 Implementation on a computer
2.3.1 The system in a box
CPU power today is still insufficient to simulate the length and time scales of real
Ion Beam Assisted Deposition (IBAD) experiments. Computers used to calculate the results
in this thesis need a few months to simulate 10000 atoms for 10
-8
s. Ten thousand atoms
form a cube only about 50 Å big. In order to give the impression of a large specimen using
only a small number of atoms, a simulation box with periodic boundary conditions is used.
This means that the system is considered to repeat itself after a certain length (so atoms at
position R are also present at positions R + L
x
, R - L
y
, R + 2L
x
 - 3L
y
etc. where L
x
and
L
y
are the length vectors of the box edges). This effectively removes unphysical boundary
effects, see fig. 5. On the other hand, the periodic lengths impose restrictions on the
phenomena that can be studied. For instance, nucleation from a melt will never occur if the
volume of the system is smaller than the critical nucleus size. Apart from the limited size
there is also the short simulation time. These and other limitations are discussed in section
2.5.3.
In the direction perpendicular to the film periodic boundary conditions are useless
because depositing atoms requires a free surface, so it appears as if only thin films not
attached to bulk material can be studied. This would in turn cause two problems: heat from
incoming particles cannot diffuse out of the system into a substrate and the momentum from
incoming atoms would eventually cause the film to drift away. For this reason atoms in the
two lowest atomic planes are given a special treatment in the simulations to overcome these
problems [8]. These ‘bottom atoms’ experience a force not only from the interactions with
other atoms but also from an artificial harmonic force which ties these atoms to their
specified equilibrium positions. This harmonic force prevents the film from drifting away.
The bottom atoms can also be used to remove excess energy. This is done by scaling their

---

17
Figure 5. Schematic representation of periodic boundary conditions. The upper black
atom, shown with its interaction range, interacts with all atoms in the hatched area,
but also with the atoms in the shaded areas through their periodic images in the
adjoining boxes. The computer calculates all interactions only once, of course.
velocities in such a way that a preset temperature T
bot
can be maintained and by temporarily
disabling the harmonic force on atoms whose velocity perpendicular to the film is directed
towards the film and whose kinetic energy exceeds 1/2*kT
bot
. The first selection
criterion ensures that the harmonic force is disabled only if atoms are moving towards the
film, ensuring that the bottom atoms can never be detached from the film. The second
selects those atoms that are ‘too hot’. In effect, these atoms are given a lower potential
energy without being accelerated by the harmonic force, thus removing energy from the
system.
The non-periodicity perpendicular to the film allows for the deposition of thicker
films using the ‘cut and shift down’ method. This method works as follows: deposition in a
simulation box is continued until the box is nearly full. Then the film is cooled to near 0 K
and the bottom part of the film is cut away. The rest of the film is shifted down. Atoms in
what are now the lowest two planes become bottom atoms. Because the film has been
cooled to near 0 K their equilibrium positions are known. The system is reheated to its
original temperature and deposition is continued. This method can be applied when
diffusion in the lower part of the film is negligible (which it usually is, see section 4.5.2).
Atoms can be introduced at the top of the box, with specified energy and angles of
incidence. The positions for introduction are selected by a random number generator.
Particles that move out of the box through the top (or sometimes through the bottom) are
assumed no longer to be of interest to the simulation and are removed from the system. A

---

18
record is kept containing information about the angles, energies, positions, and time at
which particles are introduced and removed.
2.3.2 Numerical algorithm
Eqn. (1) can usually not be solved analytically for a system of thousands of atoms.
Therefore the widely used Velocity-Verlet algorithm [9] is used to produce a numerical
approximation. From atomic positions r, velocities v , and accelerations a  at time t  the
algorithm first calculates the new positions at time t + 
∆
t
,