Dependency of Sliding Friction for Two Dimensional Systems on Electronegativity

We study the role of electronegativity in sliding friction for five different two dimensional (2D) monolayer systems using density functional theory (DFT) with van der Waals (vdW) corrections. We show that the friction between the commensurate 2D layered systems depends strongly on the electronegativity difference of the involved atoms. All the 2D layered structures exhibit almost the same magnitude of friction force when sliding along the nonpolar path, independent of the material and the surface structures. In contrast, for sliding friction along the polar path, the friction force obeys a universal linear scaling law as a function of the electronegativity difference of its constituent atoms. Further analyses demonstrate that atomic dipoles in the 2D monolayers induced by the electronegativity difference enhance the corrugation of charge distribution and increase the sliding barrier accordingly. Our studies reveal that electronegativity plays an important role in friction of low dimensional systems, and will provide a strategy for designing nanoscale devices further.


Introduction
Understanding the origin of atom-scale friction processes and controlling them for designing nanotechnological devices pose a major challenge to physicists and nanotribologists alike due to the complex energy dissipation mechanisms and intricate interfacial interaction [1][2][3].
Two-dimensional (2D) layered materials, such as graphene, hexagonal boron nitride (h-BN), and hexagonal molybdenum disulfide (MoS 2 ), due to their strong intralayer chemical bonding compared with the weak interlayer physical adsorption interaction, can serve as solid lubricants to minimize friction and wear in high local pressures and boundary contact regime in a number of applications [4][5][6]. Moreover, 2D layered structures often have novel electronic properties widely researched in literature, which can further our understanding of their tribological behavior [7][8][9][10][11].
Friction is ultimately governed by the atomistic interactions dominated by quantum mechanics [2,3,12]. Several studies have demonstrated that electronic structure, charge distribution, and even spin degree of freedom can influence the friction behaviors of low dimensional system [6,[11][12][13].
Electronegativity plays a significant role in the electronic structure, and electronegativity difference among constituent atoms has an important influence on the physical and chemical properties of materials [14]. For example, although h-BN and graphene are isoelectronic and isostructural, the former is an insulator with a wide band gap of around 6 eV, while the latter is a zero band gap semimetal [15]. This is because the large electronegativity difference between B and N atoms displaces a shared pair of electrons towards the N atom. Thus, electronegativity has an influence on the structure and hence also on the friction. The relationship between electronegativity and friction in low dimensional system was studied in Ref. [16][17][18][19][20][21]. Experiments have shown that the sliding friction in insulating multiwall BN nanotubes (BNNTs) is orders of magnitude stronger than that of semimetallic C nanotubes (CNTs), which were attributed to increased potential barrier caused by the charge localization induced by electronegativity difference between the B and N atoms in BNNTs [16]. This localization effect increases the corrugation amplitude of the interfacial potential.
Ab initio molecular dynamics calculation showed that the coefficient of friction (COF) of liquid water sliding on h-BN is about three times larger than that of on graphene, which was ascribed to the greater corrugation of the energy landscape of h-BN arising from specific electronic structure effect [17]. From DFT calculations it was found that the vdW interaction determines the interlayer binding and the electrostatic interaction mainly influences the sliding barrier of bilayer h-BN [18][19][20]. They speculated that a highly anisotropic interfacial friction should exist for the h-BN bilayer [19]. DFT calculations found that constraints on atomic motion can be employed to tune the contribution of electrostatic interactions and dispersive forces to the sliding energy profile, ultimately leading to different sliding pathways in bilayers of graphene and h-BN [20]. All the above results confirm that the polarity plays an important role in friction of low dimensional systems. However, no detailed investigations of the connection between polarity and friction properties have been published.
In this research work, we systematically investigate the role of electronegativity on friction in 2D layered systems. We show that when two monolayers slide relative to each another, the atom electronegativity difference along the sliding path strongly influence the interfacial friction properties. We show that the friction scales linearly with electronegativity difference in the polar sliding path.

Methodology
Vienna Ab-initio Simulation Package (VASP) code based on the projector augmented-wave (PAW) method was employed in the calculations [22][23][24]. The exchange-correlation interactions were treated with the generalized gradient approximation (GGA) of Perdew, Burke, and Ernzerhof (PBE) [25], with a vdW correction determined by the many-body dispersion (MBD) method [26]. An energy cut-off of 600 eV and 21×21×1 Monkhorst-Pack grids were selected for 2D irreducible Brillouin-Zone integration [27]. The convergence thresholds for total energy and Hellmann-Feynman forces are 10−5 eV and 0.01 eV/Å, respectively. A vacuum space of at least 20 Å was set to avoid the intercell interactions.
Based on the Prandtl-Tomlinson model and its extensions, several researchers have developed methods to compute friction using DFT calculations. In these methods, typical friction parameters, such as shear force, COF and potential energy surface (PES), can be obtained by calculating the energy barrier or sliding energy corrugation along the sliding path [6,11,12,18,28,29]. It should be noted that these methods only evaluate the maximum energy barrier along the sliding path, but did not consider the energy dissipation while sliding. Therefore, the friction calculated by these methods is static friction or break-loose force. In our work, the COF and PES are separately calculated by the above methods [28,29]. However, if it is assumed that all the energy needed to reach the top of the energy barrier is dissipate before climbing the next barrier, and then present calculation also gives information about the kinetic friction [30].

Results and Discussion
In order to study the influence of polarity on friction, we first discuss the polarity character in the sliding path. Fig .1(a) shows a schematic diagram of the sliding interface between two monolayers.
The atoms A, B form the upper layers, and C, D form the bottom layers arrange alternatively along the sliding direction. The total polarity of the sliding system along the sliding direction is then the sum of the polarity of each single layer in that direction, which is decided by the electronegativity differences between A and B atoms, and C and D atoms. Here, bilayer h-BN is chosen as an example to illustrate the polarity differences in different paths, as shown in Fig. 1(b). The most stable AA' stacking eclipsed with B over N atoms was chosen as an initial structure [18,19,31], and two highly symmetric directions were chosen as sliding paths. From the cross-section view, we can see that along path I, there only one kind of B (N) atom in the upper (bottom) layer, and so the polarity is zero along this direction. However, B and N atoms alternately arrange in upper and bottom layers along path II, which give rise to dipoles in both the upper and bottom layers. Therefore, we denote the path I and II as nonpolar and polar paths, respectively. respectively. Based on the lattice parameters, the stable bilayer stacking models for the three systems were obtained (see Table S1 of supporting materials), which is supported by other study [18,19,32]. Two sheets of the 2D monolayer were moved relative to each other along the path I and II to simulate the sliding process. The COF as a function of the load is shown in Fig. 2. For the identical pristine layers of graphene/graphene (Fig. 2(a)), h-BN/h-BN (Fig. 2(b)), and MoS 2 /MoS 2 (Fig. 2(c)), the COF fall in a range of 0.07-0.2, which is in agreement with earlier studies [5,11]. Among these three systems, the nonpolar graphene/graphene system exhibits isotropic friction behavior. However, in the polar h-BN and MoS 2 systems, the COF exhibit anisotropic behavior, in which the COF along polar path II is almost two times larger than for the nonpolar path. These results clearly show that the polarity plays an important role in the friction between two polar planes.
We next investigated the friction between polar and nonpolar sheets by using the isoelectronic interface of graphene on h-BN (graphene/h-BN). The computational model was simplified by enlarging the lattice constant of graphene to be equal to that of h-BN so that the influence of incommensuration and Moiré patterns on friction can be cancelled. It should be noticed that although the upper layer is nonpolar, path II is still polar due to dipoles that exist in the bottom layer along the sliding path ( Fig. 2(d)). Similar to the polar/polar interface, the COF in the nonpolar/polar system is also anisotropic, with a larger COF in polar path. These results indicate that even if dipoles exist only in one layer of the sliding interface, the polarity still has a great influence on the interfacial friction.  Fig. 2(e). It is apparent that the COF along polar path II is larger than that of the nonpolar path I. More importantly, the H-graphene/h-BN system has a larger COF than that of the graphene/h-BN, which is consistent with the increase of the polarity induced by H passivation.
The calculations provide new insight for understanding and tuning friction by surface modification of monolayers.
The variation of the interaction energy as a function of the relative lateral position of the two surfaces in contact can be represented by PES, which is relevant for the friction in the case when no external load is applied [29]. At the fundamental level, the corrugation of PES determines the intrinsic resistance to sliding and gives an indication of the maximum energy that can be dissipated during frictional processes. From the Fig. 3(a) and (b) it is observed that the h-BN system has larger PES corrugation than that of graphene system. To examine the effect of sliding direction on friction, we plotted the potential profiles along the two symmetric directions, indicated in the middle of the Fig. 3. In the graphene/graphene system, the two paths exhibit almost the same potential corrugation. However, for h-BN/h-BN system, the potential corrugation along path II is about two times larger than that of along path I, which result in different friction behaviors as observed earlier in Fig. 2. The lateral forces acting on the slider dragged along the two symmetric directions were also calculated and are shown at the bottom of Fig. 3. For the graphene system, the maximum friction force values are almost equal for the two different paths.
On the contrary, the h-BN system exhibits larger friction force along path II than that of path I.
These results further demonstrate the anisotropic friction behavior of the h-BN system, which corroborate results of COF in Fig. 2. From the above calculations, we conclude that the frictional properties exhibited by 2D systems both with and without external load, are dependent on the polarity and hence the directional sliding. Note that, the sliding potential barrier applied, the friction is decided by the E ∆ . However, for the normal loads used in our study, the major contribution to the sliding barrier comes from the work done against the external force, and the change in E ∆ accounts for only a small part (see fig. S1).
To understand the origin of the anisotropic friction behavior of the h-BN polar system, we calculated the charge density difference of h-BN system by the formula ρ and N ρ are the charge densities of h-BN, free B and N atoms, respectively. Fig. 4(a) shows that charge accumulates around N atoms from neighbor B atoms. The Mulliken population analysis estimates the transferred charge to be 0.84 e (2.18 e for Bader analysis). We further calculated the interlayer charge density difference of the stable bilayer h-BN. From the Fig. 4(b) we can see that when two sheets of h-BN approach each other, the intralayer charge transfer from B to N atoms is further enhanced. However, this effect is noticeably different in the two paths as seen in cross-sectional views of the charge density difference. Along the nonpolar path I, the net positive charge on the upper layer and the negative charge in the bottom layer separately distribute along the sliding direction, which corresponds to a small charge corrugation and smooth sliding barrier. However, along the polar path II (Fig. 4(d) [18,19,33]. The other character is that the polar path exhibits larger friction than that of nonpolar path for all systems.
To better understand the different COF along the polar direction II, we further extracted the functional relation between charge transfers and friction, as shown in Fig. 5 (b). It should be pointed out here that the charge transfer is the sum of the intralayer charge transfer in the upper and lower surfaces. As the different charge analysis method may yields some deviation, the charges transfers were separately calculated by both Mulliken population and Bader charge analysis. Although the two methods provide different amount of charge transfers, the change of the amount of charge transfer with respect to the friction force is quite similar, with errors falling within a small linear window. It is interesting that a linear functional relationship has been found between averaged frictional forces and the amount of intralayer charges transfer. It should be note that, although the COF in polar direction is larger for MoS 2 , the shear force per unit cell is small.
This because the unit of the normal load is nN/unit-cell and the MoS 2 has large cell. Compared to the nonpolar path I, the added friction in the polar path II can now be safely attributed to the enhanced charge density corrugation that appear in the polar path II. The relationship between the electronegative difference and friction is further discussed as follows. To establish a functional relation, we first quantize the electronegativity difference of a 2D interface along a sliding direction. We define the total electronegativity difference of the system along a path as