Ling Ti Kong

Phonon dispersion measured directly from molecular dynamics simulations

Abstract A method to measure the phonon dispersion of a crystal based on molecular dynamics simulation is proposed and implemented as an extension to an open source classical molecular dynamics simulation code LAMMPS. In the proposed method, the dynamical matrix is constructed by observing the displacements of atoms during molecular dynamics simulation, making use of the fluctuation–dissipation theory. The dynamical matrix can then be employed to compute the phonon spectra by evaluating its eigenvalues. It is found that the proposed method is capable of yielding the phonon dispersion accurately, while taking into account the anharmonic effect on phonons simultaneously. The implementation is done in the style of fix of LAMMPS, which is designed to run in parallel and to exploit the functions provided by LAMMPS; the measured dynamical matrices could be passed to an auxiliary postprocessing code to evaluate the phonons. Program summary Program title: FixPhonon, version 1.0 Catalogue identifier: AEJB_v1_0 Program summary URL: http://cpc.cs.qub.ac.uk/summaries/AEJB_v1_0.html Program obtainable from: CPC Program Library, Queenʼs University, Belfast, N. Ireland Licensing provisions: GNU General Public license No. of lines in distributed program, including test data, etc.: 105 393 No. of bytes in distributed program, including test data, etc.: 3 231 800 Distribution format: tar.gz Programming language: C++ Computer: All Operating system: Linux Has the code been vectorized or parallelized?: Yes. 1 to N processors may be used RAM: Depends on problem, ≈1 kB to several MB Classification: 7.8 External routines: MPI, FFT, LAMMPS version 15, January 2010 ( http://lammps.sandia.gov/ ) Nature of problem: Atoms in solids make ceaseless vibrations about their equilibrium positions, and a collective vibration forms a wave of allowed wavelength and amplitude. The quantum of such lattice vibration is called the phonon, and the so-called “lattice dynamics” is the field of study to find the normal modes of these vibrations. In other words, lattice dynamics examines the relationship between the frequencies of phonons and the wave vectors, i.e., the phonon dispersion. The evaluation of the phonon dispersion requires the construction of the dynamical matrix. In atomic scale modeling, the dynamical matrices are usually constructed by deriving the derivatives of the force field employed, which cannot account for the effect of temperature on phonons, with an exception of the tedious “quasi-harmonic” procedure. Solution method: We propose here a method to construct the dynamical matrix directly from molecular dynamics simulations, simply by observing the displacements of atoms in the system thus making the constructing of the dynamical matrix a straightforward task. Moreover, the anharmonic effect was taken into account in molecular dynamics simulations naturally, the resultant phonons therefore reflect the finite temperature effect simultaneously. Restrictions: A well defined lattice is necessary to employ the proposed method as well as the implemented code to evaluate the phonon dispersion. In other words, the system under study should be in solid state where atoms vibrate about their equilibrium positions. Besides, no drifting of the lattice is expected. The method is best suited for periodic systems, although non-periodic system with a supercell approach is also possible, it will however become inefficient when the unit cell contains too many atoms. Additional comments: The readers are encouraged to visit http://code.google.com/p/fix-phonon for subsequent update of the code as well as the associated postprocessing code, so as to keep up with the latest version of LAMMPS. Running time: Running time depends on the system size, the numbers of processors used, and the complexity of the force field, like a typical molecular dynamics simulation. For the third example shown in this paper, it took about 2.5 hours on an Intel Xeon X3220 architecture (2.4G, quadcore). References: [1] C. Campana, M.H. Muser, Phys. Rev. B 74 (2006) 075420. [2] L.T. Kong, G. Bartels, C. Campana, C. Denniston, M.H. Muser, Comp. Phys. Commun. 180 (6) (2009) 1004–1010.

Implementation of Green's function molecular dynamics: an extension to LAMMPS

The Green’s function molecular dynamics method, which enables one to study the elastic response of a three-dimensional solid to an external stress field by taking into consideration only the surface atoms, was implemented as an extension to an open source classical molecular dynamics simulation code LAMMPS. This was done in the style of fixes. The first fix, FixGFC, measures the elastic stiffness coefficients for a (small) solid block of a given material by making use of the fluctuation-dissipation theorem. With the help of the second fix, FixGFMD, the coefficients obtained from FixGFC can then be used to compute the elastic forces for a (large) block of the same material. Both fixes are designed to be run in parallel and to exploit the functions provided by LAMMPS.

The crucial role of chemical detail for slip-boundary conditions: molecular dynamics simulations of linear oligomers between sliding aluminum surfaces

We study the slip-boundary conditions of short, linear paraffins and olefins confined between two sliding aluminum surfaces with molecular dynamics. Our simulations are based on a recently developed force field for the interaction between organic molecules and bulk aluminum. The lubricant molecules investigated all consist of six monomers but differ in the existence or location of merely one double bond. It turns out that this small change in the chemistry of the lubricant molecules can alter slip lengths quite dramatically, and is not strongly correlated with surface energies and bulk viscosity of the lubricant. For example, α and β-hexene have similar large negative slip length of Λ ≈ −8 A, even though α-hexene adheres twice as strongly to the surface as β-hexene. Eliminating the double bond in β-hexene reduces the surface energy by another factor of two, but increases Λ from −8 to 120 A. These results and those of additional simulations based on unrealistic, albeit occasionally used model potentials, make us conclude that surface energies and/or molecular geometries alone are not reliable indicators for slip-boundary conditions. Instead, it is necessary to consider the full chemical detail. As a more encouraging result, we find that the bulk viscosity appears to describe the dissipation within the sheared fluid close to the wall quite well, despite significant ordering near the boundaries. Moreover, all our systems show a relatively weak dependence of the slip length on the normal pressure.

Dynamical stability of iron under high-temperature and high-pressure conditions

The dynamical stability of iron under high-temperature and high-pressure conditions was investigated based on the phonons evaluated by using a recently developed method. It is revealed that both the fcc-Fe and the hcp-Fe are dynamically stable in a wide temperature and pressure range. The bcc-Fe phase can be stable as well, while in a limited temperature/pressure regime bounded by a dynamical stability limit and a harmonic limit. Direct evidence shows that it is the entropy term that plays a critical role in stabilizing the bcc-Fe under high-temperature and high-pressure conditions.

Non-bonded force field for the interaction between metals and organic molecules: a case study of olefins on aluminum

In this work, we parameterize an empirical potential for the interaction between organic molecules and metal surfaces via force matching. This is done by pursuing a self-consistent approach similar to the ones used for equilibrium simulations; however, special attention is paid to the suitability of the resulting potential for tribological (non-equilibrium) situations. Specifically, we study olefin molecules confined between two aluminum surfaces under realistic pressures and shear rates. We find that the Buckingham potential produces better agreement with the first principle data than other force fields. While our training set only contains hexene molecules, we find that the standard error in the fitted olefin-aluminum interaction increases only by a factor of 1.15 when the force field is applied to butene, octene, and decene. Including mirror charges into the treatment only marginally improves fits. While olefins on aluminum is merely a special case, the proposed methodology can be used to parameterize any other interaction between polymers and metal surfaces for use in tribological simulations.

Quantitative results for square gradient models of fluids

Square gradient models for fluids are extensively used because they are believed to provide a good qualitative understanding of the essential physics. However, unlike elasticity theory for solids, there are few quantitative results for specific (as opposed to generic) fluids. Indeed the only numerical value of the square gradient coefficients for specific fluids have been inferred from attempts to match macroscopic properties such as surface tensions rather than from direct measurement. We employ all-atom molecular dynamics, using the TIP3P and OPLS force fields, to directly measure the coefficients of the density gradient expansion for several real fluids. For all liquids measured, including water, we find that the square gradient coefficient is negative, suggesting the need for some regularization of a model including only the square gradient, but only at wavelengths comparable to the molecular separation of molecules. The implications for liquid-gas interfaces are also examined. Remarkably, the square gradient model is found to give a reasonably accurate description of density fluctuations in the liquid state down to wavelengths close to atomic size.

An improved version of the Green's function molecular dynamics method ✩,✩✩

This work presents an improved version of the Greens function molecular dynamics method (Kong et al., 2009; Campana and Muser, 2004 [1,2]), which enables one to study the elastic response of a three-dimensional solid to an external stress field by taking into consideration only atoms near the surface. In the previous implementation, the effective elastic coefficients measured at the Γ-point were altered to reduce finite size effects: their eigenvalues corresponding to the acoustic modes were set to zero. This scheme was found to work well for simple Bravais lattices as long as only atoms within the last layer were treated as Greens function atoms. However, it failed to function as expected in all other cases. It turns out that a violation of the acoustic sum rule for the effective elastic coefficients at Γ (Kong, 2010 [3]) was responsible for this behavior. In the new version, the acoustic sum rule is enforced by adopting an iterative procedure, which is found to be physically more meaningful than the previous one. In addition, the new algorithm allows one to treat lattices with bases and the Greens function slab is no longer confined to one layer.