ReaxFF Reactive Force-Field Study of Molybdenum Disulfide (MoS2

Letter pubs.acs.org/JPCL

ReaxFF Reactive Force-Field Study of Molybdenum Disulfide (MoS2) Alireza Ostadhossein,† Ali Rahnamoun,‡ Yuanxi Wang,§ Peng Zhao,† Sulin Zhang,† Vincent H. Crespi,§,∥,⊥ and Adri C. T. van Duin*,‡ †

Department of Engineering Science and Mechanics, Pennsylvania State University, University Park, Pennsylvania 16802, United States ‡ Department of Mechanical and Nuclear Engineering, Pennsylvania State University, University Park, Pennsylvania 16802, United States § Department of Physics, Pennsylvania State University, University Park, Pennsylvania 16802, United States ∥ Department of Chemistry, Pennsylvania State University, University Park, Pennsylvania 16802, United States ⊥ Department of Materials Science and Engineering, Pennsylvania State University, University Park, Pennsylvania 16802, United States S Supporting Information *

ABSTRACT: Two-dimensional layers of molybdenum disulfide, MoS2, have been recognized as promising materials for nanoelectronics due to their exceptional electronic and optical properties. Here we develop a new ReaxFF reactive potential that can accurately describe the thermodynamic and structural properties of MoS2 sheets, guided by extensive density functional theory simulations. This potential is then applied to the formation energies of five different types of vacancies, various vacancy migration barriers, and the transition barrier between the semiconducting 2H and metallic 1T phases. The energetics of ripplocations, a recently observed defect in van der Waals layers, is examined, and the interplay between these defects and sulfur vacancies is studied. As strain engineering of MoS2 sheets is an effective way to manipulate the sheets’ electronic and optical properties, the new ReaxFF description can provide valuable insights into morphological changes that occur under various loading conditions and defect distributions, thus allowing one to tailor the electronic properties of these 2D crystals. first-principles computational approaches such as density functional theory (DFT) provide indispensable data about the electronic structure of 2D systems at a quantum mechanical level, they are limited to relatively small systems (several hundreds of atoms) and short time scales (picoseconds). Although the dynamics of 2D materials have been widely investigated with ab initio molecular dynamics (AIMD),24,25 to enable longer time (nanoseconds and beyond) and larger scale (≫1000 atoms) MD simulations an accurate force field describing the interatomic interactions is required, and to model synthesis this force field needs to be reactive. To date, several force fields have been proposed for MoS2. A valence force field (VFF) was first described in 197526 and used to study the dynamical properties of MoS2 lattices and nanotubes.27−29 The intralayer interactions were assumed to be associated with the stretching and bending of Mo−S bonds,

T

he isolation of a single-layer graphene from graphite1 was a pivotal moment in research on 2D materials. Because of its remarkable electronic, optoelectronic, and chemo-mechanical properties, graphene has become a central material for nanoscale transistors,2−4 photovoltaics,5−9 optical sensors, and topological field-effect transistors.10,11 Although graphene is currently the most widely studied 2D material, pristine graphene does not have a band gap, a property that is essential for many classes of electronic devices. Transition-metal dichalcogenides (TMDs) such as (Mo,W)(S,Se)2 have been recognized as promising alternatives due to their band gap of ∼1.2 eV in the bulk and ∼1.8 eV direct gap in single-layer (SL) MoS2.12 A MoS2 layer also shares some of graphene’s desirable properties, such as mechanical strength and electrical conductivity, opening many other possible applications such as photodetectors.13−20 Despite tremendous progress in the synthesis of TMD materials, the costly and relatively slow synthesis/characterization/test cycle for these devices has motivated the development of robust multiscale computational models that can help guide design and fabrication of 2D devices.21−23 While © 2017 American Chemical Society

Received: December 10, 2016 Accepted: January 19, 2017 Published: January 20, 2017 631

DOI: 10.1021/acs.jpclett.6b02902 J. Phys. Chem. Lett. 2017, 8, 631−640

Letter

The Journal of Physical Chemistry Letters

and the possibility to use these defects to alter or remove other defects is investigated. To fit the ReaxFF parameters for Mo/S interactions, we began with the previously published Mo/S/Ni/ O parameter set by Vasenkov et al.,45 replacing the sulfur (S) parameters with a newer set of atomic parameters from the recent publication by Islam et al.42 We then optimized these parameters against the DFT data set; Table S4 lists the resulting ReaxFF force-field parameters. The DFT data set describes the dissociation of Mo(S)−S bonds as well the S− Mo−S, Mo−S−Mo, S−S−Mo, H−S−Mo, and H−Mo−S angle bending in several small molecules (see Figure S1). The formation energy (Ef) (Figure S2) and the equations of state (EOS) (Figure S3) of condensed-phase single and double MoS2 crystalline structures are also included in the training set. To properly capture bond energies, we performed quantum mechanical (QM) calculations on Mo−S and S−S bond dissociation in various molecular clusters such as MoS4H2 and MoS3H2 (see Figure S1), each relaxed into a ground-state structure. As shown in Figure S1b, the ReaxFF bond dissociation energies for both S−S and Mo−S, achieved by scanning from the equilibrium distance to the dissociation limit, are within 5−10 kcal/mol of the QM results. The equilibrium bond lengths of single- and double-bonded Mo−S within ReaxFF and QM are shown in Figure S1(e,f). Good agreement is also achieved for the dissociation of the Mo−S bond in the MoS4H2 molecule. To optimize the angle parameters, S MoS and S−MoS bond-angle energetics were calculated in the selected molecules by varying the angle in question from 70−160° with the other parameters held fixed. To further enrich the training set, reaction energies and heats of formation associated with the Mo−S- and Mo−H-based compounds shown in Figure S2 were also incorporated. Figures S1 and S2 show that ReaxFF successfully reproduces QM equilibrium angles and the overall dissociation energy curves. The uniaxial and biaxial expansion and contraction of the MoS2 crystal both within the basal plane (a direction) and out-of-plane were also incorporated into the training set. The details of these simulations are presented in Supporting Information, Figures S3. The experimental structural parameters of free-standing SL MoS2 and bulk MoS2 are given in Table 1 and compared

and the interlayer interaction between S atoms was described as a sum of spherically symmetric two-body potentials30 between S atoms of neighboring layers. More recently, starting from the VFF model, Jiang et al.31 developed a Stillinger−Weber (SW) potential for SL MoS 2 by fitting parameters to the experimentally phonon spectrum obtained from inelastic neutron scattering and reproduced the results at the DFT level. In the SW potential, geometrical parameters are determined analytically based on the equilibrium state of each individual potential term, while the energy parameters are taken from the valence force-field model.32 The simplicity of this harmonic VFF model makes it suitable for the analytical derivation of the elastic quantities at low computational cost. Such potentials are useful for describing phonons and elastic deformations. However, they have difficulty describing states far from equilibrium.31 The Tersoff33 and Brenner34,35 models, also known as reactive empirical bond order (REBO) potentials, have had particular success for carbon-based materials and have also been adapted for MoS2.36 By explicitly incorporating the dependence of bond order (i.e., the strength of each bond) on local environment, the resulting potentials permit the calculation of structural properties and energetics for complex systems, leading to an improved description of covalent materials. However, these potentials are more computationally expensive than nonreactive FFs such as SW models. The second-generation bond-order-based potential model ReaxFF allows for covalent bond breaking and forming with associated changes in atomic hybridization. The ReaxFF reactive potential, which was developed to simulate large-scale reactive chemical systems,37 also eschews explicit bonds in favor of a bond-order concept. Because the connectivity-dependent interactions in this force field have been formulated based on bond order, a smooth transition can be achieved between nonbonded and bonded states while also retaining significant bond orders at a transition state. A detailed description of the ReaxFF formulation can be found in the original work.37 Contrary to the previous reactive potentials (i.e., Tersoff,33 Brenner,35 and REBO34), ReaxFF employs a dynamic charge scheme to vary the charge on atoms according to their chemical environment. Self-consistent charge equilibration based on the electronegativity of atoms (called EEM38) enables ReaxFF to determine atomic charge distributions consistent with explicit local electrostatics. ReaxFF includes van der Waals interactions, which play a central role in the modeling of multilayer MoS2. Also, ReaxFF employs a significantly longer-ranged bond order relationship compared with previous reactive force-field formalisms, which enables improved energies and structures for transition states and thus more precise reaction kinetics.39 This model has been extensively used in a wide range of systems.25,40−44 The ReaxFF framework allows the same atomic parameters for an element to be used regardless of its chemical environment. Among other advantages, this transferability permits atoms to migrate properly between phases during an MD simulation.39 Here we develop a ReaxFF parameter set for Mo and S to study energetics and reaction mechanisms in single-layer and multilayer (SL and ML) MoS2. The ReaxFF parameters are fit to a training set of energies, geometries, and charges derived from DFT calculations for both clusters and periodic systems. We then employ this potential (see Supporting Information) to calculate vacancy formation energies (VFEs) and the bending rigidity (BR) of pristine and defective SL MoS2 films. A new class of defects in van der Waals layers, ripplocations, is studied,

Table 1. Structural Properties of SL MoS2 As Depicted in Figure 1 structural property (Å)

ReaxFF

experiment

DFT

|a1| = |a2| d (Mo−S) h

3.19 2.41 3.11

3.157 2.38 3.116

3.16 2.43 3.11

against the DFT and ReaxFF versions of the same quantities. As seen in Figure 1a, SL MoS2 consists of a Mo atomic plane sandwiched between two planes of S atoms, each Mo being at the center of the triangular prism formed by six neighboring S atoms. The equilibrium lattice parameters within ReaxFF (corresponding to the minimum of the EOS curves) are a = 3.19 Å and c = 12.24 Å, in reasonable agreement with the experimental values of a = 3.16 Å and c = 12.32 Å.46 The vertical separation between Mo atoms in adjacent sheets is calculated as dMo−Mo = 6.112 Å, in good agreement with both DFT (6.1 Å) and experiment (6.2 Å47) for bulk MoS2. The lattice parameters, S− S separation (h), and Mo−S bond distance are calculated based 632


Letter


behavior of SL MoS2 can be described by a Taylor series expansion of the strain energy density U (εij) in terms of powers of strain. To determine the elastic constants of SL MoS2 using ReaxFF, starting from the equilibrium condition (U (ε = 0) = 0 and ∂U/∂εij (ε = 0) = 0), we apply uniaxial/biaxial stretch by step wisely increasing (decreasing) the computational cell size within the basal plane. Periodic boundary conditions were imposed along planar directions to eliminate edge effects. The elastic energy density (per unit area) can be expressed as U (ε) =

1 c11(εxx2 + εyy2 + 2εxy2 ) + c12(εxxεyy − εxy2 ) 2

(1)

where x and y correspond to the zigzag and armchair directions in the hexagonal lattice of MoS2, respectively. The isotropic linear elastic properties (i.e., surface elastic modulus Ys and Poisson’s ratio ν) can be obtained by the following expressions Ys =

2 2 − c12 c11 c and ν = 12 c11 c11

(2)

Table 2 summarizes and compares the calculated values of Ys and ν as well as the other elastic constants of MoS2, as defined by eq 2. The ultimate tensile strengths, σu, of SL MoS2 along the zigzag and armchair directions, as illustrated in Figure 1, were calculated to be 24.6 and 26.2 GPa, respectively. These values are the maximum points of stress−strain curves shown in Figure S4. It should be noted that the reported values in Table 1 are for the 2H phase of MoS2. These calculated Young’s moduli are in good agreement with the atomic force microscopy and nanoindentation experiments.57,61−63 A comprehensive study on the mechanical properties of pristine (2H and 1T phase) and parallel heterostructures of MoS2 using the present ReaxFF potential has been presented in ref 64. Phase Transition in SL MoS2 (2H ↔ 1T). It is well understood that the 2H-phase of MoS2 is the most stable phase exhibiting semiconducting behavior, while the metallic 1T-phase is only metastable. In a recent in situ scanning TEM study, Lin and coworkers65 showed that phase transitions induced by electron beam irradiation from the 2H phase to the 1T phase may produce coherent planar heterostructures of MoS2 with both phases present. Here we use ReaxFF to investigate the energetics of such polymorphic SL MoS2 heterostructures by computing the minimum energy pathways (MEPs), displayed in Figure 2, for the 1T ↔ 2H transformation in a homogeneous SL MoS2. The initial configuration, labeled by A, and final configuration, denoted by D, are shown in Figure 2. The transition state is mapped out by rigidly sliding the top S atomic layer with respect to the rest of the

Figure 1. Top and side views of SL MoS2. The lattice vectors (a1 and a2) defining the unit cell and Mo−S bond distance in 2H-MoS2 obtained from ReaxFF calculations.

on ReaxFF, and their corresponding values are compared with the results from DFT46 and experiment48 in Table 1. The extensively trained and tested ReaxFF/DFT bond and angle parameters, together with the diversity of the condensedphase structures in the training set, demonstrate the accuracy of the new force field to model other Mo−S containing systems. Force-field parameters are provided in the Supporting Information (Table S4). Structural and Mechanical Properties. While the fracture of graphene has been intensively studied,49−54 the fracture behavior of other 2D layered systems is yet to be explored.55,56 This is partly due to the lack of accurate force fields, which can reliably model bond breaking upon fracture. Because of its trilayer structure and binary atomic composition, the fracture mechanism of MoS2 is expected to be more complicated than that of graphene. In this section, mechanical and structural properties of SL MoS2 are examined through ReaxFF-based MD simulations. Comparison of the simulation results with experimental data57 and previous DFT simulations57−60 further validates the new ReaxFF description. The linear elastic

Table 2. Non-Zero Independent Components for the Elastic Tensor, Poisson’s Ratio ν, and In-Plane Stiffness Ys of 2H-MoS2 from DFT Calculationsa,b method Ys (N/m) υ c11 (N/m) c12 (N/m) c22 (N/m)

ReaxFF

ref 60

ref 59

ref 59

ref 58

experiment57

ref 46

176.32 0.39 205.08 81.59 217.13

GGA 120.1 0.254 128.4 32.6 128.4

LDA 118 0.31 130 40 130

GGA 129 0.29 140 40 140

TM 123 0.25

AFM 180 ± 60 0.27 (bulk)

HSE06-D2 134.60 145.0 82.90

a

Nominal thickness of MoS2 layers is calculated as 6.09 Å. bAFM: atomic force microscopy; TM: Trouiller-Martins; LDA: local density approximation; GGA: generalized gradient approximation; HSE: Heyd, Scuseria, and Ernzerhof. 633


Letter


Figure 2. (a) Minimum energy paths (MEPs) for the 1T → 2H phase transformation of single-layer MoS2 are calculated within both DFT (red circles) and ReaxFF (black squares), with insets showing top and side views of the two end-state structures. The system transitions from the 1T phase (configuration A) through intermediate structures B and C before ending at the 2H phase (configuration D).

Figure 3. Atomic structures and formation energies of ripplocations with different numbers n of extra units added to the top layer.

2H phase is ∼0.5 eV per formula unit more favorable than the 1T phase, while our DFT calculations yield an energy difference of ∼0.75 eV. (In both cases, the 1T-MoS2 to 2H-MoS2 transition is exothermic.) As shown in Figure 2a, the reverse energy barrier for 2H-to-1T is ∼0.7 eV within DFT, as compared with 0.45 eV with ReaxFF. The discrepancy (∼0.25 eV) between the ReaxFF and DFT computed barriers arises from some combination of different samplings of the transition path and different computational methodologies. Overall, the ReaxFF energetics of this transition are down-scaled from the DFT results by a constant factor of roughly ∼0.65. Although the discrepancy necessitates careful attention to the accuracy of the potential when modeling the absolute rates of activated processes, the errors introduced into relative rates of different kinetic processes may be more modest, corresponding to ReaxFF-based simulations at appropriately rescaled temperatures. Formation Energies of Ripplocations on the Surface of MoS2 Layers. A recent observation by Kushima et al.68 using highresolution TEM (HRTEM) revealed the existence of localized ripples in 2D crystals. The formation of these ripples involves slippage of an upper 2D layer against a lower layer by an inplane Bravais vector without breaking covalent bonds. The

crystal through climbing image nudged-elastic band (CI-NEB) sampling at the DFT level using VASP and NEB at the forcefield level using the new ReaxFF potential in LAMMPS. Our DFT results for the 1T-to-2H transition barrier in SL MoS2 are in good agreement with previously reported DFT results (1.04 eV per formula unit) by He et al.66 For the purposes of comparing the reaction pathways, both ReaxFF and DFT energetics have been referenced to the energy of the 2H phase in Figure 2. As depicted in Figure 2, the 1T-MoS2 phase in which Mo atoms have distorted octahedral coordination by sulfur atoms is energetically less favorable than the 2H-phase with a trigonal prismatic coordination of Mo. The gliding motion of the upper layer of S atoms along the transition path is consistent with the previously reported DFT results and TEM images of Guo et al.67 Note that the plotted reaction barrier is an energy per formula unit for a spatially homogeneous transformation. The actual reaction barrier in the physical transformation will depend in a complex fashion on the spatial extent of the transition region between the two phases; nevertheless, the modest size of the barrier per formula unit does suggest that the actual transition can be thermally activated under relatively modest thermal conditions. Within our ReaxFF treatments, the 634


Letter


Figure 4 compares the calculated defect formation energies of these five point defects within ReaxFF with their DFT

insertion of an extra row of unit cells into one layer relative to another creates ripples to relax compressive stress, where contact between the two layers is then lost in a domain of finite width that accommodates the extra strip of material.68 The formation mechanism of these ripples on the surface of van der Waals structures is morphologically similar to conventional edge dislocations in bulk systems. Because the newly developed ReaxFF potential is trained against the energies of highly curved MoS2 sheets, we expect it to be able to accurately describe the behavior of the ripplocations. As shown in Figure 3, ripplocations were generated by inserting n extra rows of MoS2 unit cells into the top layer, corresponding to a Burger’s vector of b = nb0 (n = 1, 2, 3, ...) and then structurally relaxing the resulting system. The formation energies of these ripples, Ef (n), are calculated using both ReaxFF and nonreactive SW potentials and compared against the DFT-vdW simulation results presented in figure 2 of ref 68. All of the initial bent structures shown in Figure 3a consist of 20 units of MoS2 in the bottom layer, that is, u = 20, with the number of units in the top layer being u + n = 21, 22, and 23. Consistent with the corresponding DFT values, ReaxFF predicts that the formation energy of these ripples varies sublinearly with respect to the magnitude of the Burger’s vector. This sublinear behavior suggests that ripplocations tend to attract and merge whenever they are in close vicinity, similar to the behavior seen in membrane mechanics.68 This is in clear contrast with the repulsive interactions between two same-sign edge dislocations. Figure 3b quantifies this sublinearity, revealing a dependence of n0.31 for ReaxFF, as compared with n0.43 from the DFT simulations. The nonlinearity for the SW potential is somewhat weaker, Ef (n) ≈ n0.5. A continuum mechanics analysis approximates this exponent to be 0.33, in close agreement with the ReaxFF result.68 Vacancy Formation Energies and Dif f usion Barriers in SL MoS2. Because the structure of SL MoS2 consists of three layers of atoms, point defects are more complex than those in SL graphene. Six types of point defects (including vacancies and antisite defects) have been identified and characterized by imaging at atomic resolution.69 Along with DFT simulations, TEM studies revealed that the electronic properties of SL MoS2 strongly depend on the structure of these defects. However, the lack of an accurate reactive potential that can reliably predict the structural and mechanical properties of defective MoS2 hinders further progress, particularly for phenomena that occur on large length and time scales. Thus we next examine the energetics of defect formation and transformation as a further verification of the new ReaxFF force field. The stability of MoS2 layers containing point defects of various types is examined through the associated formation energies. Ef defect = E defected − Epristine −

∑ i = S,Mo

(μi Ni)

Figure 4. Formation energies of point defects in monolayer MoS2, comparing ReaxFF with DFT.

counterparts from Zhou et al.,69 showing excellent agreement. Also, comparison between the experimental TEM images of ref 69 and the minimized structures from our ReaxFF simulations indicate a reasonable agreement. For the assumed values of chemical potentials, VS is the most energetically favorable defect, with the double S vacancy costing about twice as much in terms of formation energy. The ReaxFF calculated chemical potentials for S and Mo atoms are 2.82 and 6.90 eV/atom, in close agreement with the experimental values of 2.85 eV/atom for S and 6.82 eV/atom for Mo.70 To further validate the new ReaxFF potential, we next examined the migration barriers against motion of VS within SL MoS2. Recently, Le et al.71 have calculated within DFT the vacancy migration energy barriers for hopping of a S atom, shown in Figure 5 with orange, along a line of two adjacent S vacancies. Following this study, we considered two different paths (Paths 1 and 2 in Figure 5) for moving a S atom along the line of two vacancies. Path 1 considers motion of the adjacent sulfur that is coaligned with the vacancy axis, while Path 2 moves a sulfur atom that is oriented 60° off this axis. In both cases, we consider motion just between two adjacent sulfur lattice sites (i.e., a single “hop”). For Path 1, the total energies of initial state (panel I) and final state after migration (panel V) are very similar. To ensure convergence, we examined different numbers of replicas in the NEB transition path, obtaining very similar energy barriers. Table 3 compares the results of our calculations with the DFT barriers of Le et al.71 The ReaxFF-calculated diffusion barriers are within 0.13 eV of the DFT values, and in both methodologies diffusion along Path 2 has a barrier about ∼1.0 eV lower than that along Path 1. Bending Rigidity and Kinetics of Single Vacancy (VS) Defect Formation in SL and ML MoS2. The unique morphology of ripplocations, along with their high mobility and thermal stability,68 suggests the idea of using these structural perturbations to mobilize atomic defects such as vacancies and adatoms by exploiting the dependence of defect formation energy on sheet curvature. Here the BR or flexural rigidity of freestanding MoS2 layers is a critical parameter in determining the population and dynamics of defects in buckled structures. Lu et al.72 determined the elastic energy of curvature in SL graphene by calculating the formation energies of single-walled

(3)

where Ni terms are the differences in the number of S and Mo atoms in the defected and pristine structures and μi terms (i = S or Mo) represent the chemical potentials of S and Mo atoms referenced to alpha-S8 and the Mo bcc lattice, respectively. We study the VFEs of a monosulfur vacancy (VS), disulfur vacancy (VS2), the absence of one Mo atom and its three adjacent sulfur atoms (VMoS3), the same with six missing sulfur atoms (VMoS6), and the antisite defect that replaces a S2 dimer with a Mo atom (S2Mo). 635


Letter


Figure 5. (a) (I−V) Five sequential atomic configurations for an S atom diffusing along a line of two adjacent S vacancies. The hollow red circles in panels a and b indicate the missing S atoms on the top layer of SL MoS2. (b,c) Hopping of the marked (i.e., orange) S atom along the horizontal “Path 1” is presented along with the energy barrier for a 9-replica path. Yellow and cyan balls denote S and Mo atoms, respectively. The Mo atoms underneath the S vacancies row are dark blue. Path 2 is also considered in the main text.

the axial direction is free to optimize. With pure bending deformation, the strain energy density of the tube depends only on the tube curvature, as U(κ) = 1/2Dκ2. The ReaxFF treatment yields D ≈ 9.91 eV, reasonably consistent with results from a SW potential, 9.6 eV.73 To probe the effect of surface functionalization or vacancies on the BR of MoS2 (or equivalently, the dependence of defect formation energy on sheet curvature) we applied the same procedure to MoS2 tubes whose outer surface was fully covered by hydrogen (HMoS2) or weakened by a line of S-vacancies (on either the outer or inner outer surface). The strain energy density of pristine and hydrogenated SL MoS2 (U) as a function of bending curvature κ is plotted in Figure 6. Covalent functionalization of the outer surface breaks local reflection symmetry about the sheet midplane and thus shifts of the minimum energy toward nonzero curvature. Figure 7 shows the VFE for a line of VS defects on the inner or outer surfaces of a MoS2 nanotube as the curvature varies from 0.01 to 0.1 Å−1. Contrary to the case of SWCNT where vacancy

Table 3. Diffusion Barriers of Sulfur Atom near Vacancies in SL MoS2 path

ReaxFF (eV)

DFT (eV)71

Path 1 Path 2

2.44 1.48

2.32 1.35

carbon nanotubes (SWCNTs) with various diameters and found a BR of D = 0.83 eV. Following their method, a series of ReaxFF MD simulations were carried out for MoS2 monolayers rolled into cylinders of radii 1−10 nm to determine the strain energy, as shown in Figure 6. The BR of SL MoS2 can be

Figure 6. Strain energy density (strain energy per unit area) of SL MoS2 tubes as a function of bending curvature.

expressed as the second derivative of the strain energy density with respect to tube curvature, that is, ∂2U/∂κ2. In each case, the minimum energy configuration is obtained from a combination of conjugate gradient (CG) minimization and dynamically damped minimization, using the so-called “fire” method, at T = 0.1 K. To mimic pure bending of a MoS2 layer, the total potential energy is minimized when the length of the periodic box along

Figure 7. Curvature-mediated vacancy formation in SL MoS2. Whereas the formation energy of a sulfur monovacancy in the inner surface decreases with increasing curvature, similar to the behavior seen for vacancies in carbon nanotubes, the formation energy of a sulfur vacancy in the outer surface increases at higher curvature. 636


Letter


Figure 8. Formation energy and diffusion barrier for a sulfur vacancy on the upper or lower surfaces of an n = 3 ripplocation, as a function of defect location. (a) Side and top views of the ripplocation with vacancy on the upper surface. The VFE is computed for the vacancy sites shown with red circles. (b) Comparison of formation energies for upper and lower vacancies to the case of a flat SL of MoS2 as a function of distance away from the center of the ripplocation. (c,d) Diffusion barriers for straight motion of a vacancy on the upper surface (Path 1) and zigzag motion on either upper or lower surfaces (Paths 2 and 3). The red dashed line shows the vacancy diffusion barrier in flat SL MoS2.

path on upper- and lower-surfaces, respectively. Figure S5 in the Supporting Information shows the variation in energy along the full reaction path as a function of the distance from the center of the ripplocation (x). Vacancy hopping is clearly easier along the zigzag paths and is further facilitated when the vacancy is situated within a local concavity: near the center of the ripplocation for the lower vacancy (∼1.2 eV barrier) and near the edge for the upper vacancy (∼1 eV barrier). The largest diffusion barriers encountered are 5.4, 1.89, and 1.94 eV for Paths 1, 2, and 3 respectively; these occur when the vacancy sits in a local convexity. Our interest here is not in the kinetics of a vacancy diffusing across a stationary ripplocation, but the kinetics of a vacancy that has been entrained by a moving ripplocation. Consider a ripplocation that slides at a constant speed v along a MoS2 layer that hosts some areal density of sulfur vacancy defects. A vacancy on the upper surface of the layer will prefer to sit toward the leading edge of the ripplocation, that is, the site of lowest VFE. If the time for the ripplocation to slide by one sulfur lattice site is comparable to the time scale for the sulfur vacancy to hop between these sites (more precisely, the drift rate for the biased random walk of the vacancy down the local slope of VFE), then the vacancy can “surf” the ripplocation and thus be transported along the layer. Such a process has obvious utility in creating more nearly defect-free regions for use in high mobility devices. One important time scale here is the inverse of the hopping rate τ−1 = ωe−(ΔE/kT), where ω ≈ 1012 s−1 is a characteristic vibrational frequency and ΔE is the relevant hopping barrier. Fortuitously, the hopping barriers tend to be smallest (1 to 1.2 eV) where the VFEs are lowest (see Figure 8), that is, near the locations where the vacancies would be most likely entrained. For example, at a temperature of 1000 K one obtains τ ≈ 10−6 s; that is, the ripplocation can move at up to roughly one million vacancy hop lengths per second (i.e., ∼0.1 mm/s) while keeping the vacancies entrained. This crude analysis assumes perfect bias in the vacancy hopping: The actual maximum “surfing speed” is certainly lower than this, but that reduction in

formation is always easier on surfaces with higher curvature (and the vacancy does not break reflection symmetry about the midplane), the formation energy of VS in MoS2 tubes depends not only on the magnitude of the curvature but also on its sign, that is, whether the VS is on the convex or concave side of the layer. Vacancies on the outer (convex) surface are more expensive at higher curvatures, while vacancies on the inner (concave) surface have lower formation energy at higher curvature. In other words, vacancy formation on the inner surface is facilitated by sheet curvature, while that on the outer surface is suppressed. Both systems approach the flat-sheet result as the radius goes to infinity. This is consistent with the physical intuition that the vacancy “gap” in the layer is partially compensated by the bending inward of the neighboring sulfur atoms for the inner-surface vacancy. Defect formation regulated by sheet curvature might potentially be used in conjunction with the high mobility of ripplocations to sweep defects along surfaces, expelling them from selected regions by analogy to the zone refining of semiconductors. Figure 8b demonstrates the energetic interplay between sulfur vacancies and ripples. As a function of position along the n = 3 ripplocation, we calculate the formation energy for both upper-surface and lower-surface vacancies. As anticipated by the nanotube results, vacancies are most favorable in regions of high curvature, with upper and lower surface vacancies favored at the edge or center of the ripple, respectively (i.e., where the vacancy is situated in a local concavity). The energy required to move a vacancy from one side of the sheet to the other thus varies with both the magnitude and sign of the curvature. To examine the mobility of the vacancies along the ripple, we calculated the migration energy landscape for lower- and upper-surface vacancies, as shown in Figure 8c,d. Three possible diffusion pathways, each consisting of series of sulfur vacancy jumps between neighboring local minima, were considered. Path 1, Figure 8c, follows a straight line along the ripple, while Paths 2 and 3 (Figure 8d) correspond to vacancy migration along a zigzag 637


Letter

The Journal of Physical Chemistry Letters efficacy can be compensated for by multiple passes of ripplocations over the area of interest. At lower temperatures the maximum speed is lower, but the bias in hopping is larger. The passage time of a ripplocation past a given location is its width (∼100 Å) divided by its speed; this could plausibly be in the kilohertz range or above. Thus it appears plausible that induction of a repeated series of propagating ripplocations over a given region of a MoS2 layer could potentially act to sweep out sulfur vacancies from that region, given a somewhat elevated temperature and a sufficient length of treatment. Because these effects are geometrical in nature, similar effects could be anticipated for other TMD systems and potentially other defect types. In summary, we have developed a new ReaxFF reactive potential. The favorable comparison of the ReaxFF results with available DFT and experimental values confirms the capability and transferability of the new potential to investigate the properties of MoS2 nanostructures, including mechanical failure, the dynamics of rippolocations, and defect formation energies. The interaction of ripplocations and defects suggests an intriguing possibility to sweep out undesirable defects from regions of interest with a pulse train of propagating ripplocations.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. ORCID

Alireza Ostadhossein: 0000-0003-2820-9749 Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This work was supported by the U.S. National Science Foundation (NSF) grant no. CMMI-1462980. Y.W. and V.H.C. acknowledge support from U.S. Army Research Office MURI grant W911NF-11-1-0362 and the National Science Foundation Materials Innovation Platform under DMR1539916.

■

COMPUTATIONAL METHODS Ab Initio Calculations. DFT calculations were carried out by full optimization of the molecules listed in Table S1, using the Jaguar program with B3LYP exchange-correlation functional and LAV3P** basis sets for all atoms. Further DFT calculations were performed on periodic MoS2 structures using the Vienna Ab initio Simulation Package (VASP) in the generalized gradient approximation as parametrized by Perdew−Burke− Ernzerhof (GGA-PBE), using the projector augmented waves (PAW) method and Brillouin-zone integration on a Γ-centered Monkhorst−Pack grid. In all simulations, the plane-wave energy cutoff was 400 eV and atomic coordinates were relaxed until the forces on all atoms were below 0.02 eV/Å. Structural relaxation of bilayer MoS2 at different interlayer distances d was achieved by constraining the Mo atoms in two adjacent layers to a c-axis separation of d and only allowing the S atoms to relax. AA′ stacking (where Mo and S in neighboring layers eclipse each other when viewed along the c axis) is preserved in all bilayer structures. Activation barriers were calculated by the nudged elastic band (NEB) method. Despite the bandgap underestimation problem inherent to standard DFT, this error should not affect structural energies discussed in the present study where no electronic excitation is present. ReaxFF Potential for MoS2. The energy in ReaxFF is written as the sum of the bond energy (Ebond), valence-angle (Eval), torsion-angle energy (Etor), lone pair (Elp), and over-coordinate (Eover) and under-coordinate (Eunder) energy penalties plus the nonbonded van der Waals (Evdw) and Coulomb (ECoulombic) contributions

■

REFERENCES

(1) Novoselov, K. S.; Geim, A. K.; Morozov, S. V.; Jiang, D.; Zhang, Y.; Dubonos, S.; Grigorieva, I. V.; Firsov, A. A. Electric Field Effect in Atomically Thin Carbon Films. Science 2004, 306, 666−669. (2) Britnell, L.; Gorbachev, R. V.; Jalil, R.; Belle, B. D.; Schedin, F.; Mishchenko, A.; Georgiou, T.; Katsnelson, M. I.; Eaves, L.; Morozov, S. V.; et al. Field-Effect Tunneling Transistor Based on Vertical Graphene Heterostructures. Science 2012, 335, 947−950. (3) Schwierz, F. Graphene Transistors. Nat. Nanotechnol. 2010, 5, 487−496. (4) Choi, W.; Cho, M. Y.; Konar, A.; Lee, J. H.; Cha, G.-B.; Hong, S. C.; Kim, S.; Kim, J.; Jena, D.; Joo, J.; et al. High-detectivity Multilayer MoS2 Phototransistors with Spectral Response from Ultraviolet to Infrared. Adv. Mater. 2012, 24, 5832−5836. (5) Feng, J.; Qian, X.; Huang, C.-W.; Li, J. Strain-Engineered Artificial Atom as a Broad-Spectrum Solar Energy Funnel. Nat. Photonics 2012, 6, 866−872. (6) Bernardi, M.; Palummo, M.; Grossman, J. C. Extraordinary Sunlight Absorption and One Nanometer Thick Photovoltaics Using Two-Dimensional Monolayer Materials. Nano Lett. 2013, 13, 3664− 3670. (7) Pospischil, A.; Furchi, M. M.; Mueller, T. Solar-Energy Conversion and Light Emission in an Atomic Monolayer Pn Diode. Nat. Nanotechnol. 2014, 9, 257−261. (8) Shanmugam, M.; Jacobs-Gedrim, R.; Song, E. S.; Yu, B. TwoDimensional Layered Semiconductor/graphene Heterostructures for Solar Photovoltaic Applications. Nanoscale 2014, 6, 12682−12689. (9) Furchi, M. M.; Pospischil, A.; Libisch, F.; Burgdörfer, J.; Mueller, T. Photovoltaic Effect in an Electrically Tunable van Der Waals Heterojunction. Nano Lett. 2014, 14, 4785−4791. (10) Liu, J.; Hsieh, T. H.; Wei, P.; Duan, W.; Moodera, J.; Fu, L. Spin-Filtered Edge States with an Electrically Tunable Gap in a TwoDimensional Topological Crystalline Insulator. Nat. Mater. 2013, 13, 178−183. (11) Qian, X.; Liu, J.; Fu, L.; Li, J. Quantum Spin Hall Effect in TwoDimensional Transition Metal Dichalcogenides. Science 2014, 346, 1344−1347.

Esystem = E bond + Eval + E tor+ Eover + Eunder + E lp + E vdW + Ecoulombic

■

dissociation and valence angle distortion; Tables of reaction counts for Mo−S bond cleavage and formation reactions used for force-field optimization; Figures with examples of ReaxFF vs DFT energetics; uniaxial stress− strain responses of pristine 2H MoS2 sheet at room temperature; and evolution of vacancy migration energy and the optimized ReaxFF reactive force-field parameters. (PDF) initial configurations used for optimization of the potential in bgf format (PDF) trainset.in: training-set of the potential (PDF)

(4)

ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpclett.6b02902. List of parameters that were optimized within the FFtraining steps; DFT and ReaxFF calculations for bond 638


Letter


Stillinger-Weber Parametrization, Mechanical Properties, and Thermal Conductivity. J. Appl. Phys. 2013, 114, 064307. (32) Jiang, J.-W. Parametrization of Stillinger−Weber Potential Based on Valence Force Field Model: Application to Single-Layer MoS2 and Black Phosphorus. Nanotechnology 2015, 26, 315706. (33) Tersoff, J. Empirical Interatomic Potential for Carbon, with Applications to Amorphous Carbon. Phys. Rev. Lett. 1988, 61, 2879. (34) Brenner, D. W.; Shenderova, O. A.; Harrison, J. A.; Stuart, S. J.; Ni, B.; Sinnott, S. B. A Second-Generation Reactive Empirical Bond Order (REBO) Potential Energy Expression for Hydrocarbons. J. Phys.: Condens. Matter 2002, 14, 783. (35) Brenner, D. W. Empirical Potential for Hydrocarbons for Use in Simulating the Chemical Vapor Deposition of Diamond Films. Phys. Rev. B: Condens. Matter Mater. Phys. 1990, 42, 9458. (36) Liang, T.; Phillpot, S. R.; Sinnott, S. B. Parametrization of a Reactive Many-Body Potential for Mo−S Systems. Phys. Rev. B: Condens. Matter Mater. Phys. 2009, 79, 245110. (37) Van Duin, A. C.; Dasgupta, S.; Lorant, F.; Goddard, W. A. ReaxFF: A Reactive Force Field for Hydrocarbons. J. Phys. Chem. A 2001, 105, 9396−9409. (38) Mortier, W. J.; Ghosh, S. K.; Shankar, S. ElectronegativityEqualization Method for the Calculation of Atomic Charges in Molecules. J. Am. Chem. Soc. 1986, 108, 4315−4320. (39) Senftle, T. P.; Hong, S.; Islam, M. M.; Kylasa, S. B.; Zheng, Y.; Shin, Y. K.; Junkermeier, C.; Engel-Herbert, R.; Janik, M. J.; Aktulga, H. M.; et al. The ReaxFF Reactive Force-Field: Development, Applications and Future Directions. Npj Comput. Mater. 2016, 2, 15011. (40) Ostadhossein, A.; Kim, S.-Y.; Cubuk, E. D.; Qi, Y.; van Duin, A. C. Atomic Insight into the Lithium Storage and Diffusion Mechanism of SiO2/Al2O3 Electrodes of Li-Ion Batteries: ReaxFF Reactive Force Field Modeling. J. Phys. Chem. A 2016, 120, 2114. (41) Yoon, K.; Ostadhossein, A.; van Duin, A. C. Atomistic-Scale Simulations of the Chemomechanical Behavior of Graphene under Nanoprojectile Impact. Carbon 2016, 99, 58−64. (42) Islam, M. M.; Ostadhossein, A.; Borodin, O.; Yeates, A. T.; Tipton, W. W.; Hennig, R. G.; Kumar, N.; van Duin, A. C. ReaxFF Molecular Dynamics Simulations on Lithiated Sulfur Cathode Materials. Phys. Chem. Chem. Phys. 2015, 17, 3383−3393. (43) Kim, S.-Y.; Ostadhossein, A.; van Duin, A. C.; Xiao, X.; Gao, H.; Qi, Y. Self-Generated Concentration and Modulus Gradient Coating Design to Protect Si Nano-Wire Electrodes during Lithiation. Phys. Chem. Chem. Phys. 2016, 18, 3706. (44) Ostadhossein, A.; Cubuk, E. D.; Tritsaris, G. A.; Kaxiras, E.; Zhang, S.; van Duin, A. C. Stress Effects on the Initial Lithiation of Crystalline Silicon Nanowires: Reactive Molecular Dynamics Simulations Using ReaxFF. Phys. Chem. Chem. Phys. 2015, 17, 3832−3840. (45) Vasenkov, A.; Newsome, D.; Verners, O.; Russo, M. F., Jr; Zaharieva, R.; van Duin, A. C. Reactive Molecular Dynamics Study of Mo-Based Alloys under High-Pressure, High-Temperature Conditions. J. Appl. Phys. 2012, 112, 013511. (46) Peelaers, H.; Van de Walle, C. G. Elastic Constants and Pressure-Induced Effects in MoS2. J. Phys. Chem. C 2014, 118, 12073− 12076. (47) Xie, J.; Zhang, J.; Li, S.; Grote, F.; Zhang, X.; Zhang, H.; Wang, R.; Lei, Y.; Pan, B.; Xie, Y. Controllable Disorder Engineering in Oxygen-Incorporated MoS2 Ultrathin Nanosheets for Efficient Hydrogen Evolution. J. Am. Chem. Soc. 2013, 135, 17881−17888. (48) Wildervanck, J. C.; Jellinek, F. Preparation and Crystallinity of Molybdenum and Tungsten Sulfides. Z. Anorg. Allg. Chem. 1964, 328, 309−318. (49) Khare, R.; Mielke, S. L.; Paci, J. T.; Zhang, S.; Ballarini, R.; Schatz, G. C.; Belytschko, T. Coupled Quantum Mechanical/ molecular Mechanical Modeling of the Fracture of Defective Carbon Nanotubes and Graphene Sheets. Phys. Rev. B: Condens. Matter Mater. Phys. 2007, 75, 75412. (50) Pereira, V. M.; Neto, A. C.; Peres, N. M. R. Tight-Binding Approach to Uniaxial Strain in Graphene. Phys. Rev. B: Condens. Matter Mater. Phys. 2009, 80, 45401.

(12) Splendiani, A.; Sun, L.; Zhang, Y.; Li, T.; Kim, J.; Chim, C.-Y.; Galli, G.; Wang, F. Emerging Photoluminescence in Monolayer MoS2. Nano Lett. 2010, 10, 1271−1275. (13) Ataca, C.; Sahin, H.; Akturk, E.; Ciraci, S. Mechanical and Electronic Properties of MoS2 Nanoribbons and Their Defects. J. Phys. Chem. C 2011, 115, 3934−3941. (14) Lee, C.; Yan, H.; Brus, L. E.; Heinz, T. F.; Hone, J.; Ryu, S. Anomalous Lattice Vibrations of Single-and Few-Layer MoS2. ACS Nano 2010, 4, 2695−2700. (15) Molina-Sanchez, A.; Wirtz, L. Phonons in Single-Layer and FewLayer MoS2 and WS2. Phys. Rev. B: Condens. Matter Mater. Phys. 2011, 84, 155413. (16) Park, W.; Baik, J.; Kim, T.-Y.; Cho, K.; Hong, W.-K.; Shin, H.-J.; Lee, T. Photoelectron Spectroscopic Imaging and Device Applications of Large-Area Patternable Single-Layer MoS2 Synthesized by Chemical Vapor Deposition. ACS Nano 2014, 8, 4961−4968. (17) Popov, I.; Seifert, G.; Tománek, D. Designing Electrical Contacts to MoS 2 Monolayers: A Computational Study. Phys. Rev. Lett. 2012, 108, 156802. (18) Rao, C. N. R.; Gopalakrishnan, K.; Maitra, U. Comparative Study of Potential Applications of Graphene, MoS2, and Other TwoDimensional Materials in Energy Devices, Sensors, and Related Areas. ACS Appl. Mater. Interfaces 2015, 7, 7809−7832. (19) Windom, B. C.; Sawyer, W. G.; Hahn, D. W. A Raman Spectroscopic Study of MoS2 and MoO3: Applications to Tribological Systems. Tribol. Lett. 2011, 42, 301−310. (20) Yin, Z.; Li, H.; Li, H.; Jiang, L.; Shi, Y.; Sun, Y.; Lu, G.; Zhang, Q.; Chen, X.; Zhang, H. Single-Layer MoS2 Phototransistors. ACS Nano 2012, 6, 74−80. (21) Kumar, H.; Er, D.; Dong, L.; Li, J.; Shenoy, V. B. Elastic Deformations in 2D van Der Waals Heterostructures and Their Impact on Optoelectronic Properties: Predictions from a Multiscale Computational Approach. Sci. Rep. 2015, 5, 10872. (22) Kumar, H.; Dong, L.; Shenoy, V. B. Limits of Coherency and Strain Transfer in Flexible 2D van Der Waals Heterostructures: Formation of Strain Solitons and Interlayer Debonding. Sci. Rep. 2016, 5, 21516. (23) Khorshidi, A.; Peterson, A. Amp: A Modular Approach to Machine Learning in Atomistic Simulations. Comput. Phys. Commun. 2016, 207, 310−324. (24) Xie, Y.; Naguib, M.; Mochalin, V. N.; Barsoum, M. W.; Gogotsi, Y.; Yu, X.; Nam, K.-W.; Yang, X.-Q.; Kolesnikov, A. I.; Kent, P. R. Role of Surface Structure on Li-Ion Energy Storage Capacity of TwoDimensional Transition-Metal Carbides. J. Am. Chem. Soc. 2014, 136, 6385−6394. (25) Osti, N. C.; Naguib, M.; Ostadhossein, A.; Xie, Y.; Kent, P. R.; Dyatkin, B.; Rother, G.; Heller, W. T.; van Duin, A. C.; Gogotsi, Y.; et al. Effect of Metal Ion Intercalation on the Structure of MXene and Water Dynamics on Its Internal Surfaces. ACS Appl. Mater. Interfaces 2016, 8, 8859−8863. (26) Wakabayashi, N.; Smith, H. G.; Nicklow, R. M. Lattice Dynamics of Hexagonal Mo S 2 Studied by Neutron Scattering. Phys. Rev. B 1975, 12, 659. (27) Dobardžić, E.; Milošević, I.; Dakić, B.; Damnjanović, M. Raman and Infrared-Active Modes in MS 2 Nanotubes (M= Mo, W). Phys. Rev. B: Condens. Matter Mater. Phys. 2006, 74, 33403. (28) Jiang, J.-W.; Wang, B.-S.; Rabczuk, T. Phonon Modes in SingleWalled Molybdenum Disulphide Nanotubes: Lattice Dynamics Calculation and Molecular Dynamics Simulation. Nanotechnology 2014, 25, 105706. (29) Sandoval, S. J.; Yang, D.; Frindt, R. F.; Irwin, J. C. Raman Study and Lattice Dynamics of Single Molecular Layers of MoS 2. Phys. Rev. B: Condens. Matter Mater. Phys. 1991, 44, 3955. (30) Lehman, G. W.; Wolfram, T.; De Wames, R. E. Axially Symmetric Model for Lattice Dynamics of Metals with Application to Cu, Al, and Zr H 2. Phys. Rev. 1962, 128, 1593. (31) Jiang, J.-W.; Park, H. S.; Rabczuk, T. Molecular Dynamics Simulations of Single-Layer Molybdenum Disulphide (MoS2): 639


Letter

The Journal of Physical Chemistry Letters (51) Yin, H.; Qi, H. J.; Fan, F.; Zhu, T.; Wang, B.; Wei, Y. Griffith Criterion for Brittle Fracture in Graphene. Nano Lett. 2015, 15, 1918− 1924. (52) Pei, Q. X.; Zhang, Y. W.; Shenoy, V. B. A Molecular Dynamics Study of the Mechanical Properties of Hydrogen Functionalized Graphene. Carbon 2010, 48, 898−904. (53) Grantab, R.; Shenoy, V. B.; Ruoff, R. S. Anomalous Strength Characteristics of Tilt Grain Boundaries in Graphene. Science 2010, 330, 946−948. (54) Mahmoudi, N.; Ostadhossein, F.; Simchi, A. Physicochemical and Antibacterial Properties of Chitosan-polyvinylpyrrolidone Films Containing Self-organized Graphene Oxide Nanolayers. J. Appl. Polym. Sci. 2016, 133, n/a. (55) Dang, K. Q.; Spearot, D. E. Effect of Point and Grain Boundary Defects on the Mechanical Behavior of Monolayer MoS2 under Tension via Atomistic Simulations. J. Appl. Phys. 2014, 116, 013508. (56) Dang, K. Q.; Simpson, J. P.; Spearot, D. E. Phase Transformation in Monolayer Molybdenum Disulphide (MoS2) under Tension Predicted by Molecular Dynamics Simulations. Scr. Mater. 2014, 76, 41−44. (57) Bertolazzi, S.; Brivio, J.; Kis, A. Stretching and Breaking of Ultrathin MoS2. ACS Nano 2011, 5, 9703−9709. (58) Li, T. Ideal Strength and Phonon Instability in Single-Layer MoS2. Phys. Rev. B: Condens. Matter Mater. Phys. 2012, 85, 235407. (59) Cooper, R. C.; Lee, C.; Marianetti, C. A.; Wei, X.; Hone, J.; Kysar, J. W. Nonlinear Elastic Behavior of Two-Dimensional Molybdenum Disulfide. Phys. Rev. B: Condens. Matter Mater. Phys. 2013, 87, 35423. (60) Peng, Q.; Liang, C.; Ji, W.; De, S. A Theoretical Analysis of the Effect of the Hydrogenation of Graphene to Graphane on Its Mechanical Properties. Phys. Chem. Chem. Phys. 2013, 15, 2003−2011. (61) Malekmotiei, L.; Samadi-Dooki, A.; Voyiadjis, G. Z. Nanoindentation Study of Yielding and Plasticity of Poly (Methyl Methacrylate). Macromolecules 2015, 48, 5348−5357. (62) Samadi-Dooki, A.; Malekmotiei, L.; Voyiadjis, G. Z. Characterizing Shear Transformation Zones in Polycarbonate Using Nanoindentation. Polymer 2016, 82, 238−245. (63) Voyiadjis, G. Z.; Samadi-Dooki, A. Constitutive Modeling of Large Inelastic Deformation of Amorphous Polymers: Free Volume and Shear Transformation Zone Dynamics. J. Appl. Phys. 2016, 119, 225104. (64) Mortazavi, B.; Ostadhossein, A.; Rabczuk, T.; van Duin, A. C. Mechanical Response of All-MoS2 Single-Layer Heterostructures: A ReaxFF Investigation. Phys. Chem. Chem. Phys. 2016, 18, 23695− 23701. (65) Lin, Y.-C.; Dumcenco, D. O.; Huang, Y.-S.; Suenaga, K. Atomic Mechanism of the Semiconducting-to-Metallic Phase Transition in Single-Layered MoS2. Nat. Nanotechnol. 2014, 9, 391−396. (66) He, H.; Lu, P.; Wu, L.; Zhang, C.; Song, Y.; Guan, P.; Wang, S. Structural Properties and Phase Transition of Na Adsorption on Monolayer MoS2. Nanoscale Res. Lett. 2016, 11, 1−8. (67) Guo, Y.; Sun, D.; Ouyang, B.; Raja, A.; Song, J.; Heinz, T. F.; Brus, L. E. Probing the Dynamics of the Metallic-to-Semiconducting Structural Phase Transformation in MoS2 Crystals. Nano Lett. 2015, 15, 5081−5088. (68) Kushima, A.; Qian, X.; Zhao, P.; Zhang, S.; Li, J. Ripplocations in van Der Waals Layers. Nano Lett. 2015, 15, 1302−1308. (69) Zhou, W.; Zou, X.; Najmaei, S.; Liu, Z.; Shi, Y.; Kong, J.; Lou, J.; Ajayan, P. M.; Yakobson, B. I.; Idrobo, J.-C. Intrinsic Structural Defects in Monolayer Molybdenum Disulfide. Nano Lett. 2013, 13, 2615− 2622. (70) Kittel, C. Introduction to Solid State; John Wiley & Sons, 1966. (71) Le, D.; Rawal, T. B.; Rahman, T. S. Single-Layer MoS2 with Sulfur Vacancies: Structure and Catalytic Application. J. Phys. Chem. C 2014, 118, 5346−5351. (72) Lu, Q.; Arroyo, M.; Huang, R. Elastic Bending Modulus of Monolayer Graphene. J. Phys. D: Appl. Phys. 2009, 42, 102002.

(73) Jiang, J.-W.; Qi, Z.; Park, H. S.; Rabczuk, T. Elastic Bending Modulus of Single-Layer Molybdenum Disulfide (MoS2): Finite Thickness Effect. Nanotechnology 2013, 24, 435705.

640


ReaxFF Reactive Force-Field Study of Molybdenum Disulfide (MoS2

Recommend Documents