Tersoff Potential Research Articles

Structure-properties relations in graphene derivatives and metamaterials obtained by atomic-scale modeling

ABSTRACTRecent findings of atomic-scale modelling studies are reviewed on graphene derivatives and metamaterials fabricated through chemical functionalization and/or defect engineering of graphene sheets. Results of molecular-statics and molecular-dynamics simulations according to a reliable bond-order potential, as well as first-principles density functional theory calculations are reviewed that have established useful structure-properties relations in two-dimensional materials, such as graphene nanomeshes (GNMs), electron-irradiated graphene, and interlayer-bonded twisted bilayer graphene. Quantitative relationships are established for the elastic moduli, mechanical properties, and thermal conductivity of GNMs as a function of the nanomesh porosity and the mechanical response of GNMs to uniaxial tensile straining is explored over the range of nanomesh porosities. The dependence of structural, mechanical, and thermal transport properties of electron-irradiated graphene sheets on the density of irradiation-induced defects is reviewed, highlighting an amorphization transition accompanied by a brittle-to-ductile transition and a transition in thermal transport mechanism beyond a critical defect concentration. The tunability of the electronic band structure, mechanical properties, and structural response to mechanical loading of graphene-diamond nanocomposite superstructures consisting of nanodiamond superlattices in interlayer-bonded twisted bilayer graphene also is demonstrated by precise control of the density and distribution of covalent interlayer C–C bonds.

Molecular dynamics study of the material property changes induced by accumulated point defects in graphite

Molecular dynamics simulations were performed on pristine graphite using an adaptive intermolecular reactive empirical bond-order (AIREBO) potential, which introduced various concentrations of point defects and different sizes in defect clusters. Induced by the accumulation of mono-vacancy and self-interstitial atoms (SIA), the changes to the geometry, density, Young’s modulus, fracture strength, and coefficient of thermal expansion (CTE) of graphite were simulated over a range of defect concentrations and cluster sizes. We found that accumulated Frenkel pairs lead to clearly-observed changes to the geometric and elastic properties in graphite crystal, but a relatively negligible effect on CTE. As expected, changing size of defect clusters with a given concentration relates to the variation of material properties. The property changes simulated in the present work are qualitatively consistent with the experimental observations. The results are expected to provide some insights into the link between cascade-induced point defects and experimentally irradiation-assisted material property changes.

A molecular dynamics study of the mechanical properties of h-BCN monolayer using a modified Tersoff interatomic potential

Using molecular dynamics (MD) simulations, we investigate the mechanical properties of hexagonal BCN monolayer, a newly synthesized two-dimensional material with an atom ratio of B/C/N = 1:1:1. The Tersoff potential is modified to get good agreement between predicted and measured fracture strengths of graphene. With this modified Tersoff potential, we perform extensive MD simulations to study the effect of temperature, strain rate and vacancy defect on the mechanical properties of h-BCN. It is found that h-BCN is a strong material with fracture strength of 81.4–93.5 GPa, albeit ∼35% lower than that of graphene. Similar to graphene, temperature has strong effect on the mechanical properties of h-BCN. As the temperature increases from 10 K to 1300 K, the fracture strength and strain of h-BCN drops by 55% and 62%, respectively. The strain rate is found to have a moderate effect. When the strain rate increases from 0.00002 to 0.0125 ps−1, the fracture strength and strain of h-BCN increases 6.1% and 12%, respectively. As for the atomic defect, a very small concentration (0.028%) of vacancy in h-BCN is able to cause a 28% reduction in fracture strength and a 35.5% reduction in fracture strain. These findings have significance for its future applications in nanodevices.

Mechanical and fracture behavior of water submerged graphene

In this article, molecular dynamics (MD)-based simulations were performed to study the structural stability, tensile strength, and fracture behavior of graphene submerged in water. MD-based atomistic simulations were performed in conjunction with adaptive intermolecular reactive empirical bond order and TIP3P potentials for capturing the interatomic interaction in graphene and water molecules, respectively, whereas nonbonded interactions between the two were estimated with the help of Lennard-Jones potential. It was predicted from the simulations that water has a significant effect on the structural stability and the fracture behavior of graphene. The structural stability of water submerged graphene decreases as compared to dry graphene, whereas the fracture behavior of water submerged graphene shows dependency on the chirality and the passivation of crack edge atoms with hydrogen. It was also predicted from the simulations that hydrogen passivation helps in reducing the reactivity of crack edge atoms with the water molecules. Results presented in this article will help in exploring the full potential of graphene for water submerged applications.

Physically\xa0informed artificial neural networks for atomistic modeling of materials

Large-scale atomistic computer simulations of materials heavily rely on interatomic potentials predicting the energy and Newtonian forces on atoms. Traditional interatomic potentials are based on physical intuition but contain few adjustable parameters and are usually not accurate. The emerging machine-learning (ML) potentials achieve highly accurate interpolation within a large DFT database but, being purely mathematical constructions, suffer from poor transferability to unknown structures. We propose a new approach that can drastically improve the transferability of ML potentials by informing them of the physical nature of interatomic bonding. This is achieved by combining a rather general physics-based model (analytical bond-order potential) with a neural-network regression. This approach, called the physically informed neural network (PINN) potential, is demonstrated by developing a general-purpose PINN potential for Al. We suggest that the development of physics-based ML potentials is the most effective way forward in the field of atomistic simulations.

Theoretical investigations on the structural stability and miscibility in BAlN and BGaN alloys: bond-order interatomic potential calculations

Structural stability and miscibility of BxAl1−xN and BxGa1−xN alloys including 3-fold coordinated hexagonal (Hex), 4-fold coordinated wurtzite (WZ) and zinc blende (ZB) structures are systematically investigated based on empirical bond-order potential (BOP), which incorporates electrostatic energies due to bond and ionic charges. We find that the Hex structure is stabilized for free-standing BxAl1−xN and BxGa1−xN alloys with boron composition while the ZB (WZ) structure is stabilized for () both in BxAl1−xN and BxGa1−xN alloys due to the electrostatic energy. Furthermore, the miscibility of BxAl1−xN and BxGa1−xN alloys with 4-fold coordinated structures are found to be higher than that with Hex structure owing to the chemical energy contribution to the excess energy. These results suggest that our empirical interatomic potentials approach is feasible to clarify structural stability not only for the 3-fold coordinated Hex structure but also 4-fold coordinated WZ and ZB structures in BxAl1−xN and BxGa1−xN alloys.

Exploring the Effect of Dihedral Energy on the Nonlinear Mechanics of the Carbon Nanotubes Using a Multiscale Modeling

Abstract A hierarchical multiscale finite element model is employed to investigate the effect of dihedral energy term on the numerical simulation of two-dimensional materials. The numerical examples of the carbon nanotubes and graphene sheets are studied employing a refined constitutive model in conjunction with a multiscale finite element method. The constitutive law refined with the greater accuracy on the bending modulus using second generation reactive empirical bond order potential with dihedral energy term is employed to investigate the linear and nonlinear response of the carbon nanotubes incorporating material and Green–Lagrange geometric nonlinearities. The inclusion of the dihedral energy term predicts bending modulus close to those of through first principle calculations. The deformations at the nanoscale and macroscopic scales are related through the Cauchy–Born rule. The effect of the dihedral energy term on the response of the carbon nanotubes is studied in detail. The governing equation of motion for the carbon nanotubes is formulated through Hamilton’s energy principle. The spatial approximation of the carbon nanotubes at the continuum scale is attained through the finite element method. The membrane locking in the circumferential strain is eliminated through the membrane consistent interpolation functions obtained through the least-square method.

Theoretical Study on the Thermoelectric Properties of Porous Armchair Graphene Nanoribbons

Thermoelectric properties of porous graphene nanoribbons (GNRs) have been explored for a range of pore dimensions in order to achieve a high performance one-dimensional nanoscale thermoelectric device. In this paper, study has been done to observe the effect of different nanoporous shapes and their associated positions on the thermoelectric properties of GNRs. The aim of this work is to study the effect of various circular, triangular, rectangular and rhombus shape dimensions so as to tune the pore to its optimal dimension that would enhance the overall thermoelectric efficiency. Further, the effect of passivation of pore edges has been studied for all shapes so as to observe its effect on thermoelectric performance. Also, the effect of temperature dependence on thermoelectric efficiency has been studied. Ballistic transport regime and semi empirical method using Huckel basis set is used to obtain the electrical properties while the tersoff potential is used for the phononic system.

A Theoretical Study on the Thermoelectric Properties of Porous Armchair Germanene Nanoribbons

Since the limits of conventional sources of energy are rapidly approaching, the thermoelectric devices have attracted attention for their potential of power generation directly from waste heat. In this paper, thermoelectric properties of porous armchair Germanene nanoribbons (AGeNRs) have been explored for a range of pore dimensions in order to achieve a high performance two-dimensional nanoscale thermoelectric device. The work has been done to investigate the influence of different nanopore shapes and their associated positions on the thermoelectric performance so as to tune it to the optimum pore shape and position that would enhance the overall thermoelectric efficiency. Also, the effect of passivation of the pore edges on thermoelectric parameters for all shapes has been studied. Further, the influence of temperature dependence on figure of merit has been observed. Ballistic transport regime and semi-empirical method using Huckel basis set is used to obtain the electrical properties while the Brenner potential is used for the phononic properties.

(Invited) Ion Conduction in Glassy Electrolytes for Sodium Ion Batteries – an Atomistic Perspective

Solid-state sodium ion battery is a promising technology for energy storage in electrical grids. However, the characteristic low ionic conductivity of solid-state electrolytes is a major challenge towards commercialization of solid-state sodium ion batteries. Sulfide glasses have relatively high ionic conductivity at room temperature and are therefore considered promising as potential electrolytes. However, the energy barrier for ion hops in these electrolytes is highly sensitive to the chemical composition, which provides a tremendous opportunity to tailor the ionic conductivity. In-situ characterization of Na+ ion transport through experiments is extremely challenging. Atomistic simulations offer an economical route for determining the local structure of these glasses as well as Na+ ion transport mechanisms. Sodium sulfide – silicon sulfide [xNa2S – (1-x)SiS2] based electrolytes with varying compositions was chosen as a model system and atomistic simulations of various flavors were employed to determine ion conduction mechanism and evaluate ionic conductivity, each method having distinct advantages. Experimental X-ray scattering data was used for the purpose of model validation. Ab initio molecular dynamics (MD) simulations are able to accurately model the specific chemistry based on on-the-fly quantum mechanical calculations to determine interatomic interactions. Augmented with nudged elastic band analysis, these calculations can determine activation energies for ion hops. However, the maximum system size that can be realistically modeled using ab initio MD is limited to few hundred atoms and therefore cannot reliably calculate the bulk properties. On the other hand, classical MD simulations employ empirical forcefields to define pairwise interactions between constituent ions in these glasses. This method can simulate much larger systems comprising hundreds of thousands of atoms. However, the inability of these predefined forcefields to accurately capture the complex chemistry in sulfide glasses makes it extremely challenging to accurately model these glasses. Reactive force fields, such as ReaxFF, is based on a bond order potential in conjunction with charge equilibration that is based on extensive parameterization based on quantum mechanical calculations of structure and energies. ReaxFF forcefields are highly system specific. Classical MD simulations involving ReaxFF provide an effective route to accurately capture bond breakage and formation, while simulating considerably large systems (~ 100,000 atoms) and therefore accurately evaluates ionic conductivity. Based on simulation data from each of the above methods, their merits and demerits are discussed as pertains to the model system of sulfide glasses. We determined the ionic conductivities and activation energies for sodium ion hops and related them to the structural units and hop distances. Overall, these simulations provide fundamental insights into ion conduction mechanisms as well as correlate ionic conductivity with glass structure and composition. The combined approach demonstrates excellent promise in identifying optimal glassy electrolytes for sodium ion batteries.

Computational study on surface energy of amorphous DLC with respect to hybridization state of carbon and potential functions

In this study, molecular dynamics (MD) simulations were performed to measure the surface energy of amorphous carbon (a-C) in relation to its hybridization state. Liquid quenching method was used to prepare the a-C samples. To control the sp2/sp3 ratio of a-C, two different approaches were attempted, i.e., changing the quenching rate and the density of carbon atoms. The results showed that the quenching rate did not yield any appreciable change in sp2/sp3 ratio, whereas the number of sp3 bonding was proportional to the density. With the controlled sp2/sp3 ratio in a-C, its surface energy was measured using two interactive potentials, i.e., Tersoff and Airebo, where the coordination cutoff length was set to 1.75 Å. In both potentials, the measured surface energy showed a non-monotonous correlation with the hybridization state. Even though the two potentials produced very similar value of maximum surface energy of a-C, they had different trends of surface energy with sp3% (=the percentage of sp3 hybridized atoms with respect to total number of atoms in the sample) – the measured surface energy increased with sp3% under Tersoff potential but decreased with sp3% under Airebo potential. This trend of surface energy with carbon hybridization state can be interpreted by the combined effects of non-polar interaction (or dispersive component) induced by sp3 bonds and polar interactions by dangling bonds. In case of Tersoff potential, the dominant factor for the surface energy value was either the dispersive or the polar components depending on the level of sp3% in the sample, whereas the surface energy under Airebo potential was dominantly affected by the number of triple dangling bonds (sp0).

Nanomechanical properties of single- and double-layer graphene spirals: a molecular dynamics simulation

This paper relates to a computational investigation of nanomechanical properties of graphene spirals. The molecular dynamics simulation method was used to investigate the mechanical properties including the stress–strain and force–strain diagrams under tensile tests as well as the fracture characteristics of the single- and double-layer graphene spirals. The adaptive intermolecular reactive empirical bond order potential was employed to model the covalent bonds and van der Waals interactions between the carbon atoms. Also, in the last section of the paper, the mechanical behavior of the spirals is scrutinized with respect to nitrogen and boron doping with various percentages and the Young’s moduli of the graphene spirals are presented as the functions of size and doping ratios according to the stress–strain diagrams. The results reveal three major deformation phases namely, elastic due to the van der Waals interactions, elastic due to the covalent bonds, and inelastic regimes. According to the results, the graphene spirals have superelastic characteristics in the range of 2000–3000% strains and very high strength values depending on the nanostructure size.

Thermal transport across grain boundaries in polycrystalline silicene: A multiscale modeling

During the fabrication process of large scale silicene, through common chemical vapor deposition (CVD) technique, polycrystalline films are quite likely to be produced, and the existence of Kapitza thermal resistance along grain boundaries could result in substantial changes of their thermal properties. In the present study, the thermal transport along polycrystalline silicene was evaluated by performing a multiscale method. Non-equilibrium molecular dynamics simulations (NEMD) was carried out to assess the interfacial thermal resistance of various constructed grain boundaries in silicene. The effects of tensile strain and the mean temperature on the interfacial thermal resistance were also examined. In the following stage, the effective thermal conductivity of polycrystalline silicene was investigated considering the effects of grain size and tensile strain. Our results indicate that the average values of Kapitza conductance at grain boundaries at room temperature were estimated to be nearly 2.56 × 109 W/m2 K and 2.46 × 109 W/m2 K through utilizing Tersoff and Stillinger-Weber interatomic potentials respectively. Also, in spite of the mean temperature, whose increment does not change Kapitza resistance, the interfacial thermal resistance could be controlled by applying strain. Furthermore, it was found that by tuning the grain size of polycrystalline silicene, its thermal conductivity could be modulated up to one order of magnitude.

A numerical study on subsurface quality and material removal during ultrasonic vibration assisted cutting of monocrystalline silicon by molecular dynamics simulation

Molecular dynamics (MD) simulation is used to study the subsurface quality and material removal of single crystal silicon with a diamond tool during ultrasonic elliptical vibration assisted cutting (UEVAC), 1D ultrasonic vibration assisted cutting (1D UVAC) and traditional cutting (TC) process. In the simulations, a long-range analytical bond order potential is used to describe the interaction inside the silicon specimen, providing a more accurate depiction of the atomic scale mechanisms of ductile plasticity, brittle fracture, and structural changes in silicon. The results show that UEVAC and 1D UVAC in cutting brittle material silicon causes a much smaller cutting force, much lower von Mises stress at the subsurface, larger material remove rate, lower compressive normal stress and and smaller shear stress In addition, the hydrostatic stress of subsurface for TC and 1D UVAC is much higher than that for UEVAC, which results in fewer Si-II and bct5-Si formed from the original Si-I in UEVAC. Moreover, the number of other atoms for UEVAC is overall smaller than those of using TC and 1D UVAC, which confirms that UEVAC produces a better subsurface. And atomic flow field analysis shows that the UEVAC tends to cut silicon in a more ductile mode. Besides, the temperature in front of tool edge and below the tool flank face of TC is much higher. This means that 1D UVAC and UEVAC have a positive effect on the tool life. However, the temperature in subsurface zone is overall larger, which reveals that 1D UVAC and UEVAC have a negative effect on the subsurface temperature.

Analytical interatomic bond-order potential for simulations of oxygen defects in iron

We present an analytical bond-order potential for the Fe–O system, capable of reproducing the basic properties of wüstite as well as the energetics of oxygen impurities in -iron. The potential predicts binding energies of various small oxygen-vacancy clusters in -iron in good agreement with density functional theory results, and is therefore suitable for simulations of oxygen-based defects in iron. We apply the potential in simulations of the stability and structure of Fe/FeO interfaces and FeO precipitates in iron, and observe that the shape of FeO precipitates can change due to formation of well-defined Fe/FeO interfaces. The interface with crystalline Fe also ensures that the precipitates never become fully amorphous, no matter how small they are.

A molecular dynamics study of growth of silicon from melt under stress

Molecular dynamics (MD) simulations of <111> growth of silicon from melt under compressive stress have been carried out. The effect of stress on the growth mechanism and the defects formation has been studied through the simulations. Atomic interactions were described using Tersoff potential. The results indicate that the crystal growth is restrained under uniaxial compressive stress in the initial stage, but the crystal growth velocity returns to the level in absence of stress while the stress is released due to the appearance of dislocations. The competition between the FCC and the HCP stacking of biatom-layer (111) plane is discussed when silicon grows along the [111] direction from melt. In applying stress, the chance of nucleation and growth of HCP configuration increases, and consequently stacking fault area forms. The Shockley partial dislocations appear on the boundary between the normal area and the stacking fault area. Correspondingly, the stress is released, and the crystal growth velocity is restored.

Fracture Toughnesses and Crack Growth Angles of Single-Layer Graphyne Sheets

Recently, Shang et al. (Angew Chem Int Ed 57(3):774–778, 2018) have developed a method to synthesize ultrathin (around 1.9 nm) graphyne nanosheets. We reported here the mixed-mode I–II fracture toughnesses and crack growth angles of single-layer graphyne sheets using molecular dynamics (MD) simulations and the finite element (FE) method based on the boundary layer model, respectively. The various carbon–carbon bonds of graphyne sheets in the FE method are equated with the nonlinear Timoshenko beams based on the Tersoff–Brenner potential, where all the parameters of the nonlinear beams are completely determined based on the continuum modeling. All the results from the present FE method are reasonable in comparison with those from our MD simulations using the REBO potential. The present results show that both the critical stress intensity factors (SIFs) and the crack growth angle strongly depend on the chirality and loading angle $$\varphi $$ ( $$\varphi =90^{\circ }$$ and $$\varphi =0^{\circ }$$ representing pure mode I and pure mode II, respectively). Meanwhile, the fracture properties of single-layer cyclicgraphene and supergraphene sheets are also studied in order to compare with those of the graphyne sheets. The critical equivalent SIFs are derived as $$1.55<K_{{\text {eq-cy}}}$$ (cyclic) $$<1.95$$ nN A $$^{-3/2}$$ , $$1.64<K_{{\text {eq-gy}}}$$ (graphyne) $$<2.64$$ nN A $$^{-3/2}$$ and $$0.61<K_{{\text {eq-su}}}$$ (super) $$<2.04$$ nN A $$^{-3/2}$$ in the corresponding zigzag and armchair sheets using the MD simulations, while the SIFs are $$0.32<K_{{\text {eq-cy}}}$$ (cyclic) $$<0.48$$ nN A $$^{-3/2}$$ , $$1.96<K_{{\text {eq-gy}}}$$ (graphyne) $$<2.49$$ nN A $$^{-3/2}$$ and $$1.42<K_{{\text {eq-su}}}$$ (super) $$<2.95$$ nN A $$^{-3/2}$$ using the FE method. These findings should be of great help for understanding the fracture properties of carbon allotropes and designing the carbon-based nanodevices.

Refined multiscale model based on the second generation interatomic potential for the mechanics of graphene sheets

Computationally efficient multiscale model based on the refined constitutive law, including second generation reactive empirical bond order potential with dihedral energy term, is employed to investigate the static response of the graphene sheets under transverse and in-plane compressive loads including material nonlinearity through atomic interactions and Green-Lagrange geometric nonlinearity through the strain displacement relations. The bending modulus of the graphene sheet predicted without considering the dihedral energy term in the constitutive law is almost half than those of predicted through first principle calculations. The inclusion of the dihedral energy term predicts bending modulus close to those of through first principle calculations. The atomistic and continuum deformations are coupled through the Cauchy–Born rule. In the present study, the effect of the dihedral energy term on the linear and nonlinear bending and postbuckling response of the graphene sheets under transverse and in-plane compressive loads is investigated in detail. The governing finite element equations for the graphene sheet are derived through the principle of minimum potential energy. The spatial approximation of the graphene sheet at the continuum scale is attained through the finite element method.

Accelerating spin-space sampling by auxiliary spin dynamics and temperature-dependent spin-cluster expansion

Atomistic simulations of thermodynamic properties of magnetic materials rely on an accurate modelling of magnetic interactions and an efficient sampling of the high-dimensional spin space. Recent years have seen significant progress with a clear trend from model systems to material specific simulations that are usually based on electronic-structure methods. Here we develop a Hamiltonian Monte Carlo framework that makes use of auxiliary spin-dynamics and an auxiliary effective model, the temperature-dependent spin-cluster expansion, in order to efficiently sample the spin space. Our method does not require a specific form of the model and is suitable for simulations based on electronic-structure methods. We demonstrate fast warm-up and a reasonably small dynamical critical exponent of our sampler for the classical Heisenberg model. We further present an application of our method to the magnetic phase transition in bcc iron using magnetic bond-order potentials.

A Molecular Dynamics Based Comparison of the Mechanical Properties of Three Polytypes of Cubic BC3

In this work molecular dynamics simulations are employed to compare the mechanical properties and hardness of three polytypes of cubic BC3. Firstly, two interatomic Tersoff potentials, with different parameterizations, were identified from the literature for the boron-carbon system. Based on cohesive energies and structural properties from existing ab-initio studies, the suitability of the two potentials for predicting the properties of BC3 was analyzed. Secondly, using the better interatomic potential, more detailed molecular dynamics simulations were conducted to estimate and compare the elastic, yield, post-yield behavior and hardness of the three polytypes. The elastic constants compare well with existing abinitio values and vary by at most by 15% amongst the three polytypes. Response to indentation showed considerable qualitative differences in yield and post-yield response. One of the polytypes showed lower yield strength and seemed more ductile than the other two. The hardness also showed a complex dependence on both the material and the indentation depths. A peculiar, indenter-size dependent pile-up behavior was also seen. Specifically, for lower radii, pile-up was seen on indentation. As the radius of the indenter was increased, pile-up was seen only on retracting the indenter. The higher volume occupied by the indentation-amorphized material was found to be the reason for pile-up on retracting the indenter.