# Fileset

[2024A00298G_Craig_TTG_stacking.pdf](https://mdr.nims.go.jp/filesets/9a6a803f-2ac9-440d-95d9-26f9410fa735/download)

## Creator

Isaac M. Craig, Madeline Van Winkle, Catherine Groschner, Kaidi Zhang, Nikita Dowlatshahi, Ziyan Zhu, [Takashi Taniguchi](https://orcid.org/0000-0002-1467-3105), [Kenji Watanabe](https://orcid.org/0000-0003-3701-8119), Sinéad M. Griffin, D. Kwabena Bediako

## Rights

This version of the article has been accepted for publication, after peer review (when applicable) and is subject to Springer Nature’s <a href="https://www.springernature.com/gp/open-science/policies/accepted-manuscript-terms">AM terms of use</a>, but is not the Version of Record and does not reflect post-acceptance improvements, or any corrections. The Version of Record is available online at: http://dx.doi.org/10.1038/s41563-023-01783-y[In Copyright](http://rightsstatements.org/vocab/InC/1.0/)

## Other metadata

[Local atomic stacking and symmetry in twisted graphene trilayers](https://mdr.nims.go.jp/datasets/e60d8386-012b-4b85-aaf4-38d2bbd9f300)

## Fulltext

Local atomic stacking and symmetry in twistedgraphene trilayersIsaac M. Craig1, Madeline Van Winkle1, Catherine Groschner1, Kaidi Zhang1, NikitaDowlatshahi1, Ziyan Zhu2, Takashi Taniguchi3, Kenji Watanabe4, Sinéad M. Griffin5,6, andD. Kwabena Bediako1,7,*1Department of Chemistry, University of California, Berkeley, CA 94720, USA2SLAC National Accelerator Laboratory, Stanford, CA, USA3International Center for Materials Nanoarchitectonics, National Institute for Materials Science, 1-1 Namiki,Tsukuba 305-0044, Japan4Research for Functional Materials, National Institute for Materials Science, 1-1 Namiki, Tsukuba 305-0044, Japan5Molecular Foundry, Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA6Materials Sciences Division, Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA7Chemical Sciences Division, Lawrence Berkeley National Laboratory, Berkeley, CA 94720, USA*Correspondence to: bediako@berkeley.eduAbstractMoiré superlattices formed by twisting trilayers of graphene are a useful model for study-ing correlated electron behavior and offer several advantages over their formative bilayeranalogues, including a more diverse collection of correlated phases and more robust super-conductivity. Spontaneous structural relaxation significantly alters the behavior of moirésuperlattices, and has been suggested to play an important role in the relative stabilityof superconductivity in trilayers. Here, we use an interferometric four-dimensional scan-ning transmission electron microscopy approach to directly probe the local graphene layeralignment over a wide range of trilayer graphene structures. Our results inform a thoroughunderstanding of how reconstruction modulates the local lattice symmetries crucial for es-tablishing correlated phases in twisted graphene trilayers, evincing a relaxed structure thatis markedly different from that proposed previously.1Since their initial discovery, graphene-based moiré superlattices have emerged as valuabletools for investigating the balance between key parameters in strongly correlated phases[1–5]. Their effectiveness stems, in part, from the diverse phenomena they manifest whenadjusting twist angle, carrier density, and thickness. Notably, trilayer graphene showcasesmarkedly distinct properties across its native stacking arrangements [6–9]. Bernal-stacked(ABA) trilayer graphene is a semi-metal and features poorly coupled bands resembling acomposite of monolayer and bilayer graphene [10]. Conversely, rhombohedral-stacked (ABC)trilayer graphene displays hybridization among all three layers, Mott insulating states [11],metallic behavior [12], and superconductivity [13]. The variations in local lattice symme-try significantly contribute to the realization of distinct properties in few-layer graphenesystems, both within high-symmetry structures and within the locally ordered domains ofmoiré superlattices more broadly [14, 15].Particularly noteworthy is the characterization of flat bands in twisted bilayer graphene(TBLG) as a fragile topological phase protected by space-time inversion symmetry [16, 17].This protection persists despite the atomic lattice’s inversion symmetry being limited toinstances where the carbon atoms in each layer align vertically (AA stacking) [18]. The extentof overlap between localized states within the AA regions that are associated with the flatbands, and consequently the size of AA stacking regions, is believed to significantly influencethe superconducting current and transition temperature in twisted structures [19]. Twistedgraphene multi-layers, including magic angle twisted trilayer graphene (MATTLG), rely onthe same C2zT symmetry as their bilayer counterparts [20]. Similar to twisted bilayers,twisted trilayer graphene (TTLG) structures exhibit local C2zT symmetry in pockets ofAAA-type alignment. Recent findings indicating that superconducting phases in MATTLGare more resilient to variations in twist angle and gating than those observed in bilayers[21–24] have prompted intrigue regarding the potentially contrasting impact spontaneouslattice relaxation may play in these two systems. While TBLG spontaneously undergoes areduction in the proportion of AA stacking [25, 26], it has been suggested that relaxation intrilayers may instead augment the prevalence of inversion symmetric alignments [27]. This2underscores the importance of precise structural characterization to uncover the intricaciesof spontaneous lattice relaxation in twisted graphene trilayers.Here, we employ an interferometric methodology based on four-dimensional scanning trans-mission electron microscopy (4D-STEM [28]) known as Bragg interferometry (Fig. 1A) [26,29, 30]. This technique leverages the local interference pattern in diffracted electron beams tounambiguously deduce the stacking orientation of atomic layers. Unlike scanning tunnelingmicroscopy (STM) and more conventional STEM methods, this technique allows us to probemoiré patterns within encapsulated materials. Importantly, it also facilitates the selectiveimaging of individual bilayer interfaces within complex multilayered materials. These re-sults offer a direct assessment of the local atomic stacking within twisted graphene trilayers.We find that the results of this 4D-STEM measurement suggest a picture of reconstructionthat is markedly different from that proposed by previous STM work [27], and one that isconsequential for understanding the correlated electron physics in these materials.The interferometric 4D-STEM technique we use involves scanning a converged electron beamacross the sample of interest and collecting individual diffraction patterns for each real spaceposition of the probe (Fig. 1A). Throughout this work, we use the notation shown in Fig.1B to label the twist angles, θ, within the trilayer sample. Here, θ12 and θ23 denote thetwist angles between layers 1 & 2, and layers 2 & 3 respectively such that θ13 = θ12 + θ23.We further use the labels shown in Fig. 1C to denote the various high-symmetry stackingconfigurations realized within the moiré. The converged beam electron diffraction (CBED)patterns collected at each probe location then appear as shown in Figs. 1D,E, where eachlayer of the material generates a set of Bragg disks. The overlap between Bragg disksoriginating from different layers (inset Figs. 1D,E) is then used to determine the stackingorientation of those two layers. As an example, Figs. 1F,H shows how the intensity of theoverlap between layers 1 and 3—denoted as ‘(1,3) overlap’—varies across the sample. Thismodulation in intensity directly manifests the moiré pattern between layers 1 and 3 but isinsensitive to the orientation of the second layer. Similarly, the variation in the intensity ofthe (1,2,3) overlap region (Figs. 1G,I) reveals the modulation in stacking order between all3three graphene layers.Therefore, by exploiting the relationship between stacking order and overlap region intensity(see Supplemental Information sections 4–6), we map the variation in atomic stacking andhence reconstruction within trilayer graphene samples. Results of these analyses are shownin Fig. 2 for a structure that we call ‘AtA’ in which the top and bottom graphene layersare perfectly aligned to each other and twisted with respect to the middle layer (Fig. 2B).In this structure, the average intensity of the overlap regions in the first ring of Braggreflections and the average overlap region intensity in the second ring of reflections can beused to determine the local stacking configuration. Using the bi-variate color-legend shownin Fig. 2A (in which the high symmetry stacking configurations associated with each colorare overlaid), we create a map of local atomic stacking within an AtA sample (Fig. 2C).The local atomic stacking shown in Fig. 2C indicates that this particular AtA samplefor which (θ12 = 1.05◦) relaxes to decrease the total amount of AAA stacking (white) whencompared to the stacking order distribution in a rigid AtA trilayer (Fig. 2E). This is expectedas AAA stacking is roughly 29.5 and 36.5 meV/unit cell higher in energy than A-SP-A andABA stacking respectively (See Supplemental Information section 9). Further, the histogramshown (Fig. 2D) illustrates that this sample contains considerably more ABA, BAB, and SPtype stacking than AAA type stacking (See Supplemental figure 3 for the stacking histogramexpected of a rigid sample). This reconstruction also manifests in the mean line-cut shownin Fig. 2F, which corresponds to the average over all line-cuts equivalent to that denotedby the dotted line in Fig. 2C. From this profile, it is evident that the widths of the AAAregions are smaller than expected for a rigid trilayer, which is robust the presented standarddeviation (shaded region) and noise-driven differences in normalization (See SupplementalInformation section 7) as well as and limitations in spatial resolution due to appreciablebeam-width biasing (See Supplemental Information section 6).These results, as well as results for other twist angles shown in Supplemental Informationsection 8, together illustrate that AtA trilayers show an observably reconstructed atomic4stacking distribution up to a at least a twist angle of 1.81◦, with few differences seen betweenthe 1.81◦ and 1.0◦ samples. We also observe that the AtA trilayers show a pattern ofreconstruction similar to that of a twisted bilayer within this twist angle range [25, 26],although a more quantitative comparison between the bilayer and AtA trilayer reconstructionnecessitates more detailed intensity fitting approach to map strain tensor fields [26, 29] andwill be addressed in future work.A similar analysis is carried out for samples which we refer to as ‘tAB’, in which the bot-tom two graphene layers are aligned AB and the top layer is twisted (Fig. 2H), creatinga structure sometimes referred to as a ‘monolayer-on-bilayer graphene’, which also exhibitscorrelated electron physics [31, 32]. Unlike the AtA trilayers, these tAB structures showatomic reconstruction patterns driven by a preference to decrease the relative portion ofAAB stacking, as seen by comparing Fig. 2I to the rigid stacking order distribution (Fig.2K), the histogram in Fig. 2J, and the corresponding line-cut in Fig. 2L (see Supplemen-tal Information section 7 for normalization bias and error margins), which illustrates thetAB structure reconstructs such that the portions of AAB and BAB stacking are no longerequivalent as they would be in a rigid moiré. This decrease in AAB stacking is expected, itis 17.9 meV/unit cell higher in energy than ABC and BAB (See Supplemental Informationsection 9). We note here that this manifests in the AAB regions appearing to have a lowerpeak I0110 + I1010 + I1100 intensity to the BAB regions, while these regions are expected toappear sharper but with similar maximum intensity. This is likely due to broadening from anumber of factors associated with measurement acquisition and post-processing, especiallythe beam-width biasing and the Gaussian filter used (see Supplemental Information section6).The observation that AAB and BAB regions are observably distinct even at a twist angleof 1.4◦ is nonetheless notable. This effect is more dramatic at smaller twist angles as seenin the stacking order percent area trends and maps gathered within the 0.1◦ – 1.5◦ twistangle regime. While the approach used in this work prohibits a quantitative comparison ofthe AAB and BAB domain sizes (to be addressed in future work), these results still clearly5establish that the size of ABC and BAB domains are comparable to each-other and muchlarger than the AAB domains. Moreover, the stacking order distributions seen for tABappear similar to those observed in twisted bilayer graphene [26], suggesting that atomicreconstruction in tAB trilayer samples can be explained primarily from consideration of thetwisted interface.The results discussed thus far concern a limiting case of graphene trilayers wherein two of thelayers are perfectly aligned. While these materials are conceptually simpler and display richphysical properties which merit their investigation, this interferometric 4D-STEM techniquealso permits us to study a broader array of twisted trilayer structures with two independenttwist angles. In these more complex multilayered samples, the ability to selectively probeburied bilayer interfaces allows us to independently image the larger scale moiré pattern andevaluate its effect on local stacking order.Following this approach, we extract stacking order maps associated with the larger moirépattern from double overlap (Fig. 3A) and triple overlap (Fig. 3B,E). These are comparedto the maps calculated for rigid moirés in Fig. 3C. For the the double overlap case, Fig.3A reveals the presence of large local regions in which two layers are aligned directly atopeach-other (AtA or tAA, white) and regions in which two layers are aligned AB (AtB ortAB, blue). From comparing the stacking distributions of samples with 0 < θ13 ≪ θ23 (threeleftmost panels in Fig. 3A, illustrated in Fig. 3D) and 0 < θ23 ≪ θ13 (rightmost panel inFig. 3A, illustrated in Fig. 3F), we find that the two regimes display distinct reconstructionpatterns. When θ13 ≪ θ23, the observed atomic reconstruction is driven by a slight preferencefor AtA type stacking (white) over AtB (blue) and soliton-type (grey) regions. This resultis somewhat unexpected as the energetic difference between rigid AtA and AtB domains(driven only by inter-layer coupling between the top and bottom layers) has been previouslypresumed to play a minor role in reconstruction [33]. Moreover, previous STM studies [27]concluded that trilayers with 0 < θ13 ≪ θ23 relaxed to effectively eliminate AtB domains.Fig. 3A also shows (rightmost panel) that the atomic reconstruction pattern for θ23 ≪ θ136is instead driven by a preference to minimize the high energy tAA (white) domains, withinwhich every possible stacking configuration must place two carbon atoms from neighboringlayers directly atop each other — an arrangement that is sterically unfavorable [25, 26].The extent of reconstruction in these θ23 ≪ θ13 samples is therefore much larger, since theenergy difference between rigid tAA vs tAB domains (≈ 18.2 meV/unit when consideringonly adjacent interfaces) is much larger than that between rigid AtA vs AtB domains (≈0 meV/unit when considering only adjacent interfaces) [23, 33]. This is reflected in thedifference between the structures shown in the second and fourth columns of Fig. 3, inwhich both structures have comparable twist angles, but the structure in the fourth panel isobservably more reconstructed, with the tAA domains appearing as a highly contracted spot.This spot appears orange rather than white due to beam-width and data processing effects(see Supplemental Information section 6). Although the weaker higher frequency textureobserved within the white and blue domains in Fig. 3A might arise from the smaller scalemoiré pattern imparting a modulation in these stacking distributions, this pattern likelypredominantly results from a small bleed-in of the (1,2,3) interference pattern, which is hardto completely exclude with virtual apertures while retaining sufficient signal-to-noise ratios.After extracting the local AtA/tAA and AtB/tAB domains as shown in Fig. 3A, we nowexamine the (1,2,3) overlap region, which is associated with all three graphene layers (Figs.3B), to understand how the smaller scale moiré pattern modulates local stacking orderwithin these larger domains (representative regions of these maps are magnified in Fig. 3E).Additional maps are shown in Supplemental Figs. 7-8. We note that, as noted in previouswork [27], we see only two clear periodicities in our data despite the presence of threemoiré wavelengths from each twisted bilayer interface. However, this does not necessarilyimply that only two moiré wavelengths are present; inspection of the atomic stacking mapsexpected from even rigid structures (Fig. 3C and Supplemental Figures 7-8) reveals that thesmaller and larger periodicities observed reflect only the largest and smallest twist angles,respectively.Taken together, the data in Fig. 3 allow a quantification of the area fractions in TTLG7samples and the development of a qualitative model for reconstruction in the limit of slightmisalignment (Fig. 4). Fig. 4A shows that for the larger moiré pattern, the proportionof AtA/tAA and AtB/tAB stacking domains inverts across the regimes illustrated in Figs.3D,F. Associated area fractions from our measurements and those from continuum relaxationsimulations as a function of θ13− θ23 are shown in Fig. 4A. Experiment and simulation showgood agreement in the overall trends, though the measurements show a more gradual declinein area fraction of AtA/tAA (and corresponding rise in that of AtB/tAB/SP) than thesimulation with increasing θ13 − θ23. This slight discrepancy may arise because of kineticeffects preventing the system from realizing the theoretically optimal extent of relaxationdriven by layers not immediately adjacent.For the smaller scale moiré superlattice, we find that this pattern appears relatively invari-ant within the AtA, AtSP, AtB, and tAB domains (SI Section 12). Indeed, the measuredproportion of stacking orders within the AtA regions of Fig. 3 is very similar to the pureAtA sample seen in Fig. 2C, suggesting that the larger scale moiré plays a relatively minorrole in the reconstruction of the smaller scale moiré. We do however observe some differencesbetween the smaller-scale moiré within different domains. As shown in Fig. 4B, measure-ments of local θ12 values within the AtA and AtB domains of θ13 ≪ θ23 samples displaya slightly smaller θ12 angle within AtA regions as compared to the values in adjacent AtBdomains. This tightening of the smaller-scale moiré within AtB regions might help facilitatethe overall minimization of these AtB regions.Lastly, we investigate the maps of local atomic stacking order in regions with an increasinglylarge extent of extrinsic heterostrain, ϵ. From the maps shown in Fig. 5, we find that extrinsicheterostrain acts predominantly on the larger scale moiré pattern and has a diminishingeffect on the smaller scale moiré superlattice, consistent with previous work on bilayer moirésystems [26, 29]. Notably, in the most heterostrained sample of Fig. 5, despite similarθ13, the islands of AtA are deformed into stripes. These features have also been previouslyvisualized in STM studies and attributed to heterostrain between the top and bottom layers[23]. Heterostrain is therefore a powerful tuning knob for manipulating the contiguity of AtA8domains (from islands to stripes) at the expense of AtB regions, potentially modulating theemergence of correlated phases that rely on the C2zT -symmetric AtA domains.In conclusion, the nature of atomic reconstruction unveiled here for twisted trilayers ismarkedly different than that proposed in previous work, wherein it was suggested thatslightly misaligned MATTLG samples relax to almost exclusively AtA regions, with theAtB and SP regions stretched into thin domain boundaries and/or topological defects thatcontribute insignificantly to the STM measurements [27]. While Fig. 2 shows that at lengthscales where only one moiré wavelength is apparent (when θ13 ≈ 0◦), trilayers do favor theformation of large AtA domains, and in that case the local structure of trilayers is driven pri-marily by consideration of the smaller moiré, we see a clear presence of considerable AtB typestacking down to θ13 = 0.20◦ (Fig. 3). This observation contrasts previous claims of trilayersamples containing contributions from only AtA regions at a θ13 of ≈ 0.25◦, with furtherdiscussion of our results in the context of these prior measurements provided in SupplementalInformation section 10. Taken together, our measurements highlight the particular utility ofinterferometric 4D-STEM imaging alongside other scanning probe techniques like STM forcharacterizing complex multi-layered moirés, as the ability to apply a direct structural probeselectively to separate interfaces can uncover the complex picture of atomic reconstruction.The extent of AtB stacking observed could have major implications for understanding su-perconductivity in misaligned MATTLG [22, 23] and recently discovered moiré quasicrystalsystems [34]. For instance, if θ13 ≪ θ23 configurations such as MATTLG favored entirelyAtA stacking as previously proposed, their correlated behaviors could be predominantly un-derstood by consideration of the C2zT inversion symmetric AAA stacking regions much likeTBLG. While our measurements support a relative contraction of AtB domains, it is nowhere near as dramatic, revealing that sizable AtB portions remain following reconstruc-tion. These significant AtB domains may instead suggest that the ABC, AAB, and ABBregions, which have been shown to host correlated electronic phases [13, 35–37] despite alack of inversion symmetry, may play an important role in understanding correlated electronphysics in some twisted trilayers.9AcknowledgementsWe thank P. Kim, J. Ciston, K. Bustillo, and C. Ophus for vis-a-vis and epistolary discus-sions. This material is based upon work supported by the US National Science FoundationEarly Career Development Program (CAREER), under award no. 2238196 (D.K.B). I.M.C.acknowledges a pre-doctoral fellowship award under contract FA9550-21-F-0003 through theNational Defense Science and Engineering Graduate (NDSEG) Fellowship Program, spon-sored by the Air Force Research Laboratory (AFRL), the Office of Naval Research (ONR)and the Army Research Office (ARO). C.G. was supported by a grant from the W.M. KeckFoundation (Award no. 993922). Experimental work at the Molecular Foundry, LBNL wassupported by the Office of Science, Office of Basic Energy Sciences, the U.S. Departmentof Energy under Contract no. DE-AC02-05CH11231. Confocal Raman spectroscopy wassupported by a Defense University Research Instrumentation Program grant through theOffice of Naval Research under award no. N00014-20-1-2599 (D.K.B.). Other instrumenta-tion used in this work was supported by grants from the Canadian Institute for AdvancedResearch (CIFAR–Azrieli Global Scholar, Award no. GS21-011, D.K.B), the Gordon andBetty Moore Foundation EPiQS Initiative (Award no. 10637, D.K.B), and the 3M Founda-tion through the 3M Non-Tenured Faculty Award (no. 67507585, D.K.B). Theoretical workwas supported by the Theory of Materials FWP at Lawrence Berkeley National Laboratory,funded by the U.S. Department of Energy, Office of Science, Basic Energy Sciences, MaterialsSciences and Engineering Division, under Contract No. DE-AC02- 05CH11231 (S.G). K.W.and T.T. acknowledge support from JSPS KAKENHI (Grant Numbers 19H05790, 20H00354and 21H05233).Author ContributionsI.M.C., M.V., C.G., and D.K.B. conceived the study. M.V., C.G., K.Z., and N.D. fabricatedthe samples. M.V., C.G., and I.M.C. performed the experiments. I.M.C. and Z.Z. performedthe continuum simulations using code developed by Z.Z. I.M.C. developed and implementedthe interferometry code (with assistance from C.G.) and analyzed the data. T.T. and K.W.10provided the hBN crystals. D.K.B. and S.M.G. supervised the work. I.M.C. and D.K.B.wrote the manuscript with input from all co-authors.Competing InterestsThe authors declare no competing interests.Figure CaptionsFigure 1: Interferometric 4D-STEM dark field imaging of selected interfaces. (a)Schematic of the 4D-STEM approach, wherein beam interference is used to extract stackingorder. (b) Schematic illustrating the twist angle, θ, and layer numbering conventions usedto label the graphene trilayers. (c) Illustrations of various high-symmetry stacking config-urations realized within twisted trilayer moirés. (d,e) Average convergent beam electrondiffraction patterns for trilayers with θ13 ≈ 0◦ (d) and θ13 = 0.22◦(e). Overlapping TTLGBragg disks are highlighted in the insets. Attribution of each Bragg disk to a layer is moti-vated in Supplemental Information section 3. (f,h) Virtual dark field images correspondingto the overlap of layers 1 & 3. (g,i) Virtual dark field images corresponding to the overlapof all three layers. Scale bars are 1 nm−1 and 25 nm for reciprocal (d,e) and real space (f–i),respectively.Figure 2: Reconstruction in AtA and tAB trilayers. (b, h) schematic illustrating thelayer alignment in a AtA and a tAB trilayer. (a, g) Legends illustrating how color correlateswith the average first and second order Bragg disks intensities. Overlaid points are the inten-sities of high symmetry stacking orders obtained via multislice, see Methods. (c, i) Maps oflocal stacking order for AtA and tAB trilayers, cropped to exclude biasing from sample drift.(d, j) Histograms illustrating the relative prevalence of each stacking configuration. Notethat since the intensity does not depend linearly on the stacking order, a rigid bilayer willnot display a uniform distribution of intensities (see Supplemental Information sections 4–5).(e, k) Schematics illustrating the anticipated variation in local stacking order expected fora rigid structure, obtained via the expressions provided in Supplemental Information section115. The convention is used that AAA, ABA, SP, and SP* stacking denote interlayer offsetsof u12 = (0, 0), (a0/√3, 0), (a0(2√3)−1, 0), and (0, a0/2) respectively in Cartesian coordinateswhere a0 is the lattice constant. Similarly for the tAB samples, AAA, ABC, SP, and SP*stacking denote interlayer offsets of u13 = (0, 0), (a0/√3, 0), (a0(2√3)−1, 0), and (0, a0/2) re-spectively. (f, l) Intensity line-cuts corresponding to the average of all possible line-cutsequivalent to the dotted lines shown in (c, i) are given as solid lines. Shaded region repre-sents the standard deviation and arrows denote the statistically significant contractions ofAAA and AAB domains. Intensity variation expected for a rigid structure are given as dot-ted lines. Domain sizes are calculated from the full width at half max of I0110 + I1010 + I1100(red) as highlighted, where the value of I2110 + I1210 + I1120 (black) is used to distinguishbetween different high symmetry stacking orders (see Supplemental Information Section 6for validation of this approach). All scale bars are 25 nm.Figure 3: Atomic Stacking in slightly misaligned TTLG. (a) Maps of local atomicstacking from the larger moiré pattern only, corresponding to the local in-plane offset be-tween between layers 1 and 3 in panels 1-3, and the local in-plane offset between layers 2 and3 in panel 4. Colors shown correspond to the bi-variate colormap in (g), with accompanyingexpressions and simulations motivating the attribution of intensities to stacking orders aslabeled here in Supplemental Information sections 4–6. (b) Local atomic stacking obtainedfrom considering all three graphene layers. (c) Simulated stacking order maps for rigid moirésuperlattice analogues of b, obtained from the expression given in Supplemental Informa-tion Section 5. All scale bars are 25 nm. (d,f) Schematics of layer alignment in TTLGwith slightly misaligned layers. (e) Zoom-ins of the maps above illustrating the finer localmodulation of stacking order within (1) AtA, (2) AtSP, (3) AtB, and (4) tAB regions. (g)Legend illustrating how color correlates with the average first and second order Bragg disksintensities for both the two and three layer interference patterns, with labeled locations oftwo layer high symmetry regions. (h) The expected intensity variation for an individualbilayer interface within the general trilayer structure obtained via the expression providedin Supplemental Information section 4. (i) Variation in first and second order Bragg disks12intensities expected for a rigid twisted trilayer obtained via the expression provided in Sup-plemental Information section 5. Intensity relations are verified with multi-slice simulationsin Supplemental Information section 6.Figure 4: Reconstruction patterns and trends in TTLG. (a) Area fraction of atomicstacking domains from the larger moiré pattern only as a function of θ13−θ23. As both θ13 andθ23 are slightly variable for the samples discussed, additional plots of area fraction againstθ13 and θ23 independently are provided in Supplemental Information Fig. 10. Experimental(exp) data (corresponding to the maps shown in Fig. 3 and Supplemental Information Figs.7-8) are compared with relaxation simulations (sim). Area fractions associated with theexperimental and simulated data were obtained following the procedure described in themethods (with regions of I0110 + I1010 + I1100 > 0.5 and I2110 + I1210 + I1120 < 0.5 labeledAtB/tAB, I0110+I1010+I1100 > 0.5 and I2110+I1210+I1120 > 0.5 labeled AtA/tAA, and thoseremaining SP) and by applying a threshold of 0.25 degrees to the local curl respectively (withfurther details in Supplemental Information section 9). We note that these methods resultin functionally equivalent categorizations due to the small area and large intensity variationassociated with the soliton regions where these thresholds partition the data. Error barwidths are the twist angle standard deviation from 53-112 data points for each point. (b)Local twist angle associated with the smaller moiré within AtA and AtB domains. All pointscorrespond to the regime where θ13 ≪ θ23. Twist angle determination is described in themethods. (c) Qualitative schematic illustrating the atomic reconstruction patterns (largemoiré) observed for θ13 ≪ θ23 and θ23 ≪ θ13. Error bar widths are the twist angle standarddeviation from 3-33 data points for each point.Figure 5: Heterostrain Effects (a) Maps of the modulation in local stacking order be-tween layers 1 and 3 only for samples with an increasingly large percent of extrinsic het-erostrain. (b) Corresponding maps of the local stacking order modulation obtained whenconsidering all three graphene layers. Twist angles and percent heterostrains values andbounds were determined from fitting the size and asymmetry of the moiré triangles (seemethods). All scale bars are 25 nm.13References1. Balents, L., Dean, C. R., Efetov, D. K. & Young, A. F. Superconductivity and strongcorrelations in moiré flat bands. Nature Physics 16, 725–733 (2020).2. Cao, Y. et al. Correlated insulator behaviour at half-filling in magic-angle graphenesuperlattices. Nature 556, 80–84 (2018).3. Cao, Y. et al. Unconventional superconductivity in magic-angle graphene superlattices.Nature 556, 43–50 (2018).4. Lau, C. N., Bockrath, M. W., Mak, K. F. & Zhang, F. Reproducibility in the fabricationand physics of moiré materials. Nature 602, 41–50 (2022).5. Tian, H. et al. Evidence for Dirac flat band superconductivity enabled by quantumgeometry. Nature 614, 440–444 (2023).6. Bao, W. et al. Stacking-dependent band gap and quantum transport in trilayer graphene.Nature Physics 7, 948–952 (2011).7. Kumar, A. et al. Integer quantum Hall effect in trilayer graphene. Physical ReviewLetters 107, 126806 (2011).8. Lui, C. H., Li, Z., Mak, K. F., Cappelluti, E. & Heinz, T. F. Observation of an electri-cally tunable band gap in trilayer graphene. Nature Physics 7, 944–947 (2011).9. Yuan, S., Roldán, R. & Katsnelson, M. I. Landau level spectrum of ABA-and ABC-stacked trilayer graphene. Physical Review B 84, 125455 (2011).10. Craciun, M. et al. Trilayer graphene is a semimetal with a gate-tunable band overlap.Nature nanotechnology 4, 383–388 (2009).11. Chen, G. et al. Evidence of a gate-tunable Mott insulator in a trilayer graphene moirésuperlattice. Nature Physics 15, 237–241 (2019).12. Zhou, H. et al. Half-and quarter-metals in rhombohedral trilayer graphene. Nature 598,429–433 (2021).13. Zhou, H., Xie, T., Taniguchi, T., Watanabe, K. & Young, A. F. Superconductivity inrhombohedral trilayer graphene. Nature 598, 434–438 (2021).1414. Koshino, M. & McCann, E. Gate-induced interlayer asymmetry in ABA-stacked trilayergraphene. Physical Review B 79, 125443 (2009).15. Morell, E. S., Pacheco, M., Chico, L. & Brey, L. Electronic properties of twisted trilayergraphene. Physical Review B 87, 125414 (2013).16. Ahn, J., Park, S. & Yang, B.-J. Failure of Nielsen-Ninomiya theorem and fragile topol-ogy in two-dimensional systems with space-time inversion symmetry: application totwisted bilayer graphene at magic angle. Physical Review X 9, 021013 (2019).17. Song, Z. et al. All magic angles in twisted bilayer graphene are topological. Physicalreview letters 123, 036401 (2019).18. Zou, L., Po, H. C., Vishwanath, A. & Senthil, T. Band structure of twisted bilayergraphene: Emergent symmetries, commensurate approximants, and Wannier obstruc-tions. Physical Review B 98, 085435 (2018).19. Törmä, P., Peotta, S. & Bernevig, B. A. Superfluidity and quantum geometry in twistedmultilayer systems. arXiv preprint arXiv:2111.00807 (2021).20. Mora, C., Regnault, N. & Bernevig, B. A. Flatbands and perfect metal in trilayer moirégraphene. Physical review letters 123, 026402 (2019).21. Park, J. M., Cao, Y., Watanabe, K., Taniguchi, T. & Jarillo-Herrero, P. Tunablestrongly coupled superconductivity in magic-angle twisted trilayer graphene. Nature590, 249–255 (2021).22. Hao, Z. et al. Electric field–tunable superconductivity in alternating-twist magic-angletrilayer graphene. Science 371, 1133–1138 (2021).23. Kim, H. et al. Evidence for unconventional superconductivity in twisted trilayer graphene.Nature 606, 494–500 (2022).24. Zhang, X. et al. Correlated insulating states and transport signature of superconduc-tivity in twisted trilayer graphene superlattices. Physical review letters 127, 166802(2021).25. Yoo, H. et al. Atomic and electronic reconstruction at the van der Waals interface intwisted bilayer graphene. Nat. Mater. 18, 448–453 (2019).1526. Kazmierczak, N. P. et al. Strain fields in twisted bilayer graphene. Nature materials20, 956–963 (2021).27. Turkel, S. et al. Orderly disorder in magic-angle twisted trilayer graphene. Science 376,193–199 (2022).28. Ophus, C. Four-dimensional scanning transmission electron microscopy (4D-STEM):From scanning nanodiffraction to ptychography and beyond. Microscopy and Micro-analysis 25, 563–582 (2019).29. Van Winkle, M. et al. Quantitative Imaging of Intrinsic and Extrinsic Strain in Transi-tion Metal Dichalcogenide Moir\’e Bilayers. arXiv preprint arXiv:2212.07006 (2022).30. Zachman, M. J. et al. Interferometric 4D-STEM for lattice distortion and interlayerspacing measurements of bilayer and trilayer 2D materials. Small 17, 2100388 (2021).31. Tong, L.-H. et al. Spectroscopic visualization of flat bands in magic-angle twistedmonolayer-bilayer graphene: coexistence of localization and delocalization. Physical Re-view Letters 128, 126401 (2022).32. Park, Y., Chittari, B. L. & Jung, J. Gate-tunable topological flat bands in twistedmonolayer-bilayer graphene. Physical Review B 102, 035411 (2020).33. Zhu, Z., Carr, S., Massatt, D., Luskin, M. & Kaxiras, E. Twisted trilayer graphene: Aprecisely tunable platform for correlated electrons. Physical review letters 125, 116404(2020).34. Uri, A. et al. Superconductivity and strong interactions in a tunable moiréquasicrystal.Nature 620, 762–767 (2023).35. Chen, G. et al. Tunable correlated Chern insulator and ferromagnetism in a moirésuperlattice. Nature 579, 56–61 (2020).36. Xu, S. et al. Tunable van Hove singularities and correlated states in twisted monolayer–bilayer graphene. Nature Physics 17, 619–626 (2021).37. Li, S.-y. et al. Imaging topological and correlated insulating states in twisted monolayer-bilayer graphene. Nature Communications 13, 1–7 (2022).38. Kerelsky, A. et al. Maximized electron interactions at the magic angle in twisted bilayergraphene. Nature 572, 95–100 (2019).16MethodsSample PreparationAll graphene trilayers were fabricated using the tear-and-stack technique [2, 3]. Briefly, apolybisphenol-A-carbonate/polydimethylsiloxane (PC/PDMS) stamp was first used to pickup 5–10 nm hexagonal boron nitride (hBN) off a SiO2/Si substrate. This hBN was then usedto pick up and tear graphene monolayers and/or bilayers, also on SiO2/Si. The remaininggraphene portions were then sequentially rotated and picked up to construct the desiredtrilayer structure, which was stamped onto a Norcada TEM grid (200 nm silicon nitridewith 2 µm holes).Electron Microscopy MeasurementsAll electron microscopy measurements were performed at the National Center for ElectronMicroscopy at Lawrence Berkeley National Laboratory. Low-magnification dark-field TEMimages were collected using a Gatan UltraScan camera on a Thermo Fisher Scientific Titan-class microscope at 60 kV to identify regions of interest prior to 4D-STEM acquisition. 4D-STEM data were obtained using a Gatan K3 direct detection camera and Gatan Continuumimaging filter on the TEAM I microscope (aberration-corrected Thermo Fisher ScientificTitan 80–300). We operated in energy-filtered STEM mode at 80 kV using an 10-eV energyfilter centered around the zero-loss peak, convergence angle of 1.71 mrad, and a typical beamcurrent of 45–65 pA for an overall effective probe size of 1.25 nm corresponding to the full-width at half-maximum value. The diffraction patterns collected correspond to a step sizeof 0.5-2 nm depending on the sample. We operated the K3 camera in full-frame electroncounting mode with a binning of 4 × 4 pixels, energy-filtered STEM camera length of 800 mm,and exposure time of 13 ms which was the sum of multiple counted frames. Considerationsfor choice in convergence angle and camera length are discussed in Supplemental Informationsection 2.17Virtual Aperture Selection and Stacking Order MapsTo obtain the virtual dark fields shown in Figure 1 of the main text, we first had to isolatethe regions of reciprocal space associated with diffraction of (1,2,3) and (1,3) layers. Firstthe individual Bragg disks are attributed to layers as described in Supplemental Informationsection 3, using the known order in which each graphene layer was picked up, optical micro-graphs of the individual graphene layers prior to, and after assembly into the heterostructure,and conventional dark-field electron micrographs of the samples assembled on TEM grids.Virtual apertures where then obtained by using threshold intensities to draw contours inthe averaged diffraction patterns (see Supplemental Information section 1 for motivation),which were then used to obtain the masks associated with each overlap region. The intensitywithin these masked regions was then summed for each real space pixel to yield the virtualdark fields shown. In practice we used only the pixels close to the centers of each region(each region was down-sized by 25 percent) to minimize the bleed-in of intensity modula-tions originating from interference with other disks. Colored stacking order maps were thenconstructed from the virtual dark field images associated with the first and second orderBragg reflections using the provided bi-variate color-map, with red, blue, and green channelsequal to the average intensity in the first order disks, second order disks, and the average ofthe red and blue channel values, respectively. Attribution of intensities to stacking order isrationalized in Supplemental Information sections 4–6.Twist Angle and Heterostrain MeasurementTwist angles were determined through triangulating the high-symmetry stacking orders inaccordance with previous work. [26, 29]. For AtA samples, the bright AAA stacking regionswere identified by fitting the data to a series of Gaussians, the centers of which were tri-angulated using the Delaunay algorithm. The resultant moire wavelengths were then usedto determine the twist angle and percent of heterostrain through fitting to the expressionsprovided in [29, 38] using non-linear least squares. The twist angle associated with thelarger moiré was calculated similarly from the bright tAA/AtA regions within the two-layer18overlap virtual dark fields. For tAB samples, we instead inverted the virtual dark field andtriangulated the centers of the ABC stacking regions. For samples containing a large portionof heterostrain and/or an insufficiently large field of view, the smallest viable moiré trianglewas used to create a bound for the twist angle and percent heterostrain.Analysis of Stacking Area PortionsStacking area fractions were determined by thresholding the virtual dark field intensities intothe categories illustrated in Supplemental Information Figures 6 and 9. This correspondsto region definitions of I0110 + I1010 + I1100 > 0.5 and I2110 + I1210 + I1120 < 0.5 labeledAtB/tAB, I0110 + I1010 + I1100 > 0.5 and I2110 + I1210 + I1120 > 0.5 labeled AtA/tAA, andthose remaining SP. All statistical analyses and stacking order partitioning were performedafter smoothing the virtual dark field images with a Gaussian filter of σ = 0.5 nm. For theAtA and tAB samples shown, this partitioning was applied directly to the virtual dark fieldimage obtained from the three-disk overlap region. For the remaining samples, the geometricpartitioning was first performed on the the two-layer overlap virtual dark field micrographsto obtain masks for the stacking order in the larger moire pattern. These masks were thenapplied to the three-layer virtual dark fields to isolate the regions of tAB and tAA stackingso that they could be independently analyzed as in the case of AtA and tAB samples.Data AvailabilityThe data supporting the findings of this study are available within the Article and its Supple-mentary Information files. Datasets can be found at https://doi.org/10.5281/zenodo.4459669.Code AvailabilityThe code developed for data analysis in this study is available within the TrilayerTEMsub-directory at https://github.com/bediakolab/bediakolab scripts.19https://doi.org/10.5281/zenodo.4459669https://github.com/bediakolab/bediakolab_scripts