Kinetic Lattice Monte Carlo Simulations Research Articles

Atomistic models of In and Ga diffusion in Cu(In,Ga)Se2

To improve the performance of Cu(In,Ga)Se2 thin-film photovoltaic devices, a robust understanding of the dominant diffusion pathways of the alloy species In and Ga is needed. Here, the most probable defect complexes and mechanisms for In and Ga diffusion are identified with the aid of density functional theory. The binding energies and migration barriers for these complexes are calculated in bulk CuInSe2 and CuGaSe2. Analytic models and kinetic lattice Monte Carlo simulations are employed to predict the diffusivity of In and Ga under variations in composition and temperature. We find that a model based on coulombic interactions between group III antisites and vacancies on the Cu-sublattice produces results that match well with experiment.

Development of a parallel multiscale 3D model for thrombus growth under flow

Thrombus growth is a complex and multiscale process involving interactions spanning length scales from individual micron-sized platelets to macroscopic clots at the millimeter scale. Here, we describe a 3D multiscale framework to simulate thrombus growth under flow comprising four individually parallelized and coupled modules: a data-driven Neural Network (NN) that accounts for platelet calcium signaling, a Lattice Kinetic Monte Carlo (LKMC) simulation for tracking platelet positions, a Finite Volume Method (FVM) simulator for solving convection-diffusion-reaction equations describing agonist release and transport, and a Lattice Boltzmann (LB) flow solver for computing the blood flow field over the growing thrombus. Parallelization was achieved by developing in-house parallel routines for NN and LKMC, while the open-source libraries OpenFOAM and Palabos were used for FVM and LB, respectively. Importantly, the parallel LKMC solver utilizes particle-based parallel decomposition allowing efficient use of cores over highly heterogeneous regions of the domain. The parallelized model was validated against a reference serial version for accuracy, demonstrating comparable results for both microfluidic and stenotic arterial clotting conditions. Moreover, the parallelized framework was shown to scale essentially linearly on up to 64 cores. Overall, the parallelized multiscale framework described here is demonstrated to be a promising approach for studying single-platelet resolved thrombosis at length scales that are sufficiently large to directly simulate coronary blood vessels.

A three-dimensional multiscale model for the prediction of thrombus growth under flow with single-platelet resolution.

Modeling thrombus growth in pathological flows allows evaluation of risk under patient-specific pharmacological, hematological, and hemodynamical conditions. We have developed a 3D multiscale framework for the prediction of thrombus growth under flow on a spatially resolved surface presenting collagen and tissue factor (TF). The multiscale framework is composed of four coupled modules: a Neural Network (NN) that accounts for platelet signaling, a Lattice Kinetic Monte Carlo (LKMC) simulation for tracking platelet positions, a Finite Volume Method (FVM) simulator for solving convection-diffusion-reaction equations describing agonist release and transport, and a Lattice Boltzmann (LB) flow solver for computing the blood flow field over the growing thrombus. A reduced model of the coagulation cascade was embedded into the framework to account for TF-driven thrombin production. The 3D model was first tested against in vitro microfluidics experiments of whole blood perfusion with various antiplatelet agents targeting COX-1, P2Y1, or the IP receptor. The model was able to accurately capture the evolution and morphology of the growing thrombus. Certain problems of 2D models for thrombus growth (artifactual dendritic growth) were naturally avoided with realistic trajectories of platelets in 3D flow. The generalizability of the 3D multiscale solver enabled simulations of important clinical situations, such as cylindrical blood vessels and acute flow narrowing (stenosis). Enhanced platelet-platelet bonding at pathologically high shear rates (e.g., von Willebrand factor unfolding) was required for accurately describing thrombus growth in stenotic flows. Overall, the approach allows consideration of patient-specific platelet signaling and vascular geometry for the prediction of thrombotic episodes.

Computer modeling of single-layer nanocluster formation in a thin SiO2 layer buried in Si by ion mixing and thermal phase decomposition

A single sheet of Si nanoclusters with an average diameter of about 2 nm has been formed in a 30 nm Si/7 nm SiO2/Si layer stack by 50 and 60 keV Si+ ion-beam mixing at room temperature and fluences between 8.5 ⋅ 1015 and 2.6 ⋅ 1016 ions/cm2 and by subsequent thermal annealing at a temperature above 1000 °C. Computer modeling of the process is accomplished by TRIDYN dynamic ballistic simulation of ion mixing and subsequent lattice kinetic Monte Carlo simulation of the phase decomposition of substoichiometric silicon oxide into Si nanoclusters in a SiO2 matrix. The simulation algorithms are briefly described with special emphasis on the choice of governing parameters for the present system. In comparison to the experimental results, it is concluded that the predicted ion mixing profiles overestimate the interface broadening. This discrepancy is attributed to the neglect of chemical driving forces in connection with thermal-spike induced diffusion, which tends to reconstitute the Si/SiO2 interfaces. With a corresponding correction and a suitable number of Monte Carlo steps, the experimentally obtained areal densities and average diameters of the nanoclusters are successfully reproduced.

A Practical Guide to Surface Kinetic Monte Carlo Simulations.

This review article is intended as a practical guide for newcomers to the field of kinetic Monte Carlo (KMC) simulations, and specifically to lattice KMC simulations as prevalently used for surface and interface applications. We will provide worked out examples using the kmos code, where we highlight the central approximations made in implementing a KMC model as well as possible pitfalls. This includes the mapping of the problem onto a lattice and the derivation of rate constant expressions for various elementary processes. Example KMC models will be presented within the application areas surface diffusion, crystal growth and heterogeneous catalysis, covering both transient and steady-state kinetics as well as the preparation of various initial states of the system. We highlight the sensitivity of KMC models to the elementary processes included, as well as to possible errors in the rate constants. For catalysis models in particular, a recurrent challenge is the occurrence of processes at very different timescales, e.g., fast diffusion processes and slow chemical reactions. We demonstrate how to overcome this timescale disparity problem using recently developed acceleration algorithms. Finally, we will discuss how to account for lateral interactions between the species adsorbed to the lattice, which can play an important role in all application areas covered here.

A theory of growing crystalline nanorods – Mode I

Nanorods grow in two possible modes during physical vapor deposition (PVD). In mode I, monolayer surface steps dictate the diameter of nanorods. In mode II, multiple-layer surface steps dictate the diameter, which is the smallest possible under physical vapor deposition [5,10]. This paper reports closed-form theories of terrace lengths and nanorod diameter during the growth in mode I, as a function of deposition conditions. The accompanying lattice kinetic Monte Carlo simulations verify these theories. This study reveals that (1) quasi-steady growth exists for each set of nanorod growth conditions, and (2) the characteristic length scales, including terrace lengths and nanorod diameter at the quasi-steady state, depend on the deposition conditions – deposition rate F, substrate temperature T, and incidence angle θ – only as a function of l2D/tan θ, with l2D=2(v2DFcosθ)13 as a diffusion-limited length scale and v2D as the atomic diffusion jump rate over monolayer surface steps.

A Lattice Kinetic Monte Carlo Solver for First-Principles Microkinetic Trend Studies.

Mean-field microkinetic models in combination with Brønsted-Evans-Polanyi like scaling relations have proven highly successful in identifying catalyst materials with good or promising reactivity and selectivity. Analysis of the microkinetic model by means of lattice kinetic Monte Carlo promises a faithful description of a range of atomistic features involving short-range ordering of species in the vicinity of an active site. In this paper, we use the "fruit fly" example reaction of CO oxidation on fcc(111) transition and coinage metals to motivate and develop a lattice kinetic Monte Carlo solver suitable for the numerically challenging case of vastly disparate rate constants. As a result, we show that for the case of infinitely fast diffusion and absence of adsorbate-adsorbate interaction it is, in fact, possible to match the prediction of the mean-field-theory method and the lattice kinetic Monte Carlo method. As a corollary, we conclude that lattice kinetic Monte Carlo simulations of surface chemical reactions are most likely to provide additional insight over mean-field simulations if diffusion limitations or adsorbate-adsorbate interactions have a significant influence on the mixing of the adsorbates.

Acceleration and sensitivity analysis of lattice kinetic Monte Carlo simulations using parallel processing and rate constant rescaling.

Kinetic Monte Carlo (KMC) simulation provides insights into catalytic reactions unobtainable with either experiments or mean-field microkinetic models. Sensitivity analysis of KMC models assesses the robustness of the predictions to parametric perturbations and identifies rate determining steps in a chemical reaction network. Stiffness in the chemical reaction network, a ubiquitous feature, demands lengthy run times for KMC models and renders efficient sensitivity analysis based on the likelihood ratio method unusable. We address the challenge of efficiently conducting KMC simulations and performing accurate sensitivity analysis in systems with unknown time scales by employing two acceleration techniques: rate constant rescaling and parallel processing. We develop statistical criteria that ensure sufficient sampling of non-equilibrium steady state conditions. Our approach provides the twofold benefit of accelerating the simulation itself and enabling likelihood ratio sensitivity analysis, which provides further speedup relative to finite difference sensitivity analysis. As a result, the likelihood ratio method can be applied to real chemistry. We apply our methodology to the water-gas shift reaction on Pt(111).

Lattice kinetic Monte Carlo simulation study of the early stages of epitaxial GaN(0001) growth

The early stages of GaN(0001) epitaxial growth are modeled using a lattice-based kinetic Monte Carlo simulation. The simulation explicitly includes both the constituent atoms and is performed on the wurtzite crystal structure of GaN. The deposition flux is chosen to mimic conditions in a molecular beam epitaxy chamber. The surface diffusion barriers used for the growth simulation are obtained from ab initio density functional theory based calculations. The evolution of submonolayer islands from small random clusters to ordered triangular islands is captured in these simulations. Further the submonolayer island density statistics is calculated to illustrate standard scaling behavior. The island density as a function of coverage is used to show the different regimes of submonolayer growth and the transition to multi-layer growth. The island density dependence on deposition flux, Ga:N flux ratio, and temperature are also shown to be reasonable and consistent with experiments. We further highlight the importance of N-rotation, which is a diffusion mechanism for nitrogen adatoms under moderate gallium excess.

Lattice kinetic Monte Carlo simulation of epitaxial growth of silicon thin films in H2/SiH4 chemical vapor deposition systems

The use of lattice kinetic Monte Carlo techniques for reproducing epitaxial growth of silicon using H2/SiH4 chemistry is described in this work. The implementation of the surface mechanisms is limited to silane and hydrogen adsorption, surface dissociation of higher order hydrides and recombinative hydrogen desorption. Good agreement between simulation and experimental results available in the literature was obtained for H2 desorption peak positions, equilibrium hydrogen coverages and Si growth rates. The formation of (111) and (113) planes is also reproduced, however no evidence of (131) was found even at 850° C, and (119) planes were not formed at the same position than actual grown structures.

Parametric analysis of mechanically driven compositional patterning in SiGe substrates

A recently demonstrated approach for creating structured compositional gradients in the near-surface region of SiGe substrates is studied parametrically using a multiresolution coarse-grained lattice kinetic Monte Carlo simulation method. In the “stress patterning” process, a patterned elastic stress field is generated in the SiGe substrate by pressing an array of micro-indenters into it. The stressed substrate is then thermally annealed to drive the atomic diffusion in which the larger Ge atoms are pushed away from the areas of compressive stress. By varying a subset of the parameters that characterize the high-dimensional input space of the process (e.g., indenter spacing, indenter tip shape, and indenter array symmetry) we show that technologically interesting compositional configurations may be readily generated. In particular, we show that it is theoretically possible to generate arrays of well-delineated nanoscale regions of high Ge content surrounded by essentially pure Si. Such configurations may be useful as Ge “quantum dots” that exhibit three-dimensional quantum confinement, which have otherwise been very challenging to create with high degrees of size and spatial uniformity. These simulation results will be instrumental in guiding future experimental demonstrations of stress patterning.

A practical approach to the sensitivity analysis for kinetic Monte Carlo simulation of heterogeneous catalysis.

Lattice kinetic Monte Carlo simulations have become a vital tool for predictive quality atomistic understanding of complex surface chemical reaction kinetics over a wide range of reaction conditions. In order to expand their practical value in terms of giving guidelines for the atomic level design of catalytic systems, it is very desirable to readily evaluate a sensitivity analysis for a given model. The result of such a sensitivity analysis quantitatively expresses the dependency of the turnover frequency, being the main output variable, on the rate constants entering the model. In the past, the application of sensitivity analysis, such as degree of rate control, has been hampered by its exuberant computational effort required to accurately sample numerical derivatives of a property that is obtained from a stochastic simulation method. In this study, we present an efficient and robust three-stage approach that is capable of reliably evaluating the sensitivity measures for stiff microkinetic models as we demonstrate using the CO oxidation on RuO2(110) as a prototypical reaction. In the first step, we utilize the Fisher information matrix for filtering out elementary processes which only yield negligible sensitivity. Then we employ an estimator based on the linear response theory for calculating the sensitivity measure for non-critical conditions which covers the majority of cases. Finally, we adapt a method for sampling coupled finite differences for evaluating the sensitivity measure for lattice based models. This allows for an efficient evaluation even in critical regions near a second order phase transition that are hitherto difficult to control. The combined approach leads to significant computational savings over straightforward numerical derivatives and should aid in accelerating the nano-scale design of heterogeneous catalysts.

Analysis of the lattice kinetic Monte Carlo method in systems with external fields.

The lattice kinetic Monte Carlo (LKMC) method is studied in the context of Brownian particles subjected to drift forces, here principally represented by external fluid flow. LKMC rate expressions for particle hopping are derived that satisfy detailed balance at equilibrium while also providing correct dynamical trajectories in advective-diffusive situations. Error analyses are performed for systems in which collections of particles undergo Brownian motion while also being advected by plug and parabolic flows. We demonstrate how the flow intensity, and its associated drift force, as well as its gradient, each impact the accuracy of the method in relation to reference analytical solutions and Brownian dynamics simulations. Finally, we show how a non-uniform grid that everywhere retains full microscopic detail may be employed to increase the computational efficiency of lattice kinetic Monte Carlo simulations of particles subjected to drift forces arising from the presence of external fields.

The Lattice Kinetic Monte Carlo Simulation of Atomic Diffusion and Structural Transition for Gold

For the kinetic simulation of metal nanoparticles, we developed a self-consistent coordination-averaged energies for Au atoms based on energy properties of gold bulk phases. The energy barrier of the atom pairing change is proposed and holds for the microscopic reversibility principle. By applying the lattice kinetic Monte Carlo simulation on gold films, we found that the atomic diffusion of Au on the Au(111) surface undergoes a late transition state with an energy barrier of about 0.2 eV and a prefactor between 40~50 Å2/ps. This study also investigates the structural transition from spherical to faceted gold nanoparticles upon heating. The temperatures of structural transition are in agreement with the experimental melting temperatures of gold nanoparticles with diameters ranging from 2 nm to 8 nm.

(Invited) Multiscale Modeling of Stress-Mediated Compositional Patterning in SiGe Substrates

The ability to scalably fabricate periodic, large-area assemblies of Ge quantum dots (QDs) on Si or SiGe substrates with high degree of spatial and size uniformity is expected to impact numerous technologies, including nano-/micro-electronics, optoelectronics, nanosensor arrays, and high-density patterned media for data storage. One approach to achieve this goal is to harness the Stranski-Krastonov (SK) growth mechanism [1,2]. Here, a lattice mismatched species is heteroepitaxially deposited on a substrate (e.g., Ge on Si or InAs on GaAs), creating strain that increases with film thickness, leading to the spontaneous formation of a distribution of islands. The size and shape of such islands has attracted much attention in regards to their potential use as arrays of QDs, but island evolution during SK growth is inherently stochastic, leading to broad distributions of island size with random locations. As a result, much effort has been directed towards substrate patterning in order to spatially direct island formation and narrow the size distribution [2-4]. While these approaches have proven successful in principle, they are generally process-intensive and will be difficult to scale to large-area applications. Here, we describe a process in which SiGe substrates are compositionally patterned over large areas using spatially-modulated elastic fields applied by a nano-indenter array [5]. The indenter array, which is fabricated by interferometric lithography and dry etching, is pressed against a Si0.8Ge0.2 wafer in a custom-made mechanical press. The entire assembly is then annealed at high temperatures, during which the larger Ge atoms are selectively driven away from areas of compressive stress via the process of interdiffusion. Compositional analysis of the substrates demonstrates that this approach leads to a transfer of the indenter array pattern into the near-surface compositional distribution. In particular, the process leads to sharp compositional gradients that separate well-delineated regions in the near-surface region of the SiGe substrate, suggesting the possible use of this approach to simply and robustly create structures with quantum confinement properties. Employing this approach in practice to create useful configurations, however, requires a systematic assessment of the processing parameter space. For example, the final compositional distribution in the substrate depends on the indenter tip shapes and sizes, the indenter array geometric pattern, the annealing temperature-time history, as well as the initial substrate composition and thickness. In order to fulfill this task, we have developed a multiscale computer model based on the lattice kinetic Monte Carlo (LKMC) method to simulate the ‘stress transfer’ process. The LKMC simulation is propagated using rates for atomic diffusion that depend explicitly on local values of stress, composition, and temperature [6,7]. The dependence of atomic diffusion on composition is regressed to experimental data while the stress dependence is described using the theory of activation volumes [8]. The model is used to investigate the dynamics of the compositional patterning process, as well as the impact of various parameters on the compositional configurations that result. We find that certain indenter configurations produce compositional patterns that may be favorable for engineering lateral arrays quantum-confined structures. Some comparisons also are made to preliminary experimental demonstrations. [1] F. M. Ross, R. M. Tromp and M. C. Reuter, Science, 1999, 286, 1931. [2] O. G. Schmidt et al., Appl. Phys. Lett., 2000, 77, 4139. [3] F. M. Ross, M. Kammler, M. C. Reuter and R. Hull, Philos. Mag., 2004, 84, 2687. [4] J. J. Zhang et al., Appl. Phys. Lett., 2007, 91, 173115. [5] S. Ghosh, D. Kaiser, J. Bonilla, T. Sinno and S. M. Han, Appl. Phys. Lett., 2015, 107, 072106. [6] P. Castrillo, R. Pinacho, M. Jaraiz and J. E. Rubio, J. Appl. Phys., 2011, 109, 103502. [7] D. Kaiser, S. Ghosh, S. M. Han, and T Sinno, to appear in Mol. Sys. Des. & Eng. (2016). [8] M. J. Aziz, Applied Physics Letters 70, 2810 (1997).

Simulation and characterization of stress in FinFETs using novel LKMC and nanobeam diffraction methods**Project supported by the “National S&T Major Project 02”, the Opening Project of Microelectronics Devices & Bulk Si FinFET Integrated Technology, Institute of Microelectronics, Chinese Academy of Sciences (No. 2013ZX02303007-001).

A new simulation method and test instrument has been adopted to verify the traditional stress simulation in FinFET. First, a new algorithm named lattice kinetic Monte Carlo (LKMC) is used to simulate the SiGe epitaxy in source/drain regions and the stress distribution is consequently extracted after the LKMC simulation. Systematic comparison between the traditional polyhedron method and the LKMC method is carried out. The results confirm that extracted stress from both methods is consistent, which verifies the validity of traditional polyhedron method for the purpose of simulating stress in FinFET. In the following experiment, p-type FinFETs with SiGe stressors in source/drain regions are fabricated. The nano beam diffraction (NBD) method is employed to characterize the strain in Si fin. The strain value from the NBD test agrees well with the value extracted from traditional polyhedron simulation.

Kinetic Monte Carlo simulations of vacancy evolution in graphene

Kinetic lattice Monte Carlo simulations, implementing ab initio calibrated energetics, are applied to study the evolution of vacancy systems in graphene. The evolution of the under-coordinated atomic layer proceeds, in general, by island nucleation and growth. However, the first stage (nucleation) is strongly influenced by the stability of small aggregates (e.g. di-vacancies) and by the effective coalescence/attachment barriers. Quantitative predictions of the system kinetics in terms of crystal state and defects' morphology as a function of the initial state and the temperature can be obtained and readily compared with experimental structural characterization of processed samples.

Atomistic Modeling of Strained SiGe Nanofins

We utilize a hierarchical approach to atomistic modeling with density functional theory (DFT) calculations being used together with experimental observations to determine parameters for molecular dynamics (MD) potentials and kinetic lattice Monte Carlo (KLMC) simulations. MD and KLMC are used to model epitaxial growth and impurity redistribution in strained group IV alloy nanostructures.

Atomistic Modeling of Strained SiGe Nanofins

We utilize a hierarchical approach to atomistic modeling with density functional theory (DFT) calculations being used together with MD simulations and experimental observations to determine parameters for KLMC simulations. These methods are used to model epitaxial growth of strained group IV alloys and impurity redistribution in strained alloy nanostructures. We have developed kinetic lattice Monte Carlo (KLMC) methods in which the basic processes correspond directly to transitions at the atomic level. These methods are able to predict phenomena which rely on the local atomic structure, including diffusion at extremely high donor concentrations and impact of diffusion on alloy composition and anisotropic strain.Use of group IV alloy species, particularly Ge, via selective epi and/or ion implantationMoving from a pure material, Si, to a random alloy such as SiGe (or SiGeCSn) means that every lattice site has a different local atomic arrangement and thus different atomic transition rate. We utilize kinetic lattice Monte Carlo (KLMC) methods within a hierarchical framework (DFT/MD/KLMC/FEM) in order to understand epitaxial film growth, impurity redistribution and extended defect formation processes on a fundamental level in nanoscale device structures based on group IV alloys. Within KLMC, alloy species (Si, Ge, C and Sn), substitutional dopants (e.g., B, P, As), point defects (I, V), dopant/defect pairs which control diffusion (e.g., BI, AsV), dopant complexes which control activation (e.g., B3I, As4V), and point defect aggregates (e.g., 311s, dislocation loops), as well as the transitions between these meta-stable structures can all be represented atomistically with a topology that matches the physical system. Thus, models can be made more realistic and therefore predictive than is possible with continuum or non-lattice-based kMC methods.KLMC methods provide the crucial link needed for applying ab-initio results to nondilute systems (e.g., alloys, high doping levels). In KLMC, the system is considered to evolve between a series of metastable states, enabling direct atomistic simulation of material behavior for the long time scales (e.g., seconds) relevant for nanodevice fabrication. These methods are able to predict phenomena which rely on the local atomic structure, including orientation-dependent epitaxial growth, diffusion at extremely high donor concentrations and impact of diffusion on alloy composition and anisotropic strain. KLMC can also predict the correlated statistics of dopants, which have traditionally been ignored, but have a large impact on device variability.Density functional theory (DFT) calculations show that the formation and migration energies of boron interstitial complexes in SiGe depend strongly on both the local stress tensor as well as the local arrangement of Ge atoms. Using transition/hopping rates calculated from DFT in KLMC simulations gives remarkably accurate prediction of B diffusion, segregation and solubility in strained SiGe structures. Ge enhances V-mediated diffusion due to Ge-V binding, but the lower concentration of Ge leads to tensile stress in the interdiffusion region that retards vacancy-mediated diffusion. We calculate binding energies and transition barriers for Si/Ge/Sn/B/As systems to determine KLMC parameters for simulations of impurity/alloy redistribution.For modeling of epi gtrowth, DFT calculations are used to estimate the binding energies as function of local atomic arrangement and transition barrier energies for key processes such as adatom and edge diffusion, desorption, and layer incorporation. These calculations are then used to determine rates for KLMC simulation of the epitaxial growth process. Within this framework, complex phenomena such as faceting and orientation-dependent growth rate are a natural result of differences in binding energies and transition rates. We utilize molecular dynamics simulations as an additional supporting tool to connect DFT to KLMC and continuum simulations. Of particular interest for this work are MD calculations of local stress/strain based on distributions of alloy and dopant atoms. Note that lattice deformations are induced not only due to impurity incorporation but also by crystal surfaces and edges, leading to stresses even in pure Si systems. The effects of surfaces are most pronounced in nanoscale systems with large surface to volume ratios. Induced strains depend on surface orientation, so local anisotropic stresses can favor particular facets. The surface effect causes expansion in the surface plane and induces significant tensile strain even far away from the surface. A key challenge is incorporation of the effect of stresses due to impurity distributions on the redistribution of these and other impurities. We couple the KLMC simulations to MD simulations of local strain. These calculations are fully general and include impurity, point defect, extended defect and surface/interface strain effects.

KMCLib: A general framework for lattice kinetic Monte Carlo (KMC) simulations

KMCLib is a general framework for lattice kinetic Monte Carlo (KMC) simulations. The program can handle simulations of the diffusion and reaction of millions of particles in one, two, or three dimensions, and is designed to be easily extended and customized by the user to allow for the development of complex custom KMC models for specific systems without having to modify the core functionality of the program. Analysis modules and on-the-fly elementary step diffusion rate calculations can be implemented as plugins following a well-defined API. The plugin modules are loosely coupled to the core KMCLib program via the Python scripting language. KMCLib is written as a Python module with a backend C++ library. After initial compilation of the backend library KMCLib is used as a Python module; input to the program is given as a Python script executed using a standard Python interpreter. We give a detailed description of the features and implementation of the code and demonstrate its scaling behavior and parallel performance with a simple one-dimensional A–B–C lattice KMC model and a more complex three-dimensional lattice KMC model of oxygen-vacancy diffusion in a fluorite structured metal oxide. KMCLib can keep track of individual particle movements and includes tools for mean square displacement analysis, and is therefore particularly well suited for studying diffusion processes at surfaces and in solids. Program summaryProgram title: KMCLibCatalogue identifier: AESZ_v1_0Program summary URL:http://cpc.cs.qub.ac.uk/summaries/AESZ_v1_0.htmlProgram obtainable from: CPC Program Library, Queen’s University, Belfast, N. IrelandLicensing provisions: GNU General Public License, version 3No. of lines in distributed program, including test data, etc.: 49064No. of bytes in distributed program, including test data, etc.: 1575172Distribution format: tar.gzProgramming language: Python and C++.Computer: Any computer that can run a C++ compiler and a Python interpreter.Operating system: Tested on Ubuntu 12.4 LTS, CentOS release 5.9, Mac OSX 10.5.8 and Mac OSX 10.8.2, but should run on any system that can have a C++ compiler, MPI and a Python interpreter.Has the code been vectorized or parallelized?: Yes. From one to hundreds of processors depending on the type of input and simulation.RAM: From a few megabytes to several gigabytes depending on input parameters and the size of the system to simulate.Classification: 4.13, 16.13.External routines: KMCLib uses an external Mersenne Twister pseudo random number generator that is included in the code. A Python 2.7 interpreter and a standard C++ runtime library are needed to run the serial version of the code. For running the parallel version an MPI implementation is needed, such as e.g. MPICH from http://www.mpich.org or Open-MPI from http://www.open-mpi.org. SWIG (obtainable from http://www.swig.org/) and CMake (obtainable from http://www.cmake.org/) are needed for building the backend module, Sphinx (obtainable from http://sphinx-doc.org) for building the documentation and CPPUNIT (obtainable from http://sourceforge.net/projects/cppunit/) for building the C++ unit tests.Nature of problem:Atomic scale simulation of slowly evolving dynamics is a great challenge in many areas of computational materials science and catalysis. When the rare-events dynamics of interest is orders of magnitude slower than the typical atomic vibrational frequencies a straight-forward propagation of the equations of motions for the particles in the simulation cannot reach time scales of relevance for modeling the slow dynamics.Solution method:KMCLib provides an implementation of the kinetic Monte Carlo (KMC) method that solves the slow dynamics problem by utilizing the separation of time scales between fast vibrational motion and the slowly evolving rare-events dynamics. Only the latter is treated explicitly and the system is simulated as jumping between fully equilibrated local energy minima on the slow-dynamics potential energy surface.Restrictions:KMCLib implements the lattice KMC method and is as such restricted to geometries that can be expressed on a grid in space.Unusual features:KMCLib has been designed to be easily customized, to allow for user-defined functionality and integration with other codes. The user can define her own on-the-fly rate calculator via a Python API, so that site-specific elementary process rates, or rates depending on long-range interactions or complex geometrical features can easily be included. KMCLib also allows for on-the-fly analysis with user-defined analysis modules. KMCLib can keep track of individual particle movements and includes tools for mean square displacement analysis, and is therefore particularly well suited for studying diffusion processes at surfaces and in solids.Additional comments:The full documentation of the program is distributed with the code and can also be found at http://www.github.com/leetmaa/KMCLib/manualRunning time:rom a few seconds to several days depending on the type of simulation and input parameters.