Giant and Tunable Anisotropy of Nanoscale Friction in Graphene

The nanoscale friction between an atomic force microscopy tip and graphene is investigated using friction force microscopy (FFM). During the tip movement, friction forces are observed to increase and then saturate in a highly anisotropic manner. As a result, the friction forces in graphene are highly dependent on the scanning direction: under some conditions, the energy dissipated along the armchair direction can be 80% higher than along the zigzag direction. In comparison, for highly-oriented pyrolitic graphite (HOPG), the friction anisotropy between armchair and zigzag directions is only 15%. This giant friction anisotropy in graphene results from anisotropies in the amplitudes of flexural deformations of the graphene sheet driven by the tip movement, not present in HOPG. The effect can be seen as a novel manifestation of the classical phenomenon of Euler buckling at the nanoscale, which provides the non-linear ingredients that amplify friction anisotropy. Simulations based on a novel version of the 2D Tomlinson model (modified to include the effects of flexural deformations), as well as fully atomistic molecular dynamics simulations and first-principles density-functional theory (DFT) calculations, are able to reproduce and explain the experimental observations.


Fourier Transform spectra for additional crystallographic orientations
Friction force microscopy images were acquired for several scanning directions in order to analyze the periodicity of stick and slip events as a function of the scanning angle. The Fourier Transform analysis was performed at each scan line of the friction images and the average spectrum for each direction is shown in Fig. S1. Starting from the zigzag (on top), each spectrum was obtained from scans rotated by an additional ~ 5 o from the previous ones. As the AFM tip scans along the zigzag crystallographic direction, a single peak (at 4.0 nm -1 ) is observed in the FT spectrum. As the scanning direction is varied, two prominent peaks are observed in the FT spectra. The separation between these two peaks increases as the scanning direction is deviated from the zigzag.

Modeling the angle dependence of the Fourier spectra of lateral force profiles
The angle dependence of the main peaks observed in both experimental and simulated Fourier spectra can be understood by a simple model. Let us consider first the simpler cases of zigzag (θ = 0 o ) and armchair (θ = 30 o ) directions. For the zigzag direction, the tip develops a nearly straight stick-slip motion with a spatial period equal to graphene lattice constant a (Fig. S2a). The force profile resembles a sawtooth wave, so the corresponding Fourier spectrum will be composed of main peaks at (and smaller peaks at higher harmonics). However, in general the two jumps in a given zigzag will not be symmetric with respect to the scanning trajectory (dashed lines in Fig. S5e), so the corresponding force profile will have a larger period modulation of a 3 1 = λ . This gives rise to peaks at ( ) in the corresponding Fourier spectrum. This is consistent with the fact that the full period of the tip movement along the armchair direction is a In other to analyze the tip motion for a general scan angle θ, it is instructive to consider directions that give rise to commensurate (periodic) motions in the underlying graphene lattice. It is also instructive to consider the near-zigzag (small θ) limit. In this case, it is interesting to make a correspondence between the spatial period of the scanning movement and the chiral vector in a (n,1) single-wall carbon nanotube S1 . In this particular situation, as illustrated in Fig. S2b for the case n = 4, the tip will likely move in a sequence of n jumps along the 1 a  primitive vector of the graphene lattice,   Friction in the HOPG was found to be significantly smaller than for monolayer graphene. Since the interaction between the sheets in HOPG is stronger than the interaction between a monolayer graphene and the SiO 2 substrate, when the tip makes contact with the graphite surface, no puckering around the tip occurs and no wrinkle is formed at the front face of the moving tip while scanning, resulting in non-tilted force friction loops.
Energy dissipation in HOPG was obtained by integrating the friction force over the forward and backwards scans. Figure S4 shows the statistical analysis of the dissipated energy of graphite scanned along the zigzag (a) and armchair (c) directions at different applied normal forces. The energy dissipated along the armchair direction is ∼ 15% higher than along the zigzag crystallographic direction, for different values of applied normal force.

Tomlinson model simulations for HOPG
For the Tomlinson model S2,S3 used to simulate AFM scans in HOPG, described in Equations (1) and (2) Fig. S6 is a simple model based on the geometry of the stick-slip tip dynamics at the HOPG surface. It is assumed that the dissipated energy is proportional to the number of stick-slip jumps, i.e., every jump contributes equally to the dissipated energy. It is further assumed that the tip apex jumps from an hexagon center to the next. Using these assumptions, it is clear from Figs. S5d and S5e that the ratio between armchair and zigzag energy dissipation must be 15 . 1 3 2 ≈ , amazingly close to the experimentally measured 15% friction enhancement for the armchair direction with respect to zigzag, As it turns out, using the same geometrical argument, we can determine the relative energy dissipation for any scan direction with respect to zigzag. In fact, the geometry is quite analogous to that of single-wall carbon nanotubes (SWNTs) S5

Density functional theory calculations of graphene buckling
In order to analyze the Euler buckling of graphene sheets, we perform firstprinciples calculations based on density functional theory (DFT) and pseudopotentials.
The calculations were performed using the Quantum Espresso code S6 . Wavefunctions were expanded on a plane-wave basis with an energy cutoff of 50 Ry. For the electronelectron interactions, we use a PBE-GGA exchange-correlation functional S7 and, for electron-ion interactions, we have use Vanderbilt ultrasoft pseudopotentials S8 . Within this framework, the calculated lattice constant of the graphene sheet at equilibrium is 2.47 Å. In the out-of-plane direction, we have chosen a vacuum distance of 12 Å. In order to simulate the strain effects, we have built rectangular supercells in which we rescale all coordinates and lattice vectors according to the corresponding strain tensor.
We then compare the total energies for two different geometries: one in which the sheet is kept planar throughout the calculation and another in which we start with an initial sinusoidal deformation in the direction of applied stress (Fig. S7a). In both cases, atomic coordinates are fully relaxed. The wavelength of the sinusoidal deformation is chosen to match the supercell length in the corresponding direction; hence we explore different supercell sizes in order to study deformations with different wavelengths. We   Notice that the actual experimental situation differs considerably than that of an isolated and infinite graphene sheet under uniform compressive stress, basically for two reasons: (1) In the case of an AFM scan, the stress is not uniform but it is applied locally at the tip-surface contact point; (2) The sheet is not free, but interacts with a substrate through van der Waals forces. In particular, the latter condition implies that, even in the idealized situation of an infinite sheet, a finite critical stress is needed to buckle a graphene sheet.

Molecular dynamics simulations
In our molecular dynamics simulations, a Tersoff potential S10 was used to model carbon-carbon interactions and a van der Waals force, as implemented in the LAMMPS software S11 , was used to describe graphene-tip and graphene-substrate interactions.
Integration of dynamical equations was carried out considering a constant temperature controlled by a Nose-Hoover thermostat. During all calculations, a 0.1 fs timestep was adopted. Temperatures in a range going from 20K up to 350K were considered and qualitatively equivalent results were obtained for all temperatures.
Movies of the simulations are available online.