## Thickness matters in graphene stacks

If you stack graphene monolayers on top of each other, the number of layers will affect the properties of the material. Intuitively, one would expect that as the stack becomes thicker, the results will converge as the sample starts to resemble graphite. Nam *et al.* measured the conductance of graphene multilayers of increasing thickness. Studying samples up to seven layers thick, they found that in all of them, electronic correlations caused a phase transition at a nonzero critical temperature. However, the critical temperature, as well as the nature of the low-temperature state, depended strongly on the number of layers. This unexpectedly persistent dependence showed no signs of slowing down and will motivate further theoretical and experimental work.

*Science*, this issue p. 324

## Abstract

Suspended Bernal-stacked graphene multilayers up to an unexpectedly large thickness exhibit a broken-symmetry ground state whose origin remains to be understood. We show that a finite-temperature second-order phase transition occurs in multilayers whose critical temperature (*T*_{c}) increases from 12 kelvins (K) in bilayers to 100 K in heptalayers. A comparison of the data with a phenomenological model inspired by a mean-field approach suggests that the transition is associated with the appearance of a self-consistent valley- and spin-dependent staggered potential that changes sign from one layer to the next, appearing at *T*_{c} and increasing upon cooling. The systematic evolution with thickness of several measured quantities imposes constraints on any microscopic theory aiming to analyze the nature of electronic correlations in this system.

Clean two-dimensional (2D) conductors in the presence of a large perpendicular magnetic field are strongly correlated systems. Their ground states are determined by Coulomb repulsion between electrons and are characterized by broken symmetries and nontrivial topological invariants that depend sensitively on electron density (*n*) and applied field (*B*) (*1*–*4*). Accordingly, a series of quantum phase transitions occurs upon varying *n* and *B*, which manifest themselves in the very rich phenomenology ubiquitously observed in magneto-transport measurements (*5*–*13*). These considerations hold true irrespective of the specific 2D conductor considered, because they rely almost exclusively on the formation of Landau levels that quench the electron kinetic energy and allow electron-electron (*e*-*e*) interactions to dominate. In the absence of Landau levels, however, *e*-*e* interactions play a much less prominent role.

Multiple recent experiments indicate that in graphene multilayers this is not the case (*14*–*19*). A gapped insulating state at *B* = 0, first reported in Bernal-stacked bilayers (*14*–*17*), has been recently observed in all even Bernal-stacked multilayers up to octalayer graphene (8LG) (*18*, *19*). By contrast, odd Bernal-stacked multilayers (so far mono- and trilayer have been studied) remain conducting at low temperatures (*20*–*23*). These findings defy common expectations—namely, that the behavior of graphene multilayers should approach that of graphite as thickness is increased. No direct information about the nature of the insulating state in even multilayers could be obtained before now, because the phenomenon is observed only in suspended graphene devices of the highest quality (*14*–*23*), which makes measurements other than transport challenging. Because the resistance of even multilayers increases upon cooling without showing abrupt changes at any specific temperature, it is unclear even whether the insulating state results from a quantum phase transition (with a gap present at all temperatures) or whether a phase transition occurs at a critical temperature (*T*_{c}) with the gap vanishing for temperatures *T* > *T*_{c} (*16*, *17*, *24*–*28*). To address these issues, we studied ultraclean, suspended Bernal-stacked multilayers of graphene near charge neutrality and have shown that at *B* = 0 these systems unambiguously exhibit an interaction-driven finite-temperature phase transition to a broken-symmetry state at a critical temperature *T*_{c} that depends on thickness.

To draw these conclusions, we focused on the density of electrons present in the conduction band of charge-neutral multilayers, *n*_{th}(*T*), whose temperature dependence is determined by the low-energy density of states (DOS). At sufficiently low *T*, *n*_{th}(*T*) should show an exponential increase in the presence of a gap—given that at charge neutrality electrons in the conduction band are thermally activated from the valence band—or stay constant if an overlap between the valence and conduction bands is present. If the multilayers are zero-band gap semiconductors, *n*_{th}(*T*) is expected to increase with increasing *T*, with a specific dependence determined entirely by the low-energy DOS. Quantitative information can be obtained by comparing experimental data with the dependence of *n*_{th}(*T*) calculated from a chosen theoretical expression for the low-energy DOS, as we illustrate for the case of bilayers. Under the assumption that only the nearest-neighbor in-plane (γ_{0}) and out-of-plane (γ_{1}) hopping terms are relevant (the so-called minimal tight-binding model), the low-energy quadratic dispersion relation of electrons in bilayer graphene (Fig. 1A) leads to an energy-independent DOS (Fig. 1B). The resulting density *n*_{th}(*T*) of electrons in the conduction band then increases linearly with temperature (Fig. 1C). However, if a gap Δ opens (Fig. 1D), the DOS is modified (Fig. 1E), and so is the temperature dependence of *n*_{th}(*T*). Figure 1F shows *n*_{th}(*T*) expected for a gap Δ exhibiting a mean-field temperature dependence and vanishing at *T*_{c} [see the inset and (*29*)]. A transition becomes visible at *T* = *T*_{c}, below which *n*_{th}(*T*) is pronouncedly suppressed compared with the noninteracting case.

What makes these considerations useful is that the very small charge inhomogeneity present in suspended graphene multilayers enables the density *n*_{th}(*T*) of thermally excited electrons to be determined experimentally over a broad range of temperatures. *n*_{th}(*T*) can be extracted from the density dependence of the conductance *G*(*n*) near charge neutrality. To understand why and how, it is useful to look at the double logarithmic plot of *G*(*n*) (inset of Fig. 1G). We see that there is a range of *n* over which log(*G*) is constant and that log(*G*) starts increasing substantially only when *n* becomes larger than a threshold [which we denote *n**(*T*); *n**(*T*) increases at higher *T*]. The physical reason is clear: The square conductance remains virtually unchanged as long as the density *n* of electrons accumulated by the gate voltage is much smaller than the density of electrons *n*_{th}(*T*) already present owing to thermal activation from the valence band. Finding that the square conductance increases substantially only when *n* exceeds *n**(*T*) therefore implies that *n*_{th}(*T*) *~ n**(*T*), as discussed previously for mono- and bilayer graphene (*21*, *30*–*32*). In practice, the exact value of *n**(*T*) depends on the criterion used for its definition: Here, we define *n**(*T*) as the value of *n* for which the conductivity σ(*n,T*) = 1.67σ(*n* = 0*,T*), a relation that we can validate for bilayers and that is approximately correct for thicker multilayers (*29*). We emphasize, however, that none of our key conclusions depend on the precise criterion used (*29*).

The temperature dependence of *n**(*T*) extracted from measurements on four different bilayer devices is shown in Fig. 1G. For all devices, a critical temperature *T*_{c} ≅ 12 K is found, below which *n**(*T*) is suppressed compared with the gapless case. The red continuous curve in Fig. 1G, which represents *n*_{th}(*T*) expected for a gap having a mean-field temperature dependence and Δ_{0} = Δ(*T* = 0) ≅ 1.9 meV, reproduces the data very satisfactorily. The presence of a clear transition below which a finite Δ appears, the value of *T*_{c}, and the shape of *n**(*T*) are very robust against all aspects of the data analysis. They can be extracted directly from the data without any assumption about the DOS and do not depend on the criterion used to extract *n**(*T*) (only the quantitative determination of Δ requires the DOS to be specified) [see (*29*)]. Notably, we found the same value for the critical temperature in all the different devices measured, despite the sample-to-sample fluctuations in the absolute value of the low-temperature resistance and in how pronounced the insulating state was at the lowest temperature of our measurements [because of these fluctuations, even the occurrence of an insulating state driven by *e*-*e* interaction has been questioned in some past experimental work (*14*–*17*, *30*, *33*)]. The reason for this high degree of reproducibility is that our measurements effectively probe the DOS averaged over the entire device area, whereas in the insulating state the absolute value of the resistance is strongly affected by any percolating conducting path [e.g., at edges (*34*, *35*)] that occupies a negligible fraction of the total area, giving a negligible contribution to the DOS. The highly reproducible behavior of *n**(*T*) allows us to establish unambiguously the occurrence of an insulating, gapped state near charge neutrality, which is entered through a second-order phase transition at *T*_{c} = 12 K.

The same strategy outlined for bilayers can be applied to thicker even multilayers, in which multiple conduction and valence bands are present (*29*, *36*–*39*). Recent work (specifically, the quantization of the Hall effect in low-magnetic-field and magneto-Raman experiments) has provided evidence that in suspended devices at low energy these bands are well described by including only nearest-neighbor in-plane (γ_{0}) and out-of-plane (γ_{1}) hopping terms (*18*, *19*, *40*, *41*). According to this description, all conduction and valence bands in even multilayers are quadratic at low energy (with different effective masses), so that the DOS in the noninteracting case is again energy independent [see (*29*)]. The same argument used for bilayers implies that *n*_{th}(*T*) increases linearly with temperature in the absence of interactions and that the opening of a gap causes a suppression for *T* < *T*_{c}. Figure 2, A and C, illustrates the results of experiments performed on tetralayer (4LG) and hexalayer (6LG) graphene devices: *n**(*T*) depends linearly on temperature at sufficiently high temperatures, but below a critical temperature *T*_{c} (38 K for 4LG and 90 K for 6LG), *n**(*T*) is suppressed, just as for bilayers (Fig. 1, F and G). For 4LG and 6LG, the data are excellently reproduced by assuming that a gap Δ with a mean-field temperature dependence opens simultaneously at *T*_{c} on all quadratic bands (red line; Δ_{0} ≅ 5.2 meV for 4LG and Δ_{0} ≅ 13 meV for 6LG). As in the bilayer case, the temperature dependence of the resistance at charge neutrality (Fig. 2, B and D) demonstrates the insulating nature of thicker even multilayers but provides no specific feature allowing the identification of a critical temperature.

The analysis of the measured temperature dependence of *n**(*T*) and the comparison with the calculated expression for *n*_{th}(*T*) can also be performed and interpreted in the same way for odd Bernal multilayers. Within the minimal tight-binding model mentioned above, odd multilayers contain multiple conduction and valence bands, all having a quadratic dispersion, except for one that is a linear Dirac band (*29*, *36*–*39*). The DOS associated with the linear band vanishes at charge neutrality and contributes at most a few percent of the total DOS in the temperature and density ranges relevant for our work (*29*). It is therefore negligible in practice, so that the same considerations made for even multilayers hold true for odd ones. In Fig. 3, data measured on odd Bernal-stacked multilayers confirm that this is the case for both trilayer (3LG) (Fig. 3A) and heptalayer (7LG) graphene (Fig. 3C): *n**(*T*) increases linearly with temperature at high temperatures, exhibiting a pronounced suppression for *T* < *T*_{c} (with *T*_{c} = 33 K for 3LG and 100 K for 7LG), in complete analogy to the case of even multilayers. Like the data for even multilayers, the data for odd multilayers agree quantitatively with the behavior expected if a gap Δ opens simultaneously at *T*_{c} in all quadratic bands, with a mean-field temperature dependence and *T* = 0 values Δ_{0} ≅ 5 meV for 3LG and Δ_{0} ≅ 13 meV for 7LG. For odd multilayers, the low-temperature conductivity remains finite, *σ* ~ *e*^{2}/*h* (where *h* is Planck’s constant), indicating that the linear band does not gap out (Fig. 3, B and D).

In summary, we find that a phase transition occurs in all the Bernal-stacked multilayers investigated, irrespective of whether they are even or odd, with only the even ones becoming insulating at low temperatures. Data for multilayers of all thicknesses are shown in Fig. 4A, each with the corresponding fit to the calculated temperature dependence of the density of thermally activated electrons *n*_{th}(*T*). Notably, all data collapse on top of one another if *n**(*T*)/*n**(*T*_{c}) is plotted versus *T*/*T*_{c}. (Fig. 4A, inset). *T*_{c} increases linearly with increasing thickness (Fig. 4B), and so does the gap Δ_{0} that is proportional to *T*_{c} (Fig. 4C): The best linear fit (the continuous line in Fig. 4C) is close to the expected mean-field value, Δ_{0}/*k*_{B}*T*_{c} = 1.76 (where *k*_{B} is Boltzmann’s constant) (the dashed line in Fig. 4C) (*29*).

The experimental findings consistently provide substantial evidence about the microscopic nature of the broken-symmetry state, namely, that the effect of interactions is well described by a valley- and spin-dependent mean-field staggered potential that changes sign from one layer to the next, acting as an order parameter (*18*, *19*). A scenario based on the minimal tight-binding model augmented with a self-consistent staggered potential explains that at low temperatures, even multilayers become fully insulating, whereas in odd multilayers a finite conductivity of ~*e*^{2}/*h* persists (see Figs. 2, B and D, and 3, B and D); that from 1LG to 8LG, the first quantum Hall plateau appears systematically at filling factor ν = 2*N*, with σ_{xy} = 2*Ne*^{2}/*h*, where *N* is the multilayer thickness [see (*18*) and (*19*)]; that the transition occurs in odd multilayers; why the gap in the different quadratic bands has the same magnitude; and why the gap opens simultaneously in all bands at the same critical temperature (inasmuch as the latter two points are concerned, a generic mechanism could gap each band independently from the others, resulting in multiple transitions with different values of *T*_{c} and Δ) (*42*, *43*). In the scenario that we propose, with the staggered potential acting as an order parameter, the broken symmetry is discrete, which is why the occurrence of finite-temperature phase transitions in Bernal multilayers does not conflict with the Mermin-Wagner theorem (*44*), despite all multilayers being effectively 2D.

Despite the success of the proposed model, in the absence of a comprehensive theoretical analysis, important questions remain. A key question is the validity of the minimal tight-binding model because the vast existing literature on graphite suggests that more tight-binding parameters should be included. A recent theoretical analysis of *e*-*e* interaction at the Hartree-Fock level (*43*) shows that adding hopping terms present in graphite—most notably, the parameter γ_{2}, responsible for the semimetallic behavior of graphite (γ_{2} ~ −20 meV, usually)—would make the systematic behavior observed in multilayers impossible to reproduce theoretically. At the simplest level, the same holds true for other hopping parameters, such as γ_{3}, responsible for trigonal warping. Why some hopping terms that must be included to describe graphite do not appear in thin multilayers requires an explanation. A key difference may be that our experiments probe charge-neutral multilayers (i.e., *n* ~ 0), whereas in graphite approximately *n* ~ 3 × 10^{11} to 4 × 10^{11} electrons/cm^{2} are present in each individual layer (*45*). It has been established theoretically and experimentally that, as *n* approaches charge neutrality, large renormalization effects drastically change the hopping parameters in graphene: In monolayers, for instance, the Fermi velocity (and hence, γ_{0}) is predicted to diverge as *n* → 0 (*22*, *46*–*48*). Although a thorough analysis of this renormalization process for all hopping parameters is lacking, what is known at this stage implies that there is no compelling reason why the hopping terms used to describe graphite and charge-neutral thin multilayers should be the same.

Irrespective of these details, what is most notable in the experimental findings is the occurrence of a phase transition with increasingly large *T _{c}* in graphene multilayers that, until now, may have been expected to behave as bulk graphite. Why precisely

*T*

_{c}increases with increasing thickness or at which thickness this trend breaks down and the behavior of bulk graphite is recovered remains to be understood. At this stage, a thorough microscopic theoretical analysis is required not only to understand in detail the properties of this family of electronic systems but also to disclose its full potential to reveal subtle phenomena emerging from the physics of correlated electrons.

## Supplementary Materials

www.sciencemag.org/content/362/6412/324/suppl/DC1

Materials and Methods

Supplementary Text

Figs. S1 to S11

This is an article distributed under the terms of the Science Journals Default License.

## References and Notes

**Acknowledgments:**We gratefully acknowledge A. Ferreira for continued technical support of the experiments and T. Giamarchi for useful discussions.

**Funding:**Financial support from the Swiss National Science Foundation, the NCCR QSIT, and the EU Graphene Flagship Project is also gratefully acknowledged.

**Author contributions:**Y.N., D.-K.K., and D.S.-D. fabricated the devices and performed the electrical measurements. Y.N. and D.-K.K. analyzed the data under the supervision of A.F.M. All authors discussed experimental results and their interpretation. A.F.M. wrote the manuscript with input from all authors.

**Competing interests:**There are no competing financial interests.

**Data and materials availability:**The dataset used in this paper is available at Harvard Dataverse (

*50*).