next up previous contents
Next: The conjugate gradient algorithm Up: The numerical techniques Previous: The Predictor-Corrector algorithm

Molecular dynamics at constant temperature

To simulate the annealing of the disrupted crystal, periodic boundary conditions are applied, the sample being coupled to an external bath held at a constant temperature [106]. This coupling can be accomplished by inserting stochastic and friction terms in the equations of motion, yielding a Langevin equation

\begin{displaymath}m \dot{v}_i = F_i - m \gamma v_i + R(t)
\end{displaymath} (8.14)

where

\begin{displaymath}\langle R_i(t)R_j(t+\tau)\rangle = 2m\gamma k_BT'\delta(\tau)
\delta_{ij}.
\end{displaymath} (8.15)

T' is the temperature of the external bath to which the system will tend,

\begin{displaymath}\tau_T = \frac{1}{2\gamma}, \end{displaymath}

the time constant of this coupling, vi is the velocity of atom i and Fi, the force applied on it. The stochastic term is expressed by R(t) and the friction one, by $\gamma$. This equation then imposes global coupling to the heat bath but also inserts local perturbations by random noise. It can be shown that the time variation of the temperature T of the system, due to stochastic coupling, can be expressed by

 \begin{displaymath}\left(\frac{dT}{dt}\right)_{bath} = 2\gamma (T' - T),
\end{displaymath} (8.16)

leading to a new equation

 \begin{displaymath}m \dot{v}_i = F_i - m \gamma v_i\left(\frac{T'}{T} - 1\right)
\end{displaymath} (8.17)

where only global coupling appears. Then, to conserve only the global coupling and minimize local disturbance, from eq. (8.17) it can be seen that velocity rescaling should be made each time-step $\Delta t$ by a factor

 \begin{displaymath}\lambda = \left[1 + \frac{\Delta t}{\tau_T}\left(\frac{T'}{T} - 1
\right)\right]^{1/2}
\end{displaymath} (8.18)

where $\frac{\Delta t}{\tau_T}(\frac{T'}{T} - 1) \ll 1$. The current kinetic temperature T of a sample with N atoms is calculated each time-step by the equation

 \begin{displaymath}\frac{3}{2}Nk_B T = \sum_i \frac{1}{2}mv_i^2.
\end{displaymath} (8.19)

Thus, this method is very convenient since no new solutions of the motion equations are needed, the rate of the energy dissipation can be easily controlled by the time constant of the coupling $\tau_T$ and the final sample temperature is initially determined by T'. But velocity rescaling would affect the speed and the position of the atoms, and therefore the structure of the crystal. Thus no realistic results can be obtained from such calculations. A way to impose global coupling to an external bath with fixed temperature and to avoid, at the same time, perturbation of the atomic motion, is to apply velocity rescaling to only few atoms of the sample that do not lie in the damaged region. These atoms are chosen from the three outermost layers, including the periodic boundaries, keeping the crystal density constant. Thus, the rescaling velocity will maintain the bulk atoms, in a global way, at a kinetic temperature defined by eq. (9.6).

For the simulation of defect formation in the diamond crystal by energetic atoms, at 0 K, the atoms in the boundaries are kept fixed in their ``pure'' diamond lattice sites, and the three outermost layers are coupled to a bath held at a temperature of 0 K, in the way explained above. Thus, energy reflection from the boundaries is minimized and the structure of the disrupted region is determined only by the motion of the displaced atoms.

Since the model is deterministic, the initial configuration and the speed of the atoms play a crucial role in the simulation and will determine the nature of the motion. In the present research, the initial configuration for the post implantation annealing is the disrupted sample obtained from the 0 K ion impact. For the ``hot'' implantation simulation, the initial configuration is that of the perfect diamond crystal. The initial velocities are chosen randomly and suit the annealing temperature according to eq. (9.6).


next up previous contents
Next: The conjugate gradient algorithm Up: The numerical techniques Previous: The Predictor-Corrector algorithm
David Saada
2000-06-22