⇦Chinese Physics Letters, 2022, Vol. 39, No. 4, Article code 047403 Local Density of States Modulated by Strain in Marginally Twisted Bilayer Graphene Jia-Jun Ma (马佳俊)1,2†, Zhen-Yu Wang (王振宇)2†, Shui-Gang Xu (徐水钢)3,4, Yu-Xiang Gao (高于翔)1,2, Yu-Yang Zhang (张余洋)1,2,5, Qing Dai (戴庆)6, Xiao Lin (林晓)1,2,5, Shi-Xuan Du (杜世萱)1,2,5, Jindong Ren (任金东)6*, and Hong-Jun Gao (高鸿钧)1,2,5* Affiliations 1Institute of Physics, Chinese Academy of Sciences, Beijing 100190, China 2School of Physical Sciences, University of Chinese Academy of Sciences, Chinese Academy of Sciences, Beijing 100190, China 3Key Laboratory for Quantum Materials of Zhejiang Province, School of Science, Westlake University, Hangzhou 310024, China 4Institute of Natural Sciences, Westlake Institute for Advanced Study, Hangzhou 310024, China 5CAS Center for Excellence in Topological Quantum Computation, University of Chinese Academy of Sciences, Beijing 100190, China 6CAS Key Laboratory of Nanophotonic Materials and Devices, CAS Key Laboratory of Standardization and Measurement for Nano-technology, National Center for Nanoscience and Technology, Beijing 100190, China Received 16 February 2022; accepted 17 March 2022; published online 28 March 2022 †J. J. Ma and Z. Y. Wang contributed equally to this work.
*Corresponding authors. Email: renjd@nanoctr.cn; hjgao@iphy.ac.cn Citation Text: Ma J J, Wang Z Y, Xu S G et al. 2022 Chin. Phys. Lett. 39 047403    Abstract In marginally twisted bilayer graphene, the Moiré pattern consists of the maximized AB (BA) stacking regions, minimized AA stacking regions and triangular networks of domain walls. Here we realize the strain-modulated electronic structures of marginally twisted bilayer graphene by scanning tunneling microscopy/spectroscopy and density functional theory (DFT) calculations. The experimental data show four peaks near the Fermi energy at the AA regions. DFT calculations indicate that the two new peaks closer to the Fermi level may originate from the intrinsic heterostrain and the electric field implemented by back gate is likely to account for the observed shift of the four peaks. Furthermore, the $dI/dV$ map across Moiré patterns with different strain strengths exhibits a distinct appearance of the helical edge states.
cpl-39-4-047403-fig1.png
cpl-39-4-047403-fig2.png
cpl-39-4-047403-fig3.png
cpl-39-4-047403-fig4.png
 DOI:10.1088/0256-307X/39/4/047403 © 2022 Chinese Physics Society Article Text The magic angle twisted bilayer graphene (MATBG) system, the recently emerged van der Waals structure, has attracted much interest.[1–6] The bilayer graphene shows two electronic states near the Fermi level. They are the flat bands generated from the Moiré of the MATBG, where charge carriers do not attain enough kinetic energy.[2,3] Therefore, their low energy physics would be mainly determined by the Coulomb interaction.[3,4] Through electrostatic doping from back gate, the flat bands could move closer to the $E_{\rm F}$, leading to the emergence of some electronic correlation behaviors, such as superconductivity,[7] Mott metal-insulator transition[8] and magnetism.[9] In contrast to the MATBG, the marginally ($ < $ $1.1^{\circ}$) twisted bilayer graphene (TBG) also hosts much interesting physics.[10] Huang et al.[11] observed a double line profile of the topologically protected helical edge states around the domain walls when the AB and BA stacking areas are gapped out. Verbakel et al.[12] and Yao et al.[13] reported the existence of the network in the scanning tunneling microscopy (STM) image at room temperature, and Verbakel et al.[12] proved that the transport along this one-dimensional network is valley protected. The ununiformity of the Moiré superlattice in the STM images of the TBG systems can be explained by the heterostrain introduced in the preparation process.[3] Indeed, strain has been proved to be an important parameter to tune both the topographies and electronic properties of TBG.[14] For example, it could induce the pseudo-magnetic fields in the van der Waals structure, thus leading to the appearance of the pseudo-Landau levels.[15] Moreover, as reported by Edelberg et al.,[16] the strain engineering in the transition-metal-dichalcogenide homobilayer could create soliton networks, which could act as deep trapping potentials. Although strain-induced properties have been widely discussed in the MATBG, the research focusing on the marginally TBG with strain is still missing. In this work, we investigate the electronic properties, bias-dependent topographies and the helical edge states of strain-modulated marginally TBG by density functional theory (DFT) calculations and low-temperature STM/STS (scanning tunneling spectroscopy) operated at 5 K.[17–19] It is found that at AA stacking regions, four peaks labeled as V2, V1, C1 and C2 are clearly identified in the STS. DFT calculations in the atomic model of TBG with 6$^{\circ}$ indicate that the V1 and C1 peaks may originate from the strain. The calculated peak positions of these four peaks change with the vertical electric field, consistent with the peak positions' evolution as a function of the gate voltages in STS. For another marginally TBG, STM analysis reveals that there is no bias polarity dependence of the topographic difference within the Moiré pattern. Occasionally, a domain boundary is observed. The $dI/dV$ map across it with a relatively high electric field shows the presence of the helical edge states in the smaller strain Moiré pattern, while the absence of them in the other Moiré pattern indicates that the large strain may change the band structure, thus making the edge states unable to be detected.
cpl-39-4-047403-fig1.png
Fig. 1. Topography and spectroscopy of the marginally TBG. (a) Schematic diagram of STM/STS measurement of the TBG sample on h-BN/SiO$_{2}$ with back-gate voltage. The right part of (a) is the schematic diagram of the Moiré pattern caused by lattice mismatch between two layers. (b) STM image of the TBG with a twist angle of 0.09$^{\circ}$ ($V_{\rm bias} = -1.0$ V, $I_{\rm setpoint } = 100$ pA). The right part of (c) is the schematic of the real-space Moiré pattern evolving between AA, AB, BA and saddle-points. (c) Optical image of the sample. White dashed area reveals the TBG area on h-BN. (d) The $dI/dV$ spectrum taken on the AA region indicated in (c) by white dashed circle, showing four peaks around the Fermi level ($V_{\rm bias} = -0.2$ V, $I_{\rm setpoint} = 100$ pA, modulation 5 mV/714.02 Hz).
The sample in this study was fabricated by employing a “cut & stack” technique developed recently with the assist of polypropylene carbonate (PPC)/polydimethylsiloxane (PDMS) stamp. A thick h-BN and monolayer graphene were exfoliated on SiO$_{2}$/Si substrates. The monolayer graphene was cut into two parts in advance. The h-BN flake was firstly picked up by PPC. Then, h-BN was used to pick up the first part of the graphene. After the transfer stage was rotated by approximately 0.1$^{\circ}$, the second part of the graphene was picked up. In the above process, the temperature was kept at 40 ℃. To flip the stack, the PPC film with h-BN/Gr/Gr was transferred to a clean PDMS, with the stack side touching with PDMS. The PPC film can be dissolved by acetone. Finally, the stack was transferred onto pre-patterned gold electrodes. Another graphite flake was used to connect the graphene and gold electrodes. The schematic diagram of the STM experimental setup is demonstrated in Fig. 1(a). In our configuration, the tip is grounded and the sample bias ($V_{\rm bias}$) is applied across the TBG. The back-gate voltage is applied on the p$^{++}$ doped Si below the 300 nm-thick silicon dioxide, acting as the global back gate. Figure 1(b) shows a typical STM image of the TBG. The Moiré pattern caused by lattice mismatch between twisted layers can be easily found. However, the pattern sizes along different orientations are diverse, which is believed to be caused by strain coming from the fabrication processes. Figure 1(b) shows the Moiré pattern where the isosceles triangle is given by $L_{1} = 163$ nm, $L_{2} = 143$ nm and $L_{3} = 145$ nm. A large value of 20 nm can be inferred from difference between $L_{1}$ and $L_{2}$ (or $L_{3}$). As reported in the previous work,[3] the effect of homostrain on Moiré wavelength, given sub-percent strains in the experiment, can be neglected, while the heterostrain could have significant influence on the Moiré pattern periodicity. Therefore, heterostrain will be used to explain the following experimental observations. Indeed, the strength of heterostrain at $\sim $0.02% can be easily obtained from the well-known uniaxial strain model[3] considering the heterostrain to be a uniaxial strain applied on one layer of graphene and the other layer of graphene being twisted, which means that tiny heterostrain can give rise to huge variation in superlattice periodicities in marginally TBG. Furthermore, the twist angle at $\sim$$0.09^{\circ}$ was also found using the same model. At first glance, the Moiré pattern sizes in Fig. 1(b) are diverse. However, the subsequently calculated twist angle $\theta$ and strain $\varepsilon$ remain approximately fixed, i.e., $\theta \sim 0.09^{\circ}$, $\varepsilon \sim 0.02{\%}$, in the whole area. The optical image of the TBG is shown in Fig. 1(c). The pink area in the image corresponds to the bare SiO$_{2}$ surface. The white dashed cyan area is the TBG on the h-BN substrate. Since the SiO$_{2}$ surface is too rough to figure out the Moiré pattern, we choose the flat h-BN/SiO$_{2}$ as the substrate. The thick h-BN layers appear as a bright contrast in the optical image and make it easy to distinguish the scanned device. However, the exact position of the graphene-covered areas cannot be directly identified in the optical image since the monolayer or bilayer graphene is too thin to show any contrast difference. Indeed, their relative position has been marked during the sample preparation, as shown in Fig. 1(c). The Moiré pattern in the TBG manifests local stacking orders in different regions,[1–3] see inset of Fig. 1(b). The AA points, marked by the light blue circles, correspond to the positions where the top graphene hexagonal lattices are nearly exactly stacked onto the bottom graphene hexagonal lattices. The yellow strips connecting AA regions are the saddle points, also called topological domain walls in other works.[3,20] They separate the adjacent two triangular regions denoted as AB and BA, where the two graphene layers are nearly Bernal stacked and the only difference between them is the position of the top layer graphene's sublattice relative to that of the bottom layer. The area of AA and saddle-point regions are much smaller in marginally TBG, mainly determined by the structural reconstruction.[21] For example, the right part of Fig. 1(a) is the Moiré pattern formed in a 2.2$^{\circ}$ TBG. Note that the area of AA region is comparable to those of AB and BA regions. Considering the twist angle smaller than 0.5$^{\circ}$, the van der Waals interlayer interaction will maximize AB and BA regions, leading to much smaller area of AA regions and saddle-point regions since the van der Waals interlayer interaction could favor the interlayer commensurability that will compete with intralayer lattice distortion resulting in the large structure reconstruction.[21] To figure out the local density of states within AA stacking regions, we performed STS measurements on the region marked by the white dashed circle in Fig. 1(c). Figure 1(d) exhibits the observed $dI/dV$ spectrum when the back-gate voltage $V_{\rm g}$ is zero. There exist four peaks in the bias range of $\pm 0.2$ V, labeled as V2, V1, C1 and C2, according to their relative positions to the Fermi level. The electronic structure of the TBG can be easily tuned and examined through the STS measurements when varying the gate voltage. In the MATBG system, the back gate usually plays a role of tuning the carrier concentration in the sample and providing the electric field vertically across the sample in both STM and transport measurements.[2–4,7,8] This will allow researchers to study the electronic properties of this quasi two-dimension system. For instance, with varying the gate voltage, the correlated peaks could move across the Fermi level displaying the features of flat bands in MATBG.[1–3] The spectra collected at the AA regions as a function of the back-gate voltage are plotted in Fig. 2(a), where the red color represents high density of states (DOS) and the blue color represents low DOS. The spectrum with back-gate voltage less than $-40$ V is occasionally unstable and should be caused by the random noise of the system. Such an effect leads to the relatively high peak intensity in the top section of Fig. 2(a). To further analyze these spectra, $dI/dV$ spectra with obvious features are selected from Fig. 2(a) and plotted in Fig. 2(b). The related gate voltages are labeled at the right sites of each spectrum.
cpl-39-4-047403-fig2.png
Fig. 2. Gate voltage-dependent spectroscopy of the AA region in the marginally TBG. (a) The $dI/dV$ spectra on the AA region as a function of the back-gate voltage ($V_{\rm bias} = -0.2$ V, $I_{\rm setpoint} = 100$ pA, modulation 5 mV/714.02 Hz). (b) The $dI/dV$ spectra at different gate voltages selected from (a). The curves are shifted for clarity. Each gate voltage value is marked at the right site. (c) DOS of the TBG without strain (black curve) and with strain (red curve) by DFT calculations in the atomic model of TBG with 6$^{\circ}$. (d) The four peak positions' shift near the Fermi level at different electric fields based on DFT calculations in the atomic model of TBG with 6$^{\circ}$.
As shown in Figs. 2(a) and 2(b), the evolution of the four peaks in Fig. 1(d) as a function of the back-gate voltage can be clearly identified and marked. The band structure of the marginally TBG undertakes obvious changes by varying the gate voltage in our device. Note that when the back-gate voltage increases, the C2 and V2 peak positions tend to move parallel so that the energy gap between C2 and V2 peaks keeps constant. We attribute this phenomenon to the electron-hole symmetry of the TBG system.[22] When the back-gate voltage is around 0 V, at the charge neutral point, there are large energy gaps between V1 and V2 peaks ($\varDelta_{1}$), C1 and C2 peaks ($\varDelta_{2}$). Once increasing the magnitude of the back-gate voltage, both energy gaps $\varDelta_{1}$ and $\varDelta_{2}$ reduce. If the back-gate voltage is large enough, the two energy gaps are almost eliminated, i.e., the V1 and V2 peaks, C1 and C2 peaks tend to merge with each other, as shown in the top and bottom sections of Fig. 2(a). The two gaps reduce in the symmetric manner towards positive and negative directions. What's more, when the peaks move across the Fermi level in the spectra, there is no peak splitting and no appearance of the new peaks near the Fermi level, which is the signature of the correlated states in the TBG system.[2,3,7] Notice that there are four peaks around the Fermi level in the STS. This is contradictory to the previous work,[13] where only two van Hove singularities around the Fermi level were reported. The only difference between our sample and the reported one is the existence of new larger strain in TBG. To find out whether the peaks originate from the strain, we adopt the DFT calculations to compare the electronic structures of TBG both with and without strain, as shown in Fig. 2(c). Due to the computational cost, we use small TBG models containing around 400 atoms with a large twist angle of 6$^{\circ}$. The van Hove singularities are robust in a broad range of twisted angles.[23] All the DFT calculations are performed by using the Vienna ab initio simulation package (VASP)[24,25] with the projector-augmented-wave method.[26] The Perdew–Burke–Ernzerhof[27] type of exchange correlation functional and DFT-D3[28] correction were used. The kinetic energy cutoff is 400 eV. The $k$-points sampling was $3 \times 3 \times 1$ for self-consistent calculation and $18 \times 18 \times 1$ for DOS calculations, generated by Monkhorst–Pack grids with the origin at the $\varGamma$ point. The calculated DOS of TBG without strain shows two peaks near the Fermi level, which coincides well with the previous report.[13] However, two new peaks appear closer to the Fermi level when considering the strain of 5%. This result illustrates that the strain has a great chance to be the genuine origin of the observed two new peaks.
cpl-39-4-047403-fig3.png
Fig. 3. STM/STS measurements of 0.34$^{\circ}$ TBG sample. (a) STM morphology of TBG ($\theta = 0.34^{\circ}$) under different sample bias voltages and fixed tunneling current ($I_{\rm setpoint } = 100$ pA). (b) Line profiles of the STM morphology as indicated by the black arrows in (a). (c) STS measured on AA site ($V_{\rm bias} = -0.2$ V, $I_{\rm setpoint} = 100$ pA, modulation 5 mV/714.02 Hz), which is in accordance with the spectrum in Fig. 1(d) ($V_{\rm bias} = -0.2$ V, $I_{\rm setpoint} = 100$ pA, modulation 5 mV/714.02 Hz).
In order to explore the exact reason for the STS peak shift in Fig. 2(a), we calculate the DOS under different vertical electric fields, as shown in Fig. 2(d). The electric field with the order of 0.0055 V/nm may be generated while 1.0 V gate voltage is applied.[22] We calculate the DOS with an electric field ranging from $-0.3$ V/nm to 0.3 V/nm, considering the largest gate voltage 50 V in this experiment. It is found that the four peaks near the Fermi level could be moved by the electric field. The shift of the V2 and C2 peaks has a linear relationship to the strength of the electric field. On the other hand, for V1 and C1 peaks, the relationship is not linear. The energy difference between V1 and C1 peaks as a function of electric field undergoes a special variation tendency. When the electric field gradually changes from $-0.3$ V/nm to 0.3 V/nm, the energy difference becomes smaller at the beginning and starts to increase at around 0 V/nm. This is consistent with our gate-dependent spectra. In consequence, the vertical electric field may lead to the shift of the four peaks around the Fermi level in the marginally TBG. Another TBG using the same prepared method is also investigated, as shown in Fig. 3(a). The relatively uniform periodicity can be easily noticed in the STM topographies. The lengths of the Moiré pattern in the three different directions are $L_{1}= 45$ nm, $L_{2}= 45$ nm, $L_{3}= 32$ nm. The twist angle at 0.34$^{\circ}$ and the strain strength around 0.28% are also extracted from the uniaxial strain model.[3] Both the values are larger than those of the former sample. What's more, some different features are also exhibited in Fig. 3(a). With fixing the tunneling current at 100 pA, a clear contrast between the AB and BA stacking regions at both positive and negative biases can be seen. Such contrast difference has been observed in other small angle TBG.[29] In general, the topographies are only determined by the bias value rather than polarity, e.g., the contrast difference is relatively small when setting bias at 1.0 V and $-1.0$ V in Fig. 3(a). Therefore, we analyze the topographies scanned by positive bias in Fig. 3(a). Figure 3(b) reveals the observed height evolution in the STM images indicated by the black arrows as a function of the bias voltage. When reducing bias values, the contrast between AB and BA increases gradually. However, the contrast largely changes from 0.2 Å to 2 Å when further reducing bias to 0.2 V. At the same time, a three-fold symmetry feature can also be identified in the AA region when the bias is lower than 0.2 V. Similar to Fig. 1(d), the $dI/dV$ spectrum is also performed in AA region in this 0.34$^{\circ}$ TBG. Previous four peaks were also found in the range of $\pm 0.2$ V, which means that the electronic properties of marginally TBG systems should be analogous. The contrast difference helps us identify a domain boundary between two different periodicity Moiré patterns, namely $\alpha$ and $\beta$ regions, as seen in Fig. 4(a). The boundary is marked by the white dashed line. The sizes of Moiré patterns in $\alpha$ and $\beta$ regions are given by $L_{1\alpha }= 53$ nm, $L_{2\alpha }= 58$ nm, $L_{3\alpha }= 53$ nm, $L_{1\beta }= 42$ nm, $L_{2\beta }= 49$ nm, and $L_{3\beta }= 54$ nm. Thus, the twist angles and strain strengths in each domain can be calculated as $\theta_{\alpha }=0.25\mathrm{^{\circ} }$, $\varepsilon_{\alpha }=0.05\%$; ${\theta }_{\beta}=0.30\mathrm{^{\circ}}$, ${ \varepsilon }_{\beta }=0.13\%$. The strain of $\alpha$ area is smaller than that of $\beta$ area. As discussed above, the domain walls in the marginally TBG originate from the surface reconstruction and separate the AB and BA stacking regions. Researchers have shown that the energy band of the Bernal stacking TBG will exhibit an energy gap when its inversion symmetry is broken by the electric field.[30,31] The sign of valley Chern number of the bilayer graphene, which is calculated from the integral of the Berry curvature over the valley, depends on the polarity of electric field and the stacking order.[11] Therefore, in the Bernal stacking bilayer graphene with uniform configuration, inverse electric fields will lead to the opposite signs of the valley Chern number. For our TBG with AB and BA stacking regions in the same electric field, the valley Chern numbers will also exhibit the opposite signs. In this case, the domain walls between these regions with opposite valley Chern numbers could support the helical edge states.[11] Hence, if we apply a large enough gate voltage to our TBG system, the energy gap should be opened, leading to the appearance of edge states in the domain wall, see Fig. 4(b). Clear helical edge states can be distinguished in the domain walls of the $\alpha$ region, as shown in the $dI/dV$ map in Fig. 4(d). The STM topography is also listed in Fig. 4(c) for reference. We also find that no edge states can be identified in the $\beta$ region regardless of the relatively large electric field. Between these two regions are the tiny twist angle and strain differences. Further, from Fig. 3 we can conclude that the electronic properties of the marginally TBG systems should be analogous. Therefore, we propose it is the large strain in $\beta$ region that changes the band structure,[32] leading to the absence of edge states[11] in $dI/dV$ map.
cpl-39-4-047403-fig4.png
Fig. 4. Two domains with different twist angles ($\alpha$ region: $\theta_{\alpha }= 0.25^{\circ}$, $\beta$ region: $\theta_{\beta }= 0.30^{\circ}$) separated by a domain boundary. (a) STM image of the boundary. Both sides are marginally TBG. The contrast difference between AB and BA regions can be seen. (b) Schematic electronic band structure of the twisted bilayer graphene with large symmetry-breaking electric field. [(c), (d)] STM image (c) and simultaneously acquired $dI/dV$ map (d) across the domain boundary ($V_{\rm bias} = -0.2$ V, $I_{\rm setpoint} = 100$ pA). The map voltage is set as $-57$ mV with back-gate voltage at $+$50 V. The two domain areas exhibit totally different features.
In summary, we have performed gate-tunable STM/STS on the marginally TBG. Four peaks near the Fermi level are identified in the STS at the center of the AA stacking regions. DFT calculations indicate that the V1 and C1 peaks are likely to originate from the strain which was unavoidable during the preparation process and the shift of these four peaks as a function of the back-gate voltage may derive from the vertical electric field. In another tiny twist angle sample, a clear contrast between AB and BA stacking regions is exhibited, which helps to identify a domain boundary between two different Moiré patterns. The $dI/dV$ map with electric field demonstrates the presence of the helical edge states. The large strain in the graphene lattice may change the band structure of marginally TBG, hence making the edge states unable to be detected in $dI/dV$ map. Our results indicate that the heterostrain caused by sample preparation will tune the electronic states in the TBG system, which is useful in designing strain-based van der Waals materials. Acknowledgments. This work was supported by the National Natural Science Foundation of China (Grant Nos. 61888102 and U2032206), the National Key Research and Development Program of China (Grant Nos. 2018YFA0305800 and 2019YFA0308500), Chinese Academy of Sciences (Grant Nos. XDB30000000, XDB36000000, YSBR-003, and 112111KYSB20160061), and the Fundamental Research Funds for the Central Universities, China. References Electronic correlations in twisted bilayer graphene near the magic angleCharge order and broken rotational symmetry in magic-angle twisted bilayer grapheneMaximized electron interactions at the magic angle in twisted bilayer grapheneSpectroscopic signatures of many-body correlations in magic-angle twisted bilayer grapheneGraphene bilayers with a twistMomentum Space Quantum Monte Carlo on Twisted Bilayer GrapheneUnconventional superconductivity in magic-angle graphene superlatticesCorrelated insulator behaviour at half-filling in magic-angle graphene superlatticesEmergent ferromagnetism near three-quarters filling in twisted bilayer grapheneEmergence of Chern Insulating States in Non-Magic Angle Twisted Bilayer GrapheneTopologically Protected Helical States in Minimally Twisted Bilayer GrapheneValley-protected one-dimensional states in small-angle twisted bilayer grapheneSingularities and topologically protected states in twisted bilayer grapheneFano Resonance Enabled Infrared Nano-Imaging of Local Strain in Bilayer GrapheneEvidence of flat bands and correlated states in buckled graphene superlatticesTunable strain soliton networks confine electrons in van der Waals materialsKondo Effect of Cobalt Adatoms on a Graphene Monolayer Controlled by Substrate-Induced RipplesInteratomic Spin Coupling in Manganese Clusters Registered on GrapheneExperimental Synthesis of Strained Monolayer Silver Arsenide on Ag(111) SubstratesGiant oscillations in a triangular network of one-dimensional states in marginally twisted grapheneAtomic and electronic reconstruction at the van der Waals interface in twisted bilayer grapheneSpectroscopy of a tunable moiré system with a correlated and topological flat bandUnraveling the Intrinsic and Robust Nature of van Hove Singularities in Twisted Bilayer Graphene by Scanning Tunneling Microscopy and Theoretical AnalysisEfficiency of ab-initio total energy calculations for metals and semiconductors using a plane-wave basis setEfficient iterative schemes for ab initio total-energy calculations using a plane-wave basis setProjector augmented-wave methodGeneralized Gradient Approximation Made SimpleA consistent and accurate ab initio parametrization of density functional dispersion correction (DFT-D) for the 94 elements H-PuLarge-area, periodic, and tunable intrinsic pseudo-magnetic fields in low-angle twisted bilayer grapheneHelical networks in twisted bilayer graphene under interlayer biasValley Chern numbers and boundary modes in gapped bilayer grapheneBand structure engineering of graphene by strain: First-principles calculations
[1] Choi Y, Kemmer J, Peng Y, Thomson A, Arora H, Polski R, Zhang Y, Ren H, Alicea J, Refael G, Von Oppen F, Watanabe K, Taniguchi T, and Nadj-Perge S 2019 Nat. Phys. 15 1174
[2] Jiang Y, Lai X, Watanabe K, Taniguchi T, Haule K, Mao J, and Andrei E Y 2019 Nature 573 91
[3] Kerelsky A, Mcgilly L J, Kennes D M, Xian L, Yankowitz M, Chen S, Watanabe K, Taniguchi T, Hone J, Dean C, Rubio A, and Pasupathy A N 2019 Nature 572 95
[4] Xie Y, Lian B, Jack B, Liu X, Chiu C L, Watanabe K, Taniguchi T, Bernevig B A, and Yazdani A 2019 Nature 572 101
[5] Andrei E Y and Macdonald A H 2020 Nat. Mater. 19 1265
[6] Zhang X, Pan G, Zhang Y, Kang J, and Meng Z Y 2021 Chin. Phys. Lett. 38 077305
[7] Cao Y, Fatemi V, Fang S, Watanabe K, Taniguchi T, Kaxiras E, and Jarillo-Herrero P 2018 Nature 556 43
[8] Cao Y, Fatemi V, Demir A, Fang S, Tomarken S L, Luo J Y, Sanchez-Yamagishi J D, Watanabe K, Taniguchi T, Kaxiras E, Ashoori R C, and Jarillo-Herrero P 2018 Nature 556 80
[9] Sharpe A L, Fox E J, Barnard A W, Finney J, Watanabe K, Taniguchi T, Kastner M A, and Goldhaber-Gordon D 2019 Science 365 605
[10] Shen C, Ying J, Liu L, Liu J, Li N, Wang S, Tang J, Zhao Y, Chu Y, Watanabe K, Taniguchi T, Yang R, Shi D, Qu F, Lu L, Yang W, and Zhang G 2021 Chin. Phys. Lett. 38 047301
[11] Huang S, Kim K, Efimkin D K, Lovorn T, Taniguchi T, Watanabe K, Macdonald A H, Tutuc E, and Leroy B J 2018 Phys. Rev. Lett. 121 037702
[12] Verbakel J D, Yao Q, Sotthewes K, and Zandvliet H J W 2021 Phys. Rev. B 103 165134
[13] Yao Q, Chen X, Van Bremen R, Sotthewes K, and Zandvliet H J W 2020 Appl. Phys. Lett. 116 011602
[14] Du J, Lyu B, Shan W, Chen J, Zhou X, Xie J, Deng A, Hu C, Liang Q, Xie G, Li X, Luo W, and Shi Z 2021 Chin. Phys. Lett. 38 056301
[15] Mao J, Milovanovic S P, Andelkovic M, Lai X, Cao Y, Watanabe K, Taniguchi T, Covaci L, Peeters F M, Geim A K, Jiang Y, and Andrei E Y 2020 Nature 584 215
[16] Edelberg D, Kumar H, Shenoy V, Ochoa H, and Pasupathy A N 2020 Nat. Phys. 16 1097
[17] Ren J, Guo H, Pan J, Zhang Y Y, Wu X, Luo H G, Du S, Pantelides S T, and Gao H J 2014 Nano Lett. 14 4011
[18] Ren J, Guo H, Pan J, Zhang Y F, Yang Y, Wu X, Du S, Ouyang M, and Gao H J 2017 Phys. Rev. Lett. 119 176806
[19] Zhang S, Song Y, Li H, Li J M, Qian K, Liu C, Wang J O, Qian T, Zhang Y Y, Lu J C, Ding H, Lin X, Pan J, Du S X, and Gao H J 2020 Chin. Phys. Lett. 37 068103
[20] Xu S G, Berdyugin A I, Kumaravadivel P, Guinea F, Krishna K R, Bandurin D A, Morozov S V, Kuang W, Tsim B, Liu S, Edgar J H, Grigorieva I V, Fal'ko V I, Kim M, and Geim A K 2019 Nat. Commun. 10 4008
[21] Yoo H, Engelke R, Carr S, Fang S, Zhang K, Cazeaux P, Sung S H, Hovden R, Tsen A W, Taniguchi T, Watanabe K, Yi G C, Kim M, Luskin M, Tadmor E B, Kaxiras E, and Kim P 2019 Nat. Mater. 18 448
[22] Liu X, Chiu C L, Lee J Y, Farahi G, Watanabe K, Taniguchi T, Vishwanath A, and Yazdani A 2021 Nat. Commun. 12 2732
[23] Brihuega I, Mallet P, Gonzalez-Herrero H, De Trambly L G, Ugeda M M, Magaud L, Gomez-Rodriguez J M, Yndurain F, and Veuillen J Y 2012 Phys. Rev. Lett. 109 196802
[24] Kresse G and Furthmuller J 1996 Comput. Mater. Sci. 6 15
[25] Kresse G and Furthmuller J 1996 Phys. Rev. B 54 11169
[26] Blochl P E 1994 Phys. Rev. B 50 17953
[27] Perdew J P, Burke K, and Ernzerhof M 1996 Phys. Rev. Lett. 77 3865
[28] Grimme S, Antony J, Ehrlich S, and Krieg H 2010 J. Chem. Phys. 132 154104
[29] Shi H, Zhan Z, Qi Z, Huang K, Veen E V, Silva-Guillen J A, Zhang R, Li P, Xie K, Ji H, Katsnelson M I, Yuan S, Qin S, and Zhang Z 2020 Nat. Commun. 11 371
[30] San-Jose P and Prada E 2013 Phys. Rev. B 88 121408
[31] Zhang F, Macdonald A H, and Mele E J 2013 Proc. Natl. Acad. Sci. USA 110 10546
[32] Gui G, Li J, and Zhong J 2008 Phys. Rev. B 78 075435