Jacques G. Amar

The Monte Carlo method in science and engineering

Since 1953, researchers have applied the Monte Carlo method to a wide range of areas. Specialized algorithms have also been developed to extend the methods applicability and efficiency. The author describes some of the algorithms that have been developed to perform Monte Carlo simulations in science and engineering

Chapter 4 Accelerated Molecular Dynamics Methods: Introduction and Recent Developments

Abstract Because of its unrivaled predictive power, the molecular dynamics (MD) method is widely used in theoretical chemistry, physics, biology, materials science, and engineering. However, due to computational cost, MD simulations can only be used to directly simulate dynamical processes over limited timescales (e.g., nanoseconds or at most a few microseconds), even though the simulation of nonequilibrium processes can often require significantly longer timescales, especially when they involve thermal activation. In this paper, we present an introduction to accelerated molecular dynamics, a class of methods aimed at extending the timescale range of molecular dynamics, sometimes up to seconds or more. The theoretical foundations underpinning the different methods (parallel replica dynamics, hyperdynamics, and temperature-accelerated dynamics) are first discussed. We then discuss some applications and recent advances, including super-state parallel replica dynamics, self-learning hyperdynamics, and spatially parallel temperature-accelerated dynamics.

Hybrid asynchronous algorithm for parallel kinetic Monte Carlo simulations of thin film growth

We have generalized and implemented the hybrid asynchronous algorithm, originally proposed for parallel simulations of the spin-flip Ising model, in order to carry out parallel kinetic Monte Carlo (KMC) simulations. The parallel performance has been tested using a simple model of thin-film growth in both 1D and 2D. We also briefly describe how the data collection must be modified as compared to the case of the spin-flip Ising model in order to carry out rigorous data collection. Due to the presence of a wide range of rates in the simulations, this algorithm turns out to be very inefficient. The poor parallel performance results from three factors: (1) the high probability of selecting a Metropolis Monte Carlo (MMC) move, (2) the low acceptance probability of boundary moves and (3) the high cost of communications which is required before every MMC move. We also find that the parallel efficiency in two dimensions is lower than in one-dimension due to the higher probability of selecting an MMC attempt, suggesting that this algorithm may not be suitable for KMC simulations of two-dimensional thin-film growth.

Adaptive temperature-accelerated dynamics.

We present three adaptive methods for optimizing the high temperature T(high) on-the-fly in temperature-accelerated dynamics (TAD) simulations. In all three methods, the high temperature is adjusted periodically in order to maximize the performance. While in the first two methods the adjustment depends on the number of observed events, the third method depends on the minimum activation barrier observed so far and requires an a priori knowledge of the optimal high temperature T(high)(opt)(E(a)) as a function of the activation barrier E(a) for each accepted event. In order to determine the functional form of T(high)(opt)(E(a)), we have carried out extensive simulations of submonolayer annealing on the (100) surface for a variety of metals (Ag, Cu, Ni, Pd, and Au). While the results for all five metals are different, when they are scaled with the melting temperature T(m), we find that they all lie on a single scaling curve. Similar results have also been obtained for (111) surfaces although in this case the scaling function is slightly different. In order to test the performance of all three methods, we have also carried out adaptive TAD simulations of Ag/Ag(100) annealing and growth at T = 80 K and compared with fixed high-temperature TAD simulations for different values of T(high). We find that the performance of all three adaptive methods is typically as good as or better than that obtained in fixed high-temperature TAD simulations carried out using the effective optimal fixed high temperature. In addition, we find that the final high temperatures obtained in our adaptive TAD simulations are very close to our results for T(high)(opt)(E(a)). The applicability of the adaptive methods to a variety of TAD simulations is also briefly discussed.

Nuclear size corrections to the energy levels of single-electron and -muon atoms

We formulate an analytic method which accounts for the finite size of the nucleus by treating it as a boundary value problem. The method is used to obtain solutions of the Dirac equation for a central potential that is proportional to 1/r only for values of the radial coordinate greater than a given value R. Our results are applied to a non-perturbative calculation of the nuclear size corrections to the energy levels of single-electron and single-muon atoms. For values of the nuclear charge number Z greater than 40 in the case of electronic atoms, and greater than 1 in the case of muonic atoms, we find large discrepancies between our results for the atomic energy levels and those obtained from first-order relativistic perturbation theory.

Kinetics of submonolayer epitaxial growth

Abstract Molecular beam epitaxy is an important method for growing thin-films and nanostructures. One of the scientific challenges is to understand the fundamental processes that control the evolution of thin film structure and morphology. The results of kinetic Monte Carlo simulations carried out to study the dependence of the submonolayer scaled island-size distribution on the critical island-size are presented and compared with experiments. A recently developed method which involves a self-consistent coupling of evolution equations for the capture-zone distributions with rate-equations for the island-size distribution is also described. Our method explicitly takes into account the existence of a denuded (“capture”) zone around every island and the correlations between the size of an island and the corresponding average capture zone, and has been used to develop a quantitative rate-equation approach to irreversible submonolayer growth on a two-dimensional substrate. The resulting predictions for the capture-zones, capture numbers, and island-size distributions are in excellent agreement with experimental results and kinetic Monte Carlo simulations.

Diffusion-Annihilation and the Kinetics of the Ising Model in One Dimension

The relationship between the one-dimensional kinetic Ising model at zero temperature and diffusion-annihilation in one dimension is studied. Explicit asymptotic results for the average domain size, average magnetization squared, and pair-correlation function are derived for the Ising model, for arbitrary initial magnetization. For the case of zero initial magnetization (m0=0, a number of recent exact results for diffusion-annihilation with random initial conditions are obtained. However, for the casem0 not equal to zero, the asymptotic behavior turns out to be different from diffusion-annihilation with random initial conditions and at a finite density. In addition, in contrast to the case of diffusion-annihilation, the domain-size distribution scaling functionh(x) is found to depend nontrivially on the initial magnetization. The origin of these differences is clarified and the existence of nontrivial correlations in the initial wall distribution for finite initial magnetization is found to be responsible for these differences. Results of Monte Carlo simulations for the domain size distribution function for different initial magnetizations are also presented.

Self-consistent rate-equation approach to irreversible submonolayer growth in one dimension

Abstract A self-consistent rate-equation (RE) approach to irreversible submonolayer growth in one dimension is presented. Our approach is based on a set of dynamical equations for the evolution of gaps between islands which is coupled to the island-density REs via local capture numbers and explicitly takes into account correlations between the size of an island and the corresponding capture zone. In the most simple formulation, fragmentation of capture zones is not directly included, but accounted for through a uniform rescaling, while nucleation is assumed to generate only gaps with average length. Using this approach, we have been able to accurately predict the scaled island-size, capture-number, and average-gap-size distributions in the pre-coalescence regime. Our approach also leads to a novel analytical expression for the monomer capture number σ 1 =(4/ RN 1 γ ) 1/2 where N 1 is the monomer density, γ is the fraction of the substrate covered by islands, and R is the ratio D / F of the diffusion rate to deposition flux which agrees with simulations over the entire pre-coalescence regime, and implies a novel scaling behavior for the island density at low coverage, in contrast to earlier predictions. Comparisons between our RE results and kinetic Monte Carlo simulations are presented for both point islands and extended islands.

Characterization of surface morphology in epitaxial growth

Abstract Simulated kinematic antiphase RHEED and HRLEED profiles are calculated along with the surface structure factor for a model of Fe Fe (100) deposition in order to clarify the interpretation of diffraction profiles in recent experiments on Fe Fe (100) growth. Similar calculations are also presented for a self-affine surface. While self-affine surfaces do not exhibit a characteristic RHEED peak, in the case of surfaces with a typical length scale the simulated RHEED profile exhibits a peak corresponding to the typical feature size, in agreement with recent experiments. The existence of this peak appears to be due to the large amount of shadowing present in low-angle RHEED, which limits the amount of destructive interference between layers. In contrast, simulated antiphase HRLEED patterns appear to approach an invariant profile for both self-affine and mound-like surface morphologies. For the case of small mounds, our results predict a HRLEED profile with a weak peak corresponding to the average terrace size which moves outward with increasing coverage, and eventually reaches an invariant form due to angle selection. The disappearance of the HRLEED peak for surfaces which have large mound structures is explained in terms of the antiphase condition and the range of variation of terrace sizes and provides an alternative explanation for the HRLEED results observed in Ref. [10].

Kinetics of epitaxial growth and roughening

Abstract We review some recent results on epitaxial growth and surface roughening. Particular emphasis is placed on the concept of the critical island size in submonolayer growth and on the existence of scaling in both the submonolayer and multilayer growth regimes. The use of scaling ideas as well as Monte Carlo simulations and continuum equations is shown to be effective in understanding experimental results for submonolayer growth and surface roughening.