# Rippling ultrafast dynamics of suspended 2D monolayers, graphene

^{a}Physical Biology Center for Ultrafast Science and Technology, Arthur Amos Noyes Laboratory of Chemical Physics, California Institute of Technology, Pasadena, CA 91125;^{b}Theory and Simulation of Materials and National Centre for Computational Design and Discovery of Novel Materials, École Polytechnique Fédérale de Lausanne, 1015 Lausanne, Switzerland

See allHide authors and affiliations

Edited by Jacqueline K. Barton, California Institute of Technology, Pasadena, CA, and approved September 6, 2016 (received for review July 6, 2016)

## Significance

Rippling is an intrinsic feature of 2D materials, responsible for their structural stability, transport properties, and electron–hole charge redistribution. Modulating these ripples in a controlled manner not only provides a better understanding of their structural properties, but also has potential impact for applications. Here, we examine graphene monolayer as a prototypical 2D material. An ultrafast attenuation of the ripples intrinsically present in the graphene plane is followed by a significant enhancement of the rippling effect on a longer time scale, as driven by the successive excitation of in-plane and out-of-plane phonon modes. The methodology described is of a general nature and is suitable for the investigation of other 2D materials, where we expect to observe similar rippling effects.

## Abstract

Here, using ultrafast electron crystallography (UEC), we report the observation of rippling dynamics in suspended monolayer graphene, the prototypical and most-studied 2D material. The high scattering cross-section for electron/matter interaction, the atomic-scale spatial resolution, and the ultrafast temporal resolution of UEC represent the key elements that make this technique a unique tool for the dynamic investigation of 2D materials, and nanostructures in general. We find that, at early time after the ultrafast optical excitation, graphene undergoes a lattice expansion on a time scale of 5 ps, which is due to the excitation of short-wavelength in-plane acoustic phonon modes that stretch the graphene plane. On a longer time scale, a slower thermal contraction with a time constant of 50 ps is observed and associated with the excitation of out-of-plane phonon modes, which drive the lattice toward thermal equilibrium with the well-known negative thermal expansion coefficient of graphene. From our results and first-principles lattice dynamics and out-of-equilibrium relaxation calculations, we quantitatively elucidate the deformation dynamics of the graphene unit cell.

- rippling dynamics
- monolayer graphene
- 2D materials
- ultrafast electron diffraction
- first-principles lattice dynamics

For a long time it was believed that 2D materials could not exist in nature due to a crumpling effect induced by long-wavelength thermal fluctuations (1⇓–3). Recently, however, graphene was isolated from graphite by Novoselov et al. using the simple Scotch tape idea (4), and since then more and more 2D materials have been produced from bulk layered crystals (5, 6). In a strict definition, so far no isolated 2D layers have been found to be perfectly flat due to the presence of ripples, a microscopic roughening of the 2D plane, as observed in almost all cases (7, 8), apart from those in perfect contact with a substrate (9). Accordingly, rippling can be considered as an intrinsic feature of 2D materials, and is widely accepted to be at the origin of their structural stability (7, 10, 11).

Graphene, as the first 2D material, has attracted enormous interest due to extraordinary properties originating from the quasi-relativistic character of its electronic band dispersion (12). Most impressively, graphene's excellent transport properties, such as ultrahigh carrier mobility and electrical conductivity, have found a wide range of applications in new generations of nanoelectronic devices (13, 14). Considering that ripples can strongly affect the transport properties of graphene by inducing effective magnetic fields (15⇓–17) and changing local potentials (18, 19), one of the remaining open questions is whether it is possible to modulate the rippling such that active control of the transport properties can be achieved. To address this issue, the main challenge lies in the ability to modulate ripples in a controlled manner at the mesoscopic and microscopic length scales. Bao et al. (20) have introduced a strain-engineering method to create periodic ripples on suspended graphene. Their work undoubtedly provides an improved understanding of suspended 2D layers and graphene devices. However, the lateral size of the generated ripples is one order of magnitude larger than the value generally reported for intrinsic rippling (∼10 nm), which may have a different impact on the transport. Theoretically, Smolyanitsky and Tewary (21) have proposed a dynamic way to tune the ripples by using the excitation of terahertz waves. They predicted that it is possible to significantly flatten graphene exploiting the interplay between intrinsic thermally induced and externally excited ripples.

In this work, we use femtosecond laser pulses to excite the electronic states of a suspended graphene layer and modulate the rippling structure via the excitations and the anharmonic interactions of different phonon modes. To investigate the real-time response of the atomic structure of the suspended graphene to the ultrafast optical excitation, we used ultrafast electron crystallography (UEC), which represents a unique tool for investigating energy transport processes and structural dynamics of nanomaterials (22⇓–24). Compared with other ultrafast optical and X-ray time-resolved techniques, UEC offers not only a higher spatial resolution, down to the atomic scale, but also a high sensitivity to small material volumes, such as monoatomic layers, due to a higher electron/matter scattering cross-section.

Following a femtosecond optical excitation, the temporal behavior of graphene’s diffraction pattern is monitored as a function of the delay time between laser pump and electron probe. An energy evolution process from photocarrier excitation to thermal lattice equilibration is well resolved and understood in the framework of electron–phonon and phonon–phonon couplings. Here we support the experimental observation with first-principles simulations of the phonon dynamics via Boltzmann-transport equations. The excitation of in-plane and out-of-plane phonon modes via anharmonic interactions generates very specific structural dynamics. In particular, a transient nonthermal expansion mediated by in-plane phonons takes place on the 5-ps time scale, inducing an ultrafast attenuation of the ripples present in the graphene layer through the layer’s expansion. This attenuation is followed by a slower thermal equilibrium contraction, caused by out-of-plane phonons with a time constant of 50 ps, which induces a long-term recovery or even enhancement of the rippling structure. Because ripples are ubiquitous in 2D materials, we expect these results to be of general validity and of potential applications in ultrafast photoelectronics.

## Intrinsic Rippling at Steady State

From a theoretical point of view, Landau and Lifshitz (1), Peierls (2), and Mermin (3) showed that in the harmonic approximation, long-wavelength thermal fluctuations are able to destroy the long-range order of a perfect 2D lattice, resulting in melting at finite temperature. Beyond the harmonic approximation, however, the interaction between bending and stretching phonon modes could in principle suppress long-wavelength fluctuations and thus stabilize the 2D lattice through an out-of-plane deformation with strong height modulations (11). In an elegant transmission electron microscopy (TEM) experiment, Meyer et al. (7) revealed that a suspended graphene sheet is not perfectly flat but exhibits pronounced out-of-plane corrugations (ripples) with a lateral size of ∼10 nm and height up to 1 nm. The ripples, whose formation results from the freezing of thermal out-of-plane phonon modes (generally referred to as ZA modes), would lead to a gain in elastic energy and thus a minimization of the total free energy (10).

This 3D warping is able not only to stabilize the 2D lattice, but also to influence the electronic structure and transport properties of graphene. In particular, the bending of the graphene sheet will induce two main effects (19): (*i*) a π–σ rehybridization between nearest neighbors, which shifts the π-orbital energy, resulting in a local potential variation; and (*ii*) a change in the nearest-neighbor hopping integral, which introduces an effective vector potential. As a result, due to the bending, the electrons are subject to a potential depending on the local curvature of the graphene lattice (12) and, therefore, Dirac fermions will be effectively scattered by the ripples, providing an additional contribution to the electrical resistivity. In fact, carrier mobilities higher than 200,000 cm^{2}/Vs have been calculated in the case of a perfect 2D lattice, whereas one order of magnitude smaller values are experimentally reported for monolayer graphene at low temperatures (25). Moreover, a recent work has shown that periodic nanoscale ripples in graphene can induce a significant bond stretching able to open up the band gap by values as large as 0.13 eV (26). Therefore, the ability to control in space and time the lattice distortion in graphene, which allows modulation of the ripples, represents a powerful methodology for tuning graphene's electronic and transport properties.

## Rippling Dynamics Under Optical Excitation

Femtosecond laser pulses are able to perturb the anharmonic coupling between bending and stretching modes, thus modulation of the rippling structure might be expected. To monitor such changes at the atomic scale, we adopt UEC. The UEC experiments are schematically shown in Fig. 1*A* and detailed in *Materials and Methods*. Briefly, a monolayer graphene, grown by chemical vapor deposition (CVD) and characterized by Raman spectroscopy, electron microscopy, and diffraction, is transferred via a poly(methyl methacrylate) (PMMA) layer onto a TEM copper grid. The diffraction pattern in Fig. 1 *C* and *E* and the Raman spectrum in Fig. 1*D* exhibit the typical fingerprint of a single layer of sp^{2} carbon atoms (27). Here, particular care has been taken to remove the PMMA layer, and a very low amount of residue is present on the graphene, as shown by the TEM image in Fig. 1*B*. For comparison, a commercially available monolayer graphene TEM grid (from Ted Pella, Inc.), which shows an increased amount of PMMA residue, has also been investigated. The dynamics are initiated by a 120-fs laser pulse (800 nm, 2 kHz) and probed by a 20-keV electron pulse as a function of the delay time between optical pump and electron probe pulses. The evolution of the diffraction pattern is quantitatively evaluated by measuring the temporal change of the diffraction intensity and peak position of Bragg reflections.

Shown in Fig. 2*A* are the observed transients of the diffraction intensity of the second-order Bragg reflections at two different excitation fluences. Similar results are obtained for the first-order spots, although with a smaller intensity change. Because the six Bragg reflections associated with the same diffraction order show the same behavior, the transients displayed in Fig. 2*A* have been obtained by averaging over these six spots to improve the signal-to-noise ratio. Within the time window of interest, the transients can be well described by a single-exponential decay.

As shown in Fig. 2*B*, the magnitude of the intensity change as a function of excitation fluence (i.e., the laser-induced temperature rise) is quantitatively described by the in-plane Debye–Waller (DW) factor [with in-plane Debye temperature *C*). These values are similar to the decay time of optical phonons (OPs) into thermal acoustic phonons (via nonthermal acoustic phonons), reported for graphite (33⇓–35). The fluence dependence of the time constant is mainly due to an increased phonon–phonon scattering under higher excitation condition. To achieve a simple description of the decay process, we adopt the single relaxation time approximation (36, 37) and define an average relaxation time of the phonon modes *τ* as*T* is the lattice temperature at equilibrium, *k*_{B} is the Boltzmann constant, and *i*) only decays from one high-energy phonon to two low-energy ones are considered, whereas phonon reabsorption is neglected because of the lower occurrence probability; and (*ii*) the Bose–Einstein distribution for the phonon density has been expanded to a first order considering that *C*, Eq. **1** provides a satisfactory fitting of the time constant as a function of the effective lattice temperature, from which a phonon energy

It is worth mentioning that further increase of the excitation fluence will damage the graphene sample. The highest equilibrium temperature reached is 2,400 K, which is close to the stability limit of suspended CVD graphene (38), thus further verifying the validity of the estimated temperature.

As mentioned above, the intensity change in the transmission geometry adopted can mainly probe in-plane atomic vibrations. Nevertheless, out-of-plane deformation effectively modifies the projected atomic position (or unit cell) in the basal plane. It is thus possible to study out-of-plane fluctuations by monitoring the position change of Bragg reflections. To eliminate the effect of the transient electric field (TEF) (39) on the peak position, we extract the Brillouin zone area *A*, as schematically shown in Fig. 3*A* (*Inset*) instead of the absolute peak position. This method is indeed effective because the TEF shifts each Bragg reflection of the same amount along the same direction in the detector plane. The transient behavior of *A* can then be converted into the strain dynamics of the unit cell via *A*, where *A*_{0} is the zone area before optical excitation and *ε* is the in-plane strain. Differently from the case of the diffraction intensity, the strain dynamics show a more unique behavior, composed of an initial ultrafast expansion and a following slower contraction.

The thermal contraction of graphene has been widely studied at equilibrium and is attributed to the large negative thermal expansion coefficient due to thermal fluctuation of out-of-plane (ZA) modes (20, 40). The relatively long-time contraction (time constant of 50 ± 10 ps) observed here can thus be associated with the gradual population of ZA modes. Because the out-of-plane thermal fluctuations are believed to be at the origin of the intrinsic rippling observed in graphene (11), the increased population of ZA modes is responsible for an enhanced corrugation of the graphene plane, as also attested by an increased contraction observed under higher excitation fluence. Here, the temperature range in which the thermal contraction persists is consistent with the theoretical prediction of ref. 41.

Besides the slow contraction, we also observed an ultrafast (5 ps) lattice expansion. This observation demonstrates the transient atomic sensitivity of UEC, in contrast to previous investigations of the lattice deformation of suspended graphene, which have been performed only under static or quasi-static conditions (20, 40). Because the time scale of the lattice expansion is similar to the value observed for the intensity change, which involves in-plane phonon modes, and consistent with the lifetimes of OPs (∼4–5 ps) as calculated by first-principles simulations (42, 43), we attribute the expansion to the excitation of in-plane nonthermal longitudinal and transverse acoustic phonons (LA and TA) anharmonically decayed from OPs. In fact, the nonthermal excitation of LA and TA modes would be able to induce a stretching of the lattice plane, thus reducing the height of the ripples intrinsically present in graphene.

We can thus describe the observed dynamics as a fast decay from OPs to in-plane acoustic phonons, followed by a slower excitation of ZA modes. The origin of these two characteristic relaxation times can be traced back to: (*i*) the scattering selection rule (44), and (*ii*) the conservation laws of energy and momentum. The reflection symmetry z→−z of graphene limits the possible scattering processes such that only events containing an even number of out-of-plane phonons are allowed, excluding several possible channels. In addition, considering that the energy of ZA phonons is much smaller than the energy of in-plane OPs, the constraint on energy and crystal momentum conservation further suppresses the decay from OPs to ZA phonons. These intuitive arguments are corroborated by a detailed theoretical analysis of scattering rates by Singh et al. (45).

It is worth noting that Gao and Huang (32) have theoretically predicted that anharmonic interactions among in-plane phonon modes could cause a linear and positive thermal expansion, as we have observed here (Fig. 3*B*). Nevertheless, a conceptual difference still exists between the present data and their theoretical results. In fact, in the present case the attenuation of the ripples at the early time after the ultrafast excitation is due to the excitation of only in-plane LA and TA modes, which cause the stretching of graphene, while the lattice is still far from equilibrium. Out-of-plane ZA modes take a much longer time to be populated, finally leading the system to the thermal equilibrium. The cross-over from positive to negative expansion is thus mainly due to a consequential excitation of first in-plane and then out-of-plane modes. In contrast, the theoretical simulations from Gao and Huang (32) are performed at thermal equilibrium and would predict a transition from negative to positive expansion with an increased lattice temperature as caused by the competing effects induced by in-plane and out-of-plane modes.

## Graphene Dynamics Simulations

In this section we model from first-principles the lattice dynamics of monolayer graphene following the optical excitation. An ideal, yet computationally unfeasible, approach involves the description of the full excitation chain: photons exciting electrons, which decay into phonons that eventually thermalize. Because the electron relaxation occurs on a time scale much shorter than that of the lattice dynamics (46), we simplify the problem by neglecting this electron dynamical effect, and begin directly exciting a nonequilibrium phonon population and study its relaxation to equilibrium. We thus suppose that the laser instantaneously introduces additional internal energy in the layer with respect to the initial conditions, which means that the temperature of graphene is instantaneously brought from 300 K to a higher value *T*.

The initial condition in the simulations of nonthermal phonon populations is imposed by first considering that electron–phonon couplings are much stronger for in-plane OPs. Therefore, we suppose that the electronic relaxations have excited only these phonons, while in-plane and out-of-plane (both acoustic and optical) modes, which interact weakly with electrons, are left populated according to the Bose–Einstein distribution at room temperature. To excite the in-plane OPs, we assume that they follow a Bose–Einstein distribution at an excited temperature *E* of the system after the excitation: *T*.

We let the nonthermal phonon populations evolve in time according to the linear Boltzmann transport equation (BTE):*t* of the phonon mode *μ* [a shorthand notation for ** q** and phonon branches

*s*],

*V*is a normalization volume, and

*Materials and Methods*for further details). We stress that the solution to the BTE is found by diagonalizing the complete scattering matrix, i.e., beyond the single-mode relaxation time approximation. Therefore, we obtain the time evolution of the phonon population

*Materials and Methods*for details).

The time evolution of the DW factor is reported in Fig. 4*A* for various temperatures, showing a good agreement with the experimental data of Fig. 2*A*, and thus confirming that the modulation of the laser fluence is directly tuning the thermal energy of the layer. Despite the complexity of the microscopic time evolution of phonons, the DW factor relaxes to the equilibrium value approximately as a single exponential with an effective relaxation time *B* agrees well with the experimental trend in Fig. 2*C*, and qualitatively captures the values within a factor of ∼1.5. Therefore, notwithstanding the approximations adopted in the present modeling, we retain the main features of the experimental lattice dynamics and provide an adequate description of the system in the out-of-equilibrium conditions. We attribute the remaining discrepancy mainly to the choice of initial conditions of the phonon populations, which is the only free parameter of the simulation. The inclusion of the electron dynamics and their decay, which nevertheless formidably increases the modeling complexity and goes out of the scope of this article, could provide a parameter-free description of the very first femtoseconds after the laser pulse and eliminate the arbitrariness in the initial phonon populations.

The lattice dynamics simulations also provide an insight into the ultrafast evolution of strain. In particular, the excitation of the phonon gas is able to induce a change of the lattice pressure, causing the crystal to expand or contract. Within the quasi-harmonic approximation (49) (which we extend here to the time-dependent case) the lattice pressure can be obtained as*C* as a function of time for various temperatures. At short times, pressure has larger positive values, indicating the expansion of the graphene crystal. The pressure then decreases, well below the room temperature value, leading to the layer contraction. Although we are neglecting the electron dynamics, the model is able to capture the qualitative features of the ultrafast expansion/contraction as experimentally observed by the strain dynamics in Fig. 3*A*. This behavior can easily be interpreted in terms of Eq. **3** and of **3**, the equilibrium pressure is essentially a sum over Grüneisen parameters, which is dominated by the large contribution to contraction of the ZA mode.

## Summary and Conclusions

The dynamic response of suspended graphene to ultrafast laser pulses, as experimentally observed and theoretically simulated here, is schematically pictured in Fig. 5*A*. At early times (∼5 ps) after the ultrafast excitation, the ripples are stretched out by anharmonic excitation of in-plane phonons, which should reduce their height. On a longer time scale (50 ps), the excitation of ZA phonon modes increases the out-of-plane thermal fluctuation, resulting in an increase of the ripple height. The modulation of the rippling by ultrashort laser pulses not only improves our understanding of structural stability of graphene, but also provides possibilities for new generations of ultrafast photoelectronics. As described before, because ripples are “imperfections” on an otherwise perfect 2D lattice, they have a significant impact on the properties of graphene, including electronic transport, magnetoresistance, and chemical activity. By controlling the rippling structure on an ultrafast time scale, it would thus be possible to significantly modify the physical properties within a very short time window. As an example, we estimated the modulation of the electric resistivity of graphene resulting from the nonthermal in-plane stretch at different initial temperatures based on the theory given in ref. 50, as shown in Fig. 5*B*. At the maximum expansion experimentally reached, the resistivity due to thermal out-of-plane phonons is reduced more than one order of magnitude compared with the value without tension. This reduction is obtained in only few picoseconds, and thus graphene could be used as an ultrafast switch, as conceptually shown in Fig. 5*C*. When an ultrashort laser pulse illuminates the suspended graphene, the electric resistivity decreases at least by one order of magnitude in 5 ps, generating an ultrashort electrical current.

In conclusion, we have visualized the spatiotemporal behavior of suspended monolayer graphene under femtosecond optical excitation using UEC and first-principles simulations of the real-time lattice dynamics. Because of the high electron/matter scattering cross-section, the atomic-scale spatial resolution, and the ultrafast temporal resolution of UEC, we were able to resolve the energy evolution process in graphene from the optical excitation to the lattice thermal equilibrium and reveal the deformation dynamics of graphene unit cell. An ultrafast expansion occurs in 5 ps, which is responsible for an initial reduction of the ripples intrinsically present in graphene, followed by a long-term contraction in 50 ps, which leads to a recovery or even enhancement of the original rippling. This is driven by a sequential excitation and interaction of in-plane and out-of-plane acoustic phonon modes. The results reported here demonstrate the ability to modulate rippling in graphene by using ultrashort laser pulses and, because ripples are ubiquitous in 2D materials, they are expected to be of general validity with potential applications in ultrafast photoelectronics.

## Materials and Methods

### Sample Preparation.

The monolayer graphene sample used here is grown by CVD. After spin-coating a PMMA layer on one face of the graphene plane, the sample was then placed on the surface of a soluble polymer (so-called “Trivial Transfer Graphene,” purchased from ACS Material LLC). After polymer was dissolved in deionized water, the graphene with PMMA was collected with a circular aperture 2000 mesh TEM grid (hole diameter of 6.5 μm, Ted Pella, Inc.). The grids were then dried by annealing at 50 °C for 10 min, followed by a further bake at 100 °C for 10 min, to enhance the adhesion of graphene with the copper grid and flatten out the wrinkles. The grid was then dipped into acetone for 15 min to dissolve the PMMA layer and rinsed several times using ethanol to improve the cleanness of the graphene layer. Finally, the grid with graphene was annealed in a furnace with Ar (∼500 standard cubic centimeters per minute, sccm) and H_{2} (∼500 sccm) at 400 °C for 3 h to minimize the PMMA residue. We also investigated a commercially available monolayer graphene suspended on the 2000 mesh grid (Ted Pella, Inc.), which shows a nonnegligible amount of PMMA impurities resulting in a strong background scattering in the diffraction pattern. The layer number of both samples was identified by the ratio between the intensities I_{first} and I_{second} of the first-order and second-order diffraction spots, respectively. For the case of a monolayer I_{first}/I_{second} is around 1, whereas for a bilayer I_{first}/I_{second} is around 0.5 (see ref. 27).

### UEC Experiments.

The UEC experiments have been performed in the transmission geometry. Briefly, a laser-pump/electron-probe scheme is adopted to detect the photoinduced change as a function of the delay time between the pump and probe pulses. The pump pulse (800 nm, 120 fs), generated from a Ti:sapphire amplifier at the repetition rate of 2 kHz, is used to initiate the dynamics with an incident angle of 45° with respect to the sample surface. The electron probe beam (subpicosecond), emitted by a LaB_{6} photocathode under the illumination of a 266-nm UV laser pulse and then accelerated to 20 keV, is used to produce diffraction patterns off the sample at a certain delay time. Because of the limited lateral size of monolayer graphene flakes, the electron beam is tightly focused on the sample position with a spot size of 31 × 13 μm^{2}. To minimize space-charge effects, each electron probe pulse contains only ∼300 photoelectrons and the upper limit of the electron pulse duration is estimated to increase to 2.1 ps in the present experiment. The spatial overlap between pump and probe beams is ensured with a precision of ∼20 μm by maximizing the transmission of both beams through a 150-μm aperture at the sample plane. The diffraction patterns are recorded on a microchannel plate/phosphor screen/CCD assembly working in the gate mode.

### First-Principles Calculations.

We used Quantum ESPRESSO (48) for density-functional theory and density-functional perturbation theory calculations, using the local-density approximation, norm-conserving pseudopotentials from the PSLibrary (qe-forge.org/gf/project/pslibrary/), a plane-wave cutoff of 90 Ry, and a Methfessel–Paxton smearing of 0.02 Ry. The slab geometry used to simulate graphene has a relaxed lattice parameter of a = 4.607 Bohr and a cell height c = 3a. We integrate the Brillouin zone with a Gamma-centered Monkhorst–Pack mesh of 24 × 24 × 1 points and compute second- and third-order force constants on meshes of 16 × 16 × 1 and 4 × 4 × 1 points, respectively.

### Lattice Dynamics Simulations.

The scattering matrix is built to include 3-phonon and harmonic isotopic scattering (37, 51) at natural carbon abundances (98.93% C^{12}, 1.07% C^{14}), and is constructed using the same computational parameters of refs. 36, 52 (a Gaussian smearing of 10 cm^{−1} and a mesh of 128 × 128 × 1 points for integrating the Brillouin zone). The scattering matrix is diagonalized exactly using the routine PDSYEV of the ScaLAPACK library (53).

The DW factor is computed according to

with

where ** G** is the momentum transfer vector,

*b*is the index over the basis of atoms in the primitive unit cell of the crystal,

*M*

_{b}is the mass of the atom

*b*, and

## Acknowledgments

The authors gratefully acknowledge the assistance of Dr. B. Chen and Dr. X. Fu for TEM characterization and Dr. B. Liao for helpful discussion. This work was supported by the Air Force Office of Scientific Research in the Center for Physical Biology at California Institute of Technology funded by the Gordon and Betty Moore Foundation, the Swiss National Science Foundation through Grant 200021_143636, the Max Planck–Ecole Polytechnique Fédérale de Lausanne (EPFL) Center for Molecular Nanoscience and Technology, and the Swiss National Supercomputing Center (CSCS).

## Footnotes

- ↵
^{1}To whom correspondence should be addressed. Email: jhu{at}caltech.edu. ↵

^{2}Present address: Laboratory for Ultrafast Microscopy and Electron Scattering, Institute of Physics, Lausanne Center for Ultrafast Science, École Polytechnique Fédérale de Lausanne, CH-1015 Lausanne, Switzerland.↵

^{3}Deceased August 2, 2016.

Author contributions: J.H., G.M.V., and A.H.Z. designed research; J.H., G.M.V., A.C., N.M., and A.H.Z. performed research; J.H., G.M.V., A.C., N.M., and A.H.Z. analyzed data; and J.H., G.M.V., A.C., N.M., and A.H.Z. wrote the paper.

The authors declare no conflict of interest.

This article is a PNAS Direct Submission.

## References

- ↵.
- Landau LD,
- Lifshitz EM

- ↵.
- Peierls RE

- ↵
- ↵.
- Novoselov KS, et al.

- ↵.
- Novoselov KS, et al.

- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵.
- Bonaccorso F, et al.

- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵.
- Nair RR, et al.

- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵
- ↵.
- Baroni S,
- Giannozzi P,
- Isaev E

- ↵
- ↵
- ↵
- ↵.
- Blackford LS, et al.

## Citation Manager Formats

## Article Classifications

- Physical Sciences
- Chemistry